RGB->YUV

28.4.2012 20:55 | Přečteno: 1676× | Jiné |

Vrámci prográmku zmiňovaném v minulém zápisu potřebuju konverzi RGB na YUV. V GLC bylo toto:

#define RGB_TO_YCbCrJPEG_Y(Rd, Gd, Bd) \
	(    + ((306 * (Rd) + 601 * (Gd) + 117 * (Bd)) >> 10))
#define RGB_TO_YCbCrJPEG_Cb(Rd, Gd, Bd) \
	(128 - ((173 * (Rd) + 339 * (Gd) - 512 * (Bd)) >> 10))
#define RGB_TO_YCbCrJPEG_Cr(Rd, Gd, Bd) \
	(128 + ((512 * (Rd) - 429 * (Gd) -  83 * (Bd)) >> 10))

void ycbcr_bgr_to_jpeg420(save_thread_t *th, unsigned char *from)
{
	unsigned int Ypix;
	unsigned int op1, op2, op3, op4;
	unsigned char Rd, Gd, Bd;
	unsigned int ox, oy, Yy, Yx;
	unsigned char *Y, *Cb, *Cr;

	Y = th->yuv_buffer;
	Cb = &th->yuv_buffer[th->yw * th->yh];
	Cr = &th->yuv_buffer[th->yw * th->yh + th->cw * th->ch];

	oy = (th->yh - 2) * th->stride;
	ox = 0;

	for (Yy = 0; Yy < th->yh; Yy += 2) {
		for (Yx = 0; Yx < th->yw; Yx += 2) {
			op1 = ox + oy;
			op2 = op1 + 4;
			op3 = op1 + th->stride;
			op4 = op2 + th->stride;
			Rd = (from[op1 + 2] + from[op2 + 2] + from[op3 + 2] + from[op4 + 2]) >> 2;
			Gd = (from[op1 + 1] + from[op2 + 1] + from[op3 + 1] + from[op4 + 1]) >> 2;
			Bd = (from[op1 + 0] + from[op2 + 0] + from[op3 + 0] + from[op4 + 0]) >> 2;

			/* CbCr */
			*Cb++ = RGB_TO_YCbCrJPEG_Cb(Rd, Gd, Bd);
			*Cr++ = RGB_TO_YCbCrJPEG_Cr(Rd, Gd, Bd);

			/* Y' */
			Ypix = Yx + Yy * th->yw;
			Y[Ypix] = RGB_TO_YCbCrJPEG_Y(from[op3 + 2],
						     from[op3 + 1],
						     from[op3 + 0]);
			Y[Ypix + 1] = RGB_TO_YCbCrJPEG_Y(from[op4 + 2],
							 from[op4 + 1],
							 from[op4 + 0]);
			Y[Ypix + th->yw] = RGB_TO_YCbCrJPEG_Y(from[op1 + 2],
							       from[op1 + 1],
							       from[op1 + 0]);
			Y[Ypix + 1 + th->yw] = RGB_TO_YCbCrJPEG_Y(from[op2 + 2],
								   from[op2 + 1],
								   from[op2 + 0]);
			ox += 4 * 2;
		}
		ox = 0;
		oy -= 2 * th->stride;
	}
}

Šoupnul jsem to do cyklu s počtem iterací 30*20 (tj. 20s video při 30 snímcích za sekundu), rozlišení 1280x1024. Konverze je hotová za 6.243s.

Zdálo se mi to dlouho, tak jsem zkoušel optimalizovat a upatlal jsem toto:

void ycbcr_bgr_to_jpeg420_mine(save_thread_t *th, unsigned char *from)
{
  unsigned bytes_per_row;
  unsigned extra_per_row;
  
  unsigned char *row;
  unsigned char *end_pix;
  unsigned char *end_row;
  unsigned char *Y, *Cb, *Cr;
  unsigned char R, G, B;
  
  Y = th->yuv_buffer + (th->yh - 2) * th->yw;
  
	Cb = &th->yuv_buffer[th->yw * th->yh + (th->ch - 1) * th->cw];
	Cr = &th->yuv_buffer[th->yw * th->yh + th->cw * th->ch + (th->ch - 1) * th->cw];
  
  bytes_per_row = th->width * 4;
  extra_per_row = th->stride - bytes_per_row;
  
  end_row = from + th->stride * th->height;
  
  for (; from != end_row; from += th->stride + extra_per_row)
  {
    end_pix = from + bytes_per_row;
    
    for (; from != end_pix; from += 2*4)
    {
      B = (from[0 + 0] + from[4 + 0] + from[th->stride + 0] + from[th->stride + 4 + 0]) >> 2;
      G = (from[0 + 1] + from[4 + 1] + from[th->stride + 1] + from[th->stride + 4 + 1]) >> 2;
      R = (from[0 + 2] + from[4 + 2] + from[th->stride + 2] + from[th->stride + 4 + 2]) >> 2;
      
      *Cb = RGB_TO_YCbCrJPEG_Cb(R, G, B);
      *Cr = RGB_TO_YCbCrJPEG_Cr(R, G, B);
      
      Y[0] = RGB_TO_YCbCrJPEG_Y(from[th->stride + 2], from[th->stride + 1], from[th->stride + 0]);
      Y[1] = RGB_TO_YCbCrJPEG_Y(from[th->stride + 4 + 2], from[th->stride + 4 + 1], from[th->stride + 4 + 0]);
      
      Y[th->yw + 0] = RGB_TO_YCbCrJPEG_Y(from[0 + 2], from[0 + 1], from[0 + 0]);
      Y[th->yw + 1] = RGB_TO_YCbCrJPEG_Y(from[4 + 2], from[4 + 1], from[4 + 0]);
      
      Y += 2;
      ++Cb;
      ++Cr;
    }
    
    Y -= th->yw;
    Y -= th->yw;
    Y -= th->yw;
    
    Cb -= th->cw;
    Cb -= th->cw;
    
    Cr -= th->cw;
    Cr -= th->cw;
  }
}

Prohodil jsem směr (btw řádky vrácené z opengl jsou v rgb bufferu opačně (horní řádek dole)), čtu odshora dolů a zapisuju odspoda nahoru. No nevím jestli to zlepšilo to, nebo jestli to zlepšilo použítí pointerů na některých místech, každopádně dostal jsem se na 4.854s.

Ale protože Jardíkovi se to pořád zdálo hodně a ten obraz se zdál být tmavý, tak jsem našel jiný vzorec na převod a zkoušel jsem to spatlat ... v assembleru.

global bgr_to_jpeg420_sse

; Y = ( (  66 * R + 129 * G +  25 * B + 128) >> 8) +  16
; U = ( ( -38 * R -  74 * G + 112 * B + 128) >> 8) + 128
; V = ( ( 112 * R -  94 * G -  18 * B + 128) >> 8) + 128

section .data align=16

vec_alpha_1_and: align 16
	dd 0x00FFFFFF, 0x00FFFFFF, 0x00FFFFFF, 0x00FFFFFF
vec_alpha_1_or: align 16
	dd 0x01000000, 0x01000000, 0x01000000, 0x01000000

vec_dd_16: align 16
	dd 16, 16, 16, 16

vec_dd_128: align 16
	dd 128, 128, 128, 128
	
vec_y_const: align 16
	dw 25, 129, 66, 128, 25, 129, 66, 128

vec_u_const: align 16
	dw 112, -74, -38, 128, 112, -74, -38, 128

vec_v_const: align 16
	dw -18, -94, 112, 128, -18, -94, 112, 128

section .text

; this expects these xmm registers to be set:
;   xmm4 to vec_alpha_1_and
;   xmm5 to vec_alpha_1_or
;   xmm6 to vec_y_const
;   xmm7 to vec_dd_16
; input vector in xmm0 as (argb3, argb2, argb1, argb0)
; returns result as dwords in xmm1 (y3,y2,y1,y0)
; destroys content of xmm0, xmm2
%macro bgr_to_jpeg420_sse_get_Y 0
	pand      xmm0, xmm4               ; set alpha of each pixel to 1
	por       xmm0, xmm5
	
	movaps    xmm1, xmm0               ; copy to xmm1
	pxor      xmm2, xmm2               ; zero xmm2
	
	punpckhbw xmm0, xmm2         ; interleave xmm2 with first 2 pixels,  xmm = (a3,r3,g3,b3,a2,r2,g2,b2)
	punpcklbw xmm1, xmm2         ; interleave xmm2 with second 2 pixels, xmm = (a1,r1,g1,b1,a0,r0,g0,b0)
	
	pmaddwd xmm0, xmm6           ; xmm = (128a3+66r3 = A3,129g3+25b3 = B3,128a2+66r2 = A2,129g2+25b2 = B2)
	pmaddwd xmm1, xmm6           ; xmm = (128a1+66r1 = A1,129g1+25b1 = B1,128a0+66r0 = A0,129g0+25b0 = B0)
	
	pshufd xmm0, xmm0, 0b11011000 ; (A3, A2, B3, B2)
	pshufd xmm1, xmm1, 0b11011000 ; (A1, A0, B1, B0)
	
	movaps     xmm2, xmm1 ; save xmm1
	punpckhqdq xmm1, xmm0 ; xmm1 = (A3, A2, A1, A0)
	punpcklqdq xmm2, xmm0 ; xmm2 = (B3, B2, B1, B0)
	
	paddd xmm1, xmm2      ; xmm1 = (A3+B3, A2+B2, A1+B1, A0+B0)
	psrld xmm1, 8         ; xmm1 /= 256
	paddd xmm1, xmm7      ; xmm1 += 16
	; xmm1 = ( y3, y2, y1, y0)
%endmacro

; void bgr_to_jpeg420_sse(
;   unsigned char *yuv_buffer,       // rdi, must be 16B aligned
;   const unsigned char *bgr_buffer, // rsi, must be 16B aligned
;   unsigned int width,              // rdx, must be multiple of 4
;   unsigned int height              // rcx, must be multiple of 2
; )
bgr_to_jpeg420_sse:
	push rax
	push rbx
	push rcx
	push rdx
	push rsi
	push rdi
	push rbp
	
	mov rbx, rdx    ; rbx = width, mul will overwrite rdx later
	
	mov rax, rdx
	mul rcx         ; rax = width * height
	
	mov rdx, rbx    ; move back width to rdx
	
	mov rbp, rdi
	add rbp, rax    ; rbp = pointer to U part
	
	mov rbx, rax
	shr rbx, 2      ; rbp+rbx = pointer to V part
	
	shl rax, 2      ; rax = rgb buffer size (width * height * 4)
	
	mov  r9, rdx    ; r9 = width
	shl  r9, 2      ; r9 = width*4 = rgb_buffer stride
	
	add rsi, rax    ; make rsi point to the last but one row
	sub rsi, r9     ; of rgb_buffer
	sub rsi, r9
	
	shr rcx, 1      ; we process 2 rows in one loop,so make height half
	
	movaps xmm4, [rel vec_alpha_1_and]
	movaps xmm5, [rel vec_alpha_1_or]
	movaps xmm6, [rel vec_y_const]
	movaps xmm7, [rel vec_dd_16]
	movaps xmm3, [rel vec_dd_128]
	
	movaps xmm10, [rel vec_u_const]
	movaps xmm11, [rel vec_v_const]
	
.__outer_loop:
	
	mov r8, rdx           ; r8 = width
	shr r8, 2             ; r8 = width / 4 (4 pixels in one loop)
	
.__inner_loop:
	movaps    xmm8, [rsi]              ; load the 4 pixels in first row, xmm = (argb3,argb2,argb1,argb0)
	movaps    xmm9, [rsi+r9]           ; load the 4 pixels in second row, xmm = (argb3,argb2,argb1,argb0)
	
	movaps xmm0, xmm9             ; get Y value for each pixel in the vector
	bgr_to_jpeg420_sse_get_Y
	; xmm1 = ( y3, y2, y1, y0)
	
	pextrw eax, xmm1, 0           ; store it in yuv_buffer
	mov BYTE [rdi], al
	pextrw eax, xmm1, 2
	mov BYTE [rdi+1], al
	pextrw eax, xmm1, 4
	mov BYTE [rdi+2], al
	pextrw eax, xmm1, 6
	mov BYTE [rdi+3], al
	
	movaps xmm0, xmm8             ; get Y value for each pixel in the vector
	bgr_to_jpeg420_sse_get_Y
	; xmm1 = ( y3, y2, y1, y0)
	
	pextrw eax, xmm1, 0           ; store it in yuv_buffer
	mov BYTE [rdi+rdx], al
	pextrw eax, xmm1, 2
	mov BYTE [rdi+rdx+1], al
	pextrw eax, xmm1, 4
	mov BYTE [rdi+rdx+2], al
	pextrw eax, xmm1, 6
	mov BYTE [rdi+rdx+3], al
	
	
	; Calc U and V here
	movaps    xmm0, xmm8
	pavgb     xmm0, xmm9
	pshufd    xmm1, xmm0, 0b10000000
	pshufd    xmm0, xmm0, 0b11010000
	pavgb     xmm0, xmm1
	
	pand      xmm0, xmm4               ; set alpha of each pixel to 1
	por       xmm0, xmm5
	
	pxor      xmm2, xmm2
	punpckhbw xmm0, xmm2               ; xmm0 now contains words (a1,r1,g1,b1,a0,r0,g0,b0)
	
	movaps    xmm8, xmm0               ; save it to xmm8 as U calculation will overwrite xmm0
	
	; calc U
	pmaddwd xmm0, xmm10                ; xmm = (128a1-38r1 = A1,-74g1+112b1 = B1,128a0-38r0 = A0,-74g0+112b0 = B0)
	pshufd  xmm1, xmm0, 0b10110001
	paddd xmm0, xmm1                   ; xmm = (A1+B1, A1+B1, A0+B0, A0+B0)
	psrld xmm0, 8                      ; xmm /= 256
	paddd xmm0, xmm3                   ; xmm += 128
	
	pextrw eax, xmm0, 0
	mov BYTE [rbp], al
	pextrw eax, xmm0, 4
	mov BYTE [rbp+1], al
	
	; calc V
	movaps xmm0, xmm8
	pmaddwd xmm0, xmm11                ; xmm = (128a1-38r1 = A1,-74g1+112b1 = B1,128a0-38r0 = A0,-74g0+112b0 = B0)
	pshufd  xmm1, xmm0, 0b10110001
	paddd xmm0, xmm1                   ; xmm = (A1+B1, A1+B1, A0+B0, A0+B0)
	psrld xmm0, 8                      ; xmm /= 256
	paddd xmm0, xmm3                   ; xmm += 128
	
	pextrw eax, xmm0, 0
	mov BYTE [rbp+rbx], al
	pextrw eax, xmm0, 4
	mov BYTE [rbp+rbx+1], al
	
	add rbp, 2
	
	add rsi, 16          ; rgb_buffer += 16
	add rdi, 4           ; Y += 4
	dec r8               ; remaining iterations for inner loop
	jnz .__inner_loop
	
	add rdi, rdx        ; skip the row
	
	sub rsi, r9         ; make rsi point to the 2 previous rows
	sub rsi, r9
	sub rsi, r9
	
	dec rcx              ; decrement number of remaining rows
	jnz .__outer_loop    ; if not zero, continue
	
.__end:
	pop rbp
	pop rdi
	pop rsi
	pop rdx
	pop rcx
	pop rbx
	pop rax
	ret

S tímto jsem se dostal na 3.509s. Kód je docela ošklivý a patlal jsem to s instrukcema, co jsem našel v různých seznamech instrukcí po netu se válících. No a proč to sem házím ... kdyby se našel nějaký expert co by to rád zoptimalizoval, nebo poradil, kde a co tam zlepšit, .. byl bych moc vděčný.

Používám tam x86-64 registry, s xmm registry jsem se prostě nevešel a když ty konstanty nebudu mít přednačtené a budu je pak z paměti načítat v cyklu, tak to zdržuje a algoritmus to zpomalí. A když už jsem se nevešel s xmm registry, tak už jsem rovnou použil i x86-64 GPR.

       

Komentáře

[XSSds]

Mozem poprosit x86-16 a x86-128 bitovu verziu?

[Jardík]

Můžete

[User682]

zdravim,

dam Vam jenom inspiraci na hledani, ale ne reseni.

Kdysi jsem prochazel git repositar ffmpeg projektu a je tam vcelku dost asembleru na ruzne konverze. Mozna, ze tam je neco co hledate.

gf

[Miriam]

Ahoj, ještě pořád nemáš holku? Moje nabídka stále platí, nezajdem někam na drink?

[Josef Kufner]

Řekl bych, že ta první otázka je vzhledem k obsahu zápisku celkem zbytečná

[Miriam]

Jsem slušná holka, né děvka. K tomu patří i formality ;)

[Miriam]

A vlastně bych už měla jít spát. Tak pa a brou noc

[blizzz]

Áno, určite si holka ale s kokotom :D

[XSSds]

Ba dum tsss ... LOL True Story :D

[blizzz]

čo tým chcel básnik povedať?

[munich machine]

Ha - ha - ha.

[Josef Kufner]

Však to samozřejmě. Jsem tím ani nic nevhodného nemyslel

[Jardík]

Promiň, ale já nemůžu, mám poruchu osobnosti. Leda bys chtěla pejska, kterého bys ráda krmila a zbytek dne bys si ho nevšímala aby si mohl vklidu omílat, že je život zbytečnej a že nemá smysl chodit do práce, protože za chození do práce musíš platit nadlidem u koryta.

[Vskutečnosti Saýc]

Premyslim jak rychle by o zvladla CUDA.

[Miloslav Ponkrác]

Moc tomu nerozumím a proto se předem omlouvám za blbosti.

Ale pokud vím, v tomto systém je U a V barevná složka, která se odvozuje od jasové složky.

Tedy po výpočtu jasové složky stačí:

U = konstanta * (B - Y')

V = konstanta * (R - Y')

Přijde mi zbytečné rozepisovat to násobením přes plnou transformační matici.

[Miloslav Ponkrác]

U = 0,4921 * (B - Y')

V = 0,8773 * (R - Y')

[Josef Kufner]

Jo, jo... a na Wikipedii je těch (podobně jednoduchých) vzorečků spousta.

Btw, pár desítek let staré televize tenhle převod počítají analogově a stíhají to v reálném čase

[Miloslav Ponkrác]

Analogové operace jsou mnohonásobně rychlejší a efektivnější než digitální operace. V tomto případě je to i primitivní elektronický obvod.

Ale hw obvod navržený na míru by to spočítal i digitálně v reálném čase.

Ostatně to, že analogové obvody mají vyšší rychlosti i kapacitu (ale nižší odolnost proti šumu a ztrátě informací) už je znovuobjeven. Takové SSD disky a flashky už v zásadě nezaznamenávají informaci digitálně (pokud nejsou SLC) ale v zásadě analogově (MLC). Navyšuje to kapacitu, ale zhoršuje stabilitu uložených dat a rychleji to ztrácí data.

Kromě šumu je analogová technika ve všech vlastnostech o mnoho řádů napřed. To co se dá v analogové technice vyřešit součástkami za 10 Kč, na vyřešení stejné operace v digitální technice je potřeba stroj za sto tisíc dolarů a ještě bude pokulhávat za analogem.

Analogová elektronika je dokonce schopná řešit diferenciální rovnice v reálném čase.

[Martin Mareš]

Kromě šumu je analogová technika ve všech vlastnostech o mnoho řádů napřed. To co se dá v analogové technice vyřešit součástkami za 10 Kč, na vyřešení stejné operace v digitální technice je potřeba stroj za sto tisíc dolarů a ještě bude pokulhávat za analogem.

Pokud si vyberete jednu z těch několika málo operací, které se dají analogově snadno realizovat, tak ano

Ale rád se nechám přesvědčit – stačí, když ukážete, jak jednoduše analogově spočítat RSA

[pc2005]

No já nevím, není náhodou MLC pouze soustava o jiném základu než 2?

BTW ono záleží na rychlosti. S operačníma zesilovačema bys měl docela problém pro rychlost operací frekvencí xGHz. Resp ty operáky by asi nestály 10Kč.

[ustp]

Trochu demagogie s tema ssd. MLC bunky pouzivaj vic napetovejch urovni, ale porad ukladaj digitalni data.

[Ondrej 'SanTiago' Zajicek]

Presne tak, 'digitalni' obvykle znaci, ze vyznamy jsou prirazene diskretnim urovnim, nikoliv ze tech urovni musi byt pouze 2. Pouziti vice diskretnich urovni neni neobvykle (napr. 100Mbps a 1Gbps ethernet).

[pc2005]

Vlastně i výstup z DAC by byl před low pass filtrem diskrétní.

[Ondrej 'SanTiago' Zajicek]

To vicemene ano, ale nemel diskretne definovane urovne.

[pc2005]

No leda přes dílky z nějakého Vref.

[kralyk z abclinuxu]

Hmm, tohle bylo lepší než obvyklé statě o C++

[Agent]

Možná trochu OT, ale připomělo mi to, jak v Amáru vycházel návod na videotitulkovač někdy v polovině 90 let. Bedna plná integráčů a jediný co to umělo byly titulky. Netrvalo dlouho a takováto udělátka ztratila smysl. Jediná výhoda možná byla, že to bylo levnější než tehdejší PC s MPEG2 kartou.

[Vskutečnosti Saýc]

No staci si vsimnout, ze jardik zpracuje tricetisekundove video za necele ctyri sekundy, coz na realny cas spokojene staci.

[Miloslav Ponkrác]

Nikdo se ani slovem nevyjádřil, že jardík nestíhá v reálném čase.

Možná by si stačilo všimnout, že řeč byla o tom, že jakýsi předpotopní televizáč to stíhá v reálném čase.

A možná by si stačilo všimnout, že nebyla jakákoli zmínka o jardíkovi, jeho algoritmu ani srovnávání rychlostí.

Jinak řečeno, chápání psaného textu, nic moc. Doučte se to, studente a trénujte.

[Vskutečnosti Saýc]

Mozna by si stacilo vsimnout, na ktery prispevek odpovidam. Napoveda: neni to ten Vas.

[Jardík]

Mým výsledkem má být yuv420, pro každé 4 pixely (2x2) je tam jeden U a jeden V. Takže když mám pro 4 pixely už 4 Y', tak abych mohl použít ten vzorec řekněme pro U, musel bych nejdříve zprůměrovat B hodnoty pro ty 4 pixely, zprůměrovat Y hodnoty pro ty pixely a použít vzorec.

Momentálně to dělám tak, pro každý pixel vypočítám Y, zprůměruju rgb hodnoty pro 2x2 sousední pixely a spočtu U. Pravděpodobně by mi použítí "jednodušší" vzorce v tomto případě moc nepomohlo. Leda bych U nepočítal z průměrů ale pouze z jednoho pixelu, ale pak bych asi nedostal hezký výsledek.

[Miloslav Ponkrác]

Jak jsem psal, nerozumím tomu do detailu.

Nemusíte průměrovat nic, prostě hodnoty sečtete a to bude 4× požadovaná hodnota. Konstanty podle toho upravíte.

Rovnou se přiznám, že kódy v článku jsme moc neluštil. Chybí mi další informace jako definice struktury a vůbec nějaké komentáře, takže jsem to jenom přelétl.

Něčemu málo rozumím, ale ne celku. Nerozumím přesně v jakém formátu a rozsahu chodí vstupní RGB data.

Ale jednu věc vím zcela určitě. Volba lepšího algoritmu a pohrání si s matematikou dává mnohem větší zrychlení, než ďábelsky kódovat nanosekundy pomocí asm. Možná se ukáže, že to není k ničemu.

Určitě bych si vzal základní „atom“, tedy nejmenší pravidelný rámec, který mám vyrobit. Nejdříve bych si pohrál matematicky, pak by se snažil převést do bajtových operací tak, aby se to šikovně optimalizovalo. Neoptimalizoval bych jeden bod, ale celý nejmenší blok, který se pravidelně opakuje.

Mám trochu pocit z článku, že náhodně hledáte vzorce a že v matematice plavete a tak se snažíte ďábelsky kódovat ty nanosekundy.

Pochopit ty vzorce je jednoduché. Transformační matice s převodem na int počítání a posun je dělána právě pro SIMD, kde jsou tyhle nadbytečné operace vcelku zdarma. I proto Vám v klasickém C vyšel mnohem pomalejší algoritmus, protože jste tam rval věci optimalizované pro SIMD operace.

[Jardík]

Tenhle obrázek z wiki je můj výstupní formát.

Vstup jsou za sebou jdoucí pixely B8G8R8A8 s tím, že celý obraz je vzhůru nohama (tak ho vrátí glReadPixels).

[Miloslav Ponkrác]

Snad rozumím. Pro jistotu popíšu:

V zásadě se tedy kódují po čtveřicích pixelů. U každé čtveřice pixelů (čtverec 2×2 body) se spočítá Y', tedy vyjdou čtyři bajty Y'. Dále se ale pro celou čtveřici pixelů spočítá jen jedno U a následně jedno V.

Celá čtveřice bodů, které přijdou jako 16 bajtů RGBA se zakódují do šesti bajtů Y'1, Y'2, Y'3, Y'4, U, V.

Analogová hodnota Y' se pohybuje v rozsahu od 0 do 1, a transformuje se do bajtu Y' jako hodnota od 16 do 235. Hodnoty pod 16 a nad 235 se nepoužívají a těchto hodnot nemohou uložené bajty Y' nikdy nabýt.

[Miloslav Ponkrác]

Co takhle?

Spočítat pro každou čtveřici pixelů Y':

Y' = (66R + 129G + 25B) >> 8 + 16;

Ale algoritmus bych přepsal takto:

===============================================

Varianta 1 pro extrémní šetření pamětí:

a) Udělal bych si 128 KB velkou předpočítanou tabulku, tedy pole o 64 K prvcích, každý prvek by obsahovat 2 bajtovou hodnotu.

index_1 = B + 256 * G; pole_1[index_1] = (16 * 256 + 25 * B + 129 * G);

b) Udělal bych si druhou 512 B velkou předpočítanou tabulku:

index_2 = R; pole_2[index_2] = (66 * R);

c) Výsledné Y' bych spočítal jako:

Y' = (pole_1[index_1] + pole_2[index_2]) >> 8;

===============================================

Nejrychlejší varianta 2 pro zabrání 16 MB RAM (dnes mají počítače několik GB RAM):

a) předpočítanou tabulku pro převod:

index = B + 256 * G + 65536 * R (tedy nic počítat nemusíte, pouze použijete bajty B,G,R tak jak na intel procesorech leží) pole[index] = výsledná hodnota Y' bajtu

[kutr]

Jak jsem psal, nerozumím tomu do detailu. Ale jednu věc vím zcela určitě. Volba lepšího algoritmu a pohrání si s matematikou dává mnohem větší zrychlení, než ďábelsky kódovat nanosekundy pomocí asm. Možná se ukáže, že to není k ničemu.

Jenže Jardík obvykle jde cestou předčasné optimalizace a ukazuje nám proč se nemá dělat. (viz. jeho starší blog kde ďábelsky optimalizoval v ASM, až se dostal na stejný čas jako naivní C program)

[Miloslav Ponkrác]

Mě naučili, co to je skutečná optimalizace, a jak efektivně matematické záležitosti a algoritmy mohou fungovat až na vysoké škole.

Do té doby jsem také optimalizoval asm.

Jenže díky bohu za to. Protože jsem se naučil optimalizovat algoritmy pro procesor a programovací jazyky a pak později i na úrovni algoritmů.

To je mimochodem z jeden důvodů, proč je dobré udělat vysokou školu. Naučíte se spoustu věcí, na které nikdy nepřijdete (protože Vás k tomu nikdo nebude nutit a protože to nikdo kolem vás nebude znát). Uděláte skok o desítky let dopředu oproti stavu kdy se budete učit pouze praxí.

Jardík si možná prošlápne cestu jak se to dělat nemá, ale jako vedlejší efekt se naučí být machr v postalgoritmické optimalizaci. Navíc ho to baví.

Časem na to přijde a bude machr po všech stránkách.

Nemůžete být v něčem dobří, pokud nejdete občas vedlejšími cestami, neděláte chyby a nevyzkoušíte si i co nefunguje.

Ale Jardíkův zápal se mi strašně líbí.

[anon]

Fajnovy posledni odstavce.

Kazdopadne s tou vs nemohu souhlasit - me vs co se tyka odbornych znalosti dala temer nic. Nadruhou stranu jsem mel hromadu volnyho casu ve kterem jsem mohl sam nabrat hodne cenych zkusenosti h oboru I mimo. I kdyz vlastne kurzy mimo it mi daly zabrat no it me akorat trapilo ze to vsecho jsem se naucil sam uz pred x rokama a kdy jako budem delat neco zajimavyho.

[Kvakor]

Kdyby to bylo YUV na RGB, tam to jsem před lety dělal a někde jsem našel a přepsal rutinu, která mohutně používá přetpočítané tabulky a ve vnitřením cyklu je jenom sčítaní, bitové posuvy a vytažení hodnoty z tabulky. Rozdíl mezi verzí co v cyklu počítala podle vzorečku a tou s předpočítanými tabulkami byl tenkrát asi o řád. Na druhou stranu, dnešní procesory jsou vůči pamětem tak rychlé, že násobení může vyjít levněji než načtení hodnoty z paměti. I tak budu komentovat jen céčkovou část:

Aby se vyhodilo násobení z hlavniho cyklu, stačilo by předpočítat si dopředu tabulky pro každou hodnotu Y, Cr a Cb místo násobení uvnitř maker RGB_TO_YCbCrJPEG_* a pošoupnout všechno doprava až v okamžiku, kdy se ukládá finální hodnota (popř. to procpat ještě přes nějakou tabulku, co by dělala clipping). Dá se tím i zaokrouhlovat, pokud se dobře nastaví konstanta, co se přičítá před posuvem.

Nebo možná udělat "unrolling" toho chroma subsamplingu a vecpal do něj i počítání luminance, protože mám pocit, že by to méně zaplácávalo cache - pokud už je některá hodnota z paměti načtená, měla by být využitá až "do mrtě". A to pošoupnutí ">> 2" u onoho supsamplingu se mi také moc nelíbí, protože to se může bez ztráty kytičky udělat až na konci a získámé tím trochu vyšší přenost.

No a ještě bych nejspíš držel dva ukazatele do RGB pole from, pošouplé o th->stride a jenom je inkrementoval, protože i když z toho překladač udělá nepřímé adresování, stále to bude sahat do paměti. To samé i u výsledného Y a th->yw.

[Miloslav Ponkrác]

Já bych reálně počítal Y' z RGB.

Pak bych udělal předpočítané tabulky (B,Y')->U a (R,Y')->V.

[Bedňa]

Jardíku držím palce.

[mirec]

Držím palec. Práve riešim real time renderovanie OpenGL do YUV streamu s následným enkódovaním do MPEG2 a vysielaním cez rtmp.

[Miloslav Ponkrác]

Tak jsem prostudoval ITU dokument.

Pochopil jsem, proč se výpočet počítá transformační maticí a nikoli odčítáním. Protože výsledné Y' je jen přibližným zaokrouhlením a celý výpočet se počítá na integerech jako fixed point čísla.

Použití bajtu Y' by přidalo do výpočtu zaokrouhlovací chybu (výsledek by se mohl lišit o jedničku nahoru nebo dolů od skutečné hodnoty). Výpočet by se tedy musel prováděl na 16bitovém integeru namísto 8bitovém bajtu.

Tedy výpočet Y' by pak do výpočtu U a V vstoupit jako celé fixed point číslo (16bitový integer) – tedy Y' před rotací doprava o 8 bitů.

[Duff]

Pro výpočet pomocí procesoru by se to dalo pěkně paralelizovat. Matice by se rozdělila na počet částí jako procesorů/jader.

A ještě lepší možnost je to počítat pomocí grafiky a pixel shaderů. I na smartfounech to půjde v reálném čase 30 snímků za sekundu.