1 //fibe.h	FAST IIR BLUR ENGINE
2 
3 /*
4 Copyright (C) 2011  Marko Cebokli    http://lea.hamradio.si/~s57uuu
5 
6 
7  This program is free software; you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation; either version 2 of the License, or
10  (at your option) any later version.
11 
12  This program is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with this program; if not, write to the Free Software
19  Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
20 
21 */
22 
23 /*
24    NOTE: float_rgba must be declared externally as
25 
26  typedef struct
27     {
28     float r;
29     float g;
30     float b;
31     float a;
32     } float_rgba;
33 
34     --------------------------------
35        CONTENTS OF FIBE.H FILE:
36     --------------------------------
37 
38 calcab_lp1()	auxilliary function to calculate lowpass
39         tap coefficients for FIBE-2
40 
41 young_vliet()	auxilliary function to calculate tap coefs
42         for Gauss approximation with FIBE-3
43 
44 rep()		auxilliary function to calculate wraparound
45         values for FIBE-2
46 
47 fibe1o_8()	one tap quadrilateral IIR filter
48         speed optimized C function
49         includes 8bit/float conversions
50 
51 fibe2o_8()	two tap quadrilateral IIR filter
52         speed optimized C function
53         includes 8bit/float conversions
54 
55 fibe3_8()	three tap quadrilateral IIR filter
56         includes 8bit/float conversions
57 
58 The functions work internally with floats. I have included
59 the 8bit/float conversion into the first and last
60 processing loops, to avoid two additional cache polluting
61 and therefore time consuming "walks" through memory.
62 
63 */
64 
65 
66 
67 //edge compensation average size
68 #define EDGEAVG 8
69 
70 #include <sys/types.h>
71 #include <string.h>
72 #include "frei0r_math.h"
73 
74 //---------------------------------------------------------
75 //koeficienti za biquad lowpass  iz f in q
76 // f v Nyquistih    0.0 < f < 0.5
calcab_lp1(float f,float q,float * a0,float * a1,float * a2,float * b0,float * b1,float * b2)77 void calcab_lp1(float f, float q, float *a0, float *a1, float *a2, float *b0, float *b1, float *b2)
78 {
79     float a,b;
80 
81     a=sinf(PI*f)/2.0/q;
82     b=cosf(PI*f);
83     *b0=(1.0-b)/2.0;
84     *b1=1.0-b;
85     *b2=(1.0-b)/2.0;
86     *a0=1.0+a;
87     *a1=-2.0*b;
88     *a2=1.0-a;
89 }
90 
91 //---------------------------------------------------------
92 //3tap iir coefficients for Gauss approximation according to:
93 //Ian T. Young, Lucas J. van Vliet:
94 //Recursive implementation of the Gaussian filter
95 //Signal Processing 44 (1995) 139-151
96 // s=sigma    0.5 < s < 200.0
young_vliet(float s,float * a0,float * a1,float * a2,float * a3)97 void young_vliet(float s, float *a0, float *a1, float *a2, float *a3)
98 {
99     float q;
100 
101     q=0.0;
102     if (s>2.5)
103     {
104         q = 0.98711*s - 0.96330;
105     }
106     else
107     {	//to velja za s>0.5 !!!!
108         q = 3.97156 - 4.14554*sqrtf(1.0-0.26891*s);
109     }
110 
111     *a0 = 1.57825 + 2.44413*q + 1.4281*q*q + 0.422205*q*q*q;
112     *a1 = 2.44413*q + 2.85619*q*q + 1.26661*q*q*q;
113     *a2 = -1.4281*q*q - 1.26661*q*q*q;
114     *a3 = 0.422205*q*q*q;
115 }
116 
117 //---------------------------------------------------
118 //kompenzacija na desni
119 //c=0.0 "odziv na zacetno stanje" (zunaj crno)
120 //gain ni kompenziran
rep(float v1,float v2,float c,float * i1,float * i2,int n,float a1,float a2)121 void rep(float v1, float v2, float c, float *i1, float *i2, int n,  float a1, float a2)
122 {
123     int i;
124     float lb[8192];
125 
126     lb[0]=v1;lb[1]=v2;
127     for (i=2;i<n-2;i++)
128     {
129         lb[i]=c-a1*lb[i-1]-a2*lb[i-2];
130     }
131 
132     lb[n-2]=0.0;lb[n-1]=0.0;
133     for (i=n-3;i>=0;i--)
134     {
135         lb[i]=lb[i]-a1*lb[i+1]-a2*lb[i+2];
136     }
137 
138     *i1=lb[0]; *i2=lb[1];
139 }
140 
141 //---------------------------------------------------------
142 // 1-tap IIR v 4 smereh
143 //optimized for speed
144 //loops rearanged for more locality (better cache hit ratio)
145 //outer (vertical) loop 2x unroll to break dependency chain
146 //simplified indexes
fibe1o_8(const uint32_t * inframe,uint32_t * outframe,float_rgba * s,int w,int h,float a,int ec)147 void fibe1o_8(const uint32_t* inframe, uint32_t* outframe, float_rgba *s, int w, int h, float a, int ec)
148 {
149     int i,j;
150     float b,g,g4,avg,avg1,cr,cg,cb,g4a,g4b;
151     int p,pw,pj,pwj,pww,pmw;
152 
153     avg=EDGEAVG;	//koliko vzorcev za povprecje pri edge comp
154     avg1=1.0/avg;
155 
156     g=1.0/(1.0-a);
157     g4=1.0/g/g/g/g;
158 
159     //predpostavimo, da je "zunaj" crnina (nicle)
160     b=1.0/(1.0-a)/(1.0+a);
161 
162     //prvih avg vrstic
163     for (i=0;i<avg;i++)
164     {
165         p=i*w;pw=p+w;
166         if (ec!=0)
167         {
168             cr=0.0;cg=0.0;cb=0.0;
169             for (j=0;j<avg;j++)
170             {
171                 s[p+j].r=(float)(inframe[p+j]&0xFF);
172                 s[p+j].g=(float)((inframe[p+j]&0xFF00)>>8);
173                 s[p+j].b=(float)((inframe[p+j]&0xFF0000)>>16);
174                 cr=cr+s[p+j].r;
175                 cg=cg+s[p+j].g;
176                 cb=cb+s[p+j].b;
177             }
178             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
179             s[p].r=cr*g+b*(s[p].r-cr);
180             s[p].g=cg*g+b*(s[p].g-cg);
181             s[p].b=cb*g+b*(s[p].b-cb);
182         }
183         else
184             for (j=0;j<avg;j++)
185             {
186                 s[p+j].r=(float)(inframe[p+j]&0xFF);
187                 s[p+j].g=(float)((inframe[p+j]&0xFF00)>>8);
188                 s[p+j].b=(float)((inframe[p+j]&0xFF0000)>>16);
189             }
190 
191         for (j=1;j<avg;j++)	//tja  (ze pretvorjeni)
192         {
193             s[p+j].r=s[p+j].r+a*s[p+j-1].r;
194             s[p+j].g=s[p+j].g+a*s[p+j-1].g;
195             s[p+j].b=s[p+j].b+a*s[p+j-1].b;
196         }
197         for (j=avg;j<w;j++)	//tja  (s pretvorbo)
198         {
199             s[p+j].r=(float)(inframe[p+j]&0xFF);
200             s[p+j].g=(float)((inframe[p+j]&0xFF00)>>8);
201             s[p+j].b=(float)((inframe[p+j]&0xFF0000)>>16);
202             s[p+j].r=s[p+j].r+a*s[p+j-1].r;
203             s[p+j].g=s[p+j].g+a*s[p+j-1].g;
204             s[p+j].b=s[p+j].b+a*s[p+j-1].b;
205         }
206 
207         if (ec!=0)
208         {
209             cr=0.0;cg=0.0;cb=0.0;
210             for (j=w-avg;j<w;j++)
211             {
212                 cr=cr+s[p+j].r;
213                 cg=cg+s[p+j].g;
214                 cb=cb+s[p+j].b;
215             }
216             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
217             s[pw-1].r=cr*g+b*(s[pw-1].r-cr);
218             s[pw-1].g=cg*g+b*(s[pw-1].g-cg);
219             s[pw-1].b=cb*g+b*(s[pw-1].b-cb);
220         }
221         else
222         {
223             s[pw-1].r=b*s[pw-1].r;
224             s[pw-1].g=b*s[pw-1].g;
225             s[pw-1].b=b*s[pw-1].b;
226         }
227         for (j=w-2;j>=0;j--)	//nazaj
228         {
229             s[p+j].r=a*s[p+j+1].r+s[p+j].r;
230             s[p+j].g=a*s[p+j+1].g+s[p+j].g;
231             s[p+j].b=a*s[p+j+1].b+s[p+j].b;
232         }
233     }
234 
235     //prvih avg vrstic samo navzdol (nazaj so ze)
236     for (i=0;i<w;i++)
237     {
238         if (ec!=0)
239         {
240             cr=0.0;cg=0.0;cb=0.0;
241             for (j=0;j<avg;j++)
242             {
243                 cr=cr+s[i+w*j].r;
244                 cg=cg+s[i+w*j].g;
245                 cb=cb+s[i+w*j].b;
246             }
247             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
248             s[i].r=cr*g+b*(s[i].r-cr);
249             s[i].g=cg*g+b*(s[i].g-cg);
250             s[i].b=cb*g+b*(s[i].b-cb);
251         }
252         for (j=1;j<avg;j++)	//dol
253         {
254             s[i+j*w].r=s[i+j*w].r+a*s[i+w*(j-1)].r;
255             s[i+j*w].g=s[i+j*w].g+a*s[i+w*(j-1)].g;
256             s[i+j*w].b=s[i+j*w].b+a*s[i+w*(j-1)].b;
257         }
258     }
259 
260     for (i=avg;i<h-1;i=i+2)	//po vrsticah navzdol
261     {
262         p=i*w; pw=p+w; pww=pw+w; pmw=p-w;
263         if (ec!=0)
264         {
265             cr=0.0;cg=0.0;cb=0.0;
266             for (j=0;j<avg;j++)
267             {
268                 s[p+j].r=(float)(inframe[p+j]&0xFF);
269                 s[p+j].g=(float)((inframe[p+j]&0xFF00)>>8);
270                 s[p+j].b=(float)((inframe[p+j]&0xFF0000)>>16);
271                 cr=cr+s[p+j].r;
272                 cg=cg+s[p+j].g;
273                 cb=cb+s[p+j].b;
274             }
275             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
276             s[p].r=cr*g+b*(s[p].r-cr);
277             s[p].g=cg*g+b*(s[p].g-cg);
278             s[p].b=cb*g+b*(s[p].b-cb);
279             cr=0.0;cg=0.0;cb=0.0;
280             for (j=0;j<avg;j++)
281             {
282                 s[pw+j].r=(float)(inframe[pw+j]&0xFF);
283                 s[pw+j].g=(float)((inframe[pw+j]&0xFF00)>>8);
284                 s[pw+j].b=(float)((inframe[pw+j]&0xFF0000)>>16);
285                 cr=cr+s[pw+j].r;
286                 cg=cg+s[pw+j].g;
287                 cb=cb+s[pw+j].b;
288             }
289             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
290             s[pw].r=cr*g+b*(s[pw].r-cr);
291             s[pw].g=cg*g+b*(s[pw].g-cg);
292             s[pw].b=cb*g+b*(s[pw].b-cb);
293         }
294         else
295         {
296             for (j=0;j<avg;j++)
297             {
298                 s[p+j].r=(float)(inframe[p+j]&0xFF);
299                 s[p+j].g=(float)((inframe[p+j]&0xFF00)>>8);
300                 s[p+j].b=(float)((inframe[p+j]&0xFF0000)>>16);
301             }
302             for (j=0;j<avg;j++)
303             {
304                 s[pw+j].r=(float)(inframe[pw+j]&0xFF);
305                 s[pw+j].g=(float)((inframe[pw+j]&0xFF00)>>8);
306                 s[pw+j].b=(float)((inframe[pw+j]&0xFF0000)>>16);
307             }
308         }
309         for (j=1;j<avg;j++)	//tja  (ze pretvojeni)
310         {
311             pj=p+j;pwj=pw+j;
312             s[pj].r=s[pj].r+a*s[pj-1].r;
313             s[pj].g=s[pj].g+a*s[pj-1].g;
314             s[pj].b=s[pj].b+a*s[pj-1].b;
315             s[pwj].r=s[pwj].r+a*s[pwj-1].r;
316             s[pwj].g=s[pwj].g+a*s[pwj-1].g;
317             s[pwj].b=s[pwj].b+a*s[pwj-1].b;
318         }
319         for (j=avg;j<w;j++)	//tja  (s pretvorbo)
320         {
321             pj=p+j;pwj=pw+j;
322             s[pj].r=(float)(inframe[pj]&0xFF);
323             s[pj].g=(float)((inframe[pj]&0xFF00)>>8);
324             s[pj].b=(float)((inframe[pj]&0xFF0000)>>16);
325             s[pj].r=s[pj].r+a*s[pj-1].r;
326             s[pj].g=s[pj].g+a*s[pj-1].g;
327             s[pj].b=s[pj].b+a*s[pj-1].b;
328             s[pwj].r=(float)(inframe[pwj]&0xFF);
329             s[pwj].g=(float)((inframe[pwj]&0xFF00)>>8);
330             s[pwj].b=(float)((inframe[pwj]&0xFF0000)>>16);
331             s[pwj].r=s[pwj].r+a*s[pwj-1].r;
332             s[pwj].g=s[pwj].g+a*s[pwj-1].g;
333             s[pwj].b=s[pwj].b+a*s[pwj-1].b;
334         }
335 
336         if (ec!=0)
337         {
338             cr=0.0;cg=0.0;cb=0.0;
339             for (j=w-avg;j<w;j++)
340             {
341                 cr=cr+s[p+j].r;
342                 cg=cg+s[p+j].g;
343                 cb=cb+s[p+j].b;
344             }
345             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
346             s[pw-1].r=cr*g+b*(s[pw-1].r-cr);
347             s[pw-1].g=cg*g+b*(s[pw-1].g-cg);
348             s[pw-1].b=cb*g+b*(s[pw-1].b-cb);
349             cr=0.0;cg=0.0;cb=0.0;
350             for (j=w-avg;j<w;j++)
351             {
352                 cr=cr+s[pw+j].r;
353                 cg=cg+s[pw+j].g;
354                 cb=cb+s[pw+j].b;
355             }
356             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
357             s[pww-1].r=cr*g+b*(s[pww-1].r-cr);
358             s[pww-1].g=cg*g+b*(s[pww-1].g-cg);
359             s[pww-1].b=cb*g+b*(s[pww-1].b-cb);
360         }
361         else
362         {
363             s[pw-1].r=b*s[pw-1].r;	//rep H
364             s[pw-1].g=b*s[pw-1].g;
365             s[pw-1].b=b*s[pw-1].b;
366             s[pww-1].r=b*s[pww-1].r;
367             s[pww-1].g=b*s[pww-1].g;
368             s[pww-1].b=b*s[pww-1].b;
369         }
370 
371         //zacetek na desni
372         s[pw-2].r=s[pw-2].r+a*s[pw-1].r;	//nazaj
373         s[pw-2].g=s[pw-2].g+a*s[pw-1].g;
374         s[pw-2].b=s[pw-2].b+a*s[pw-1].b;
375         s[pw-1].r=s[pw-1].r+a*s[p-1].r;		//dol
376         s[pw-1].g=s[pw-1].g+a*s[p-1].g;
377         s[pw-1].b=s[pw-1].b+a*s[p-1].b;
378 
379         for (j=w-2;j>=1;j--)	//nazaj
380         {
381             pj=p+j;pwj=pw+j;
382             s[pj-1].r=a*s[pj].r+s[pj-1].r;
383             s[pj-1].g=a*s[pj].g+s[pj-1].g;
384             s[pj-1].b=a*s[pj].b+s[pj-1].b;
385             s[pwj].r=a*s[pwj+1].r+s[pwj].r;
386             s[pwj].g=a*s[pwj+1].g+s[pwj].g;
387             s[pwj].b=a*s[pwj+1].b+s[pwj].b;
388             //zdaj naredi se en piksel vertikalno dol, za vse stolpce
389             //dva nazaj, da ne vpliva na H nazaj
390             s[pj].r=s[pj].r+a*s[pmw+j].r;
391             s[pj].g=s[pj].g+a*s[pmw+j].g;
392             s[pj].b=s[pj].b+a*s[pmw+j].b;
393             s[pwj+1].r=s[pwj+1].r+a*s[pj+1].r;
394             s[pwj+1].g=s[pwj+1].g+a*s[pj+1].g;
395             s[pwj+1].b=s[pwj+1].b+a*s[pj+1].b;
396         }
397         //konec levo
398         s[pw].r=s[pw].r+a*s[pw+1].r;	//nazaj
399         s[pw].g=s[pw].g+a*s[pw+1].g;
400         s[pw].b=s[pw].b+a*s[pw+1].b;
401         s[p].r=s[p].r+a*s[pmw].r;	//dol
402         s[p].g=s[p].g+a*s[pmw].g;
403         s[p].b=s[p].b+a*s[pmw].b;
404         s[pw+1].r=s[pw+1].r+a*s[p+1].r;	//dol
405         s[pw+1].g=s[pw+1].g+a*s[p+1].g;
406         s[pw+1].b=s[pw+1].b+a*s[p+1].b;
407         s[pw].r=s[pw].r+a*s[p].r;	//dol
408         s[pw].g=s[pw].g+a*s[p].g;
409         s[pw].b=s[pw].b+a*s[p].b;
410     }
411 
412     //ce je sodo stevilo vrstic, moras zadnjo posebej
413     if (i!=h)
414     {
415         p=i*w; pw=p+w;
416         for (j=1;j<w;j++)	//tja
417         {
418             s[p+j].r=s[p+j].r+a*s[p+j-1].r;
419             s[p+j].g=s[p+j].g+a*s[p+j-1].g;
420             s[p+j].b=s[p+j].b+a*s[p+j-1].b;
421         }
422 
423         s[pw-1].r=b*s[pw-1].r;	//rep H
424         s[pw-1].g=b*s[pw-1].g;
425         s[pw-1].b=b*s[pw-1].b;
426 
427         for (j=w-2;j>=0;j--)	//nazaj in dol
428         {
429             s[p+j].r=a*s[p+j+1].r+s[p+j].r;
430             s[p+j].g=a*s[p+j+1].g+s[p+j].g;
431             s[p+j].b=a*s[p+j+1].b+s[p+j].b;
432 
433             //zdaj naredi se en piksel vertikalno dol, za vse stolpce
434             //dva nazaj, da ne vpliva na H nazaj
435             s[p+j+1].r=s[p+j+1].r+a*s[p-w+j+1].r;
436             s[p+j+1].g=s[p+j+1].g+a*s[p-w+j+1].g;
437             s[p+j+1].b=s[p+j+1].b+a*s[p-w+j+1].b;
438         }
439         //levi piksel vert
440         s[p].r=s[p].r+a*s[p-w].r;
441         s[p].g=s[p].g+a*s[p-w].g;
442         s[p].b=s[p].b+a*s[p-w].b;
443     }
444 
445     //zadnja vrstica (h-1)
446     g4b=g4*b;
447     g4a=g4/(1.0-a);
448     p=(h-1)*w;
449     if (ec!=0)
450     {
451         for (i=0;i<w;i++)	//po stolpcih
452         {
453             cr=0.0;cg=0.0;cb=0.0;
454             for (j=h-avg;j<h;j++)
455             {
456                 cr=cr+s[i+w*j].r;
457                 cg=cg+s[i+w*j].g;
458                 cb=cb+s[i+w*j].b;
459             }
460             cr=cr*avg1; cg=cg*avg1; cb=cb*avg1;
461             s[i+p].r=g4a*cr+g4b*(s[i+p].r-cr);
462             s[i+p].g=g4a*cg+g4b*(s[i+p].g-cg);
463             s[i+p].b=g4a*cb+g4b*(s[i+p].b-cb);
464             outframe[p+i]=((uint32_t)s[p+i].r&0xFF) + (((uint32_t)s[p+i].g&0xFF)<<8) + (((uint32_t)s[p+i].b&0xFF)<<16);
465         }
466     }
467     else
468     {
469         for (j=0;j<w;j++)	//po stolpcih
470         {
471             s[j+p].r=g4b*s[j+p].r;	//rep V
472             s[j+p].g=g4b*s[j+p].g;
473             s[j+p].b=g4b*s[j+p].b;
474             outframe[p+j]=((uint32_t)s[p+j].r&0xFF) + (((uint32_t)s[p+j].g&0xFF)<<8) + (((uint32_t)s[p+j].b&0xFF)<<16);
475         }
476     }
477 
478     for (i=h-2;i>=0;i--)	//po vrsticah navzgor
479     {
480         p=i*w; pw=p+w;
481         for (j=0;j<w;j++)	//po stolpcih
482         {
483             s[p+j].r=a*s[pw+j].r+g4*s[p+j].r;
484             s[p+j].g=a*s[pw+j].g+g4*s[p+j].g;
485             s[p+j].b=a*s[pw+j].b+g4*s[p+j].b;
486             outframe[p+j]=((uint32_t)s[p+j].r&0xFF) + (((uint32_t)s[p+j].g&0xFF)<<8) + (((uint32_t)s[p+j].b&0xFF)<<16);
487         }
488     }
489 
490 }
491 
492 //-------------------------------------------------------
493 // 2-tap IIR v stirih smereh   a only verzija, a0=1.0
494 //desno kompenzacijo izracuna direktno (rdx,rsx,rcx)
495 //optimized for speed
fibe2o_8(const uint32_t * inframe,uint32_t * outframe,float_rgba s[],int w,int h,float a1,float a2,float rd1,float rd2,float rs1,float rs2,float rc1,float rc2,int ec)496 void fibe2o_8(const uint32_t* inframe, uint32_t* outframe, float_rgba s[], int w, int h, float a1, float a2,  float rd1, float rd2, float rs1, float rs2, float rc1, float rc2, int ec)
497 {
498     float cr,cg,cb,g,g4,avg,gavg,avgg,iavg;
499     float_rgba rep1,rep2;
500     int i,j;
501     int jw,jww,h1w,h2w,iw,i1w,i2w;
502 
503     g=1.0/(1.0+a1+a2);
504     g4=1.0/g/g/g/g;
505     avg=EDGEAVG;	//koliko vzorcev za povprecje pri edge comp
506     gavg=g4/avg;
507     avgg=1.0/g/avg;
508     iavg=1.0/avg;
509 
510     for (j=0;j<avg;j++)	//prvih avg vrstic tja in nazaj
511     {
512         jw=j*w; jww=jw+w;
513         cr=0.0;cg=0.0;cb=0.0;
514         if (ec!=0)
515         {	//edge comp (popvprecje prvih)
516             for (i=0;i<avg;i++)
517             {
518                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
519                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
520                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
521                 cr=cr+s[jw+i].r;
522                 cg=cg+s[jw+i].g;
523                 cb=cb+s[jw+i].b;
524             }
525             cr=cr*gavg; cg=cg*gavg; cb=cb*gavg;
526         }
527         else
528             for (i=0;i<avg;i++)
529             {
530                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
531                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
532                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
533             }
534 
535         s[jw].r=g4*s[jw].r-(a1+a2)*g*cr;
536         s[jw].g=g4*s[jw].g-(a1+a2)*g*cg;
537         s[jw].b=g4*s[jw].b-(a1+a2)*g*cb;
538         s[jw+1].r=g4*s[jw+1].r-a1*s[jw].r-a2*g*cr;
539         s[jw+1].g=g4*s[jw+1].g-a1*s[jw].g-a2*g*cg;
540         s[jw+1].b=g4*s[jw+1].b-a1*s[jw].b-a2*g*cb;
541 
542         if (ec!=0)
543         {	//edge comp za nazaj
544             cr=0.0;cg=0.0;cb=0.0;
545             for (i=w-avg;i<w;i++)
546             {
547                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
548                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
549                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
550                 cr=cr+s[jw+i].r;
551                 cg=cg+s[jw+i].g;
552                 cb=cb+s[jw+i].b;
553             }
554             cr=cr*gavg; cg=cg*gavg; cb=cb*gavg;
555         }
556         else
557             for (i=w-avg;i<w;i++)
558             {
559                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
560                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
561                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
562             }
563 
564         for (i=2;i<avg;i++)	//tja (ze pretv. levo)
565         {
566             s[jw+i].r=g4*s[jw+i].r-a1*s[jw+i-1].r-a2*s[jw+i-2].r;
567             s[jw+i].g=g4*s[jw+i].g-a1*s[jw+i-1].g-a2*s[jw+i-2].g;
568             s[jw+i].b=g4*s[jw+i].b-a1*s[jw+i-1].b-a2*s[jw+i-2].b;
569         }
570 
571         for (i=avg;i<w-avg;i++)	//tja (s pretvorbo)
572         {
573             s[jw+i].r=(float)(inframe[jw+i]&0xFF);
574             s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
575             s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
576             s[jw+i].r=g4*s[jw+i].r-a1*s[jw+i-1].r-a2*s[jw+i-2].r;
577             s[jw+i].g=g4*s[jw+i].g-a1*s[jw+i-1].g-a2*s[jw+i-2].g;
578             s[jw+i].b=g4*s[jw+i].b-a1*s[jw+i-1].b-a2*s[jw+i-2].b;
579         }
580 
581         for (i=w-avg;i<w;i++)	//tja (ze pretv. desno)
582         {
583             s[jw+i].r=g4*s[jw+i].r-a1*s[jw+i-1].r-a2*s[jw+i-2].r;
584             s[jw+i].g=g4*s[jw+i].g-a1*s[jw+i-1].g-a2*s[jw+i-2].g;
585             s[jw+i].b=g4*s[jw+i].b-a1*s[jw+i-1].b-a2*s[jw+i-2].b;
586         }
587 
588         rep1.r=(s[jww-1].r+s[jww-2].r)*0.5*rs1+(s[jww-1].r-s[jww-2].r)*rd1;
589         rep1.g=(s[jww-1].g+s[jww-2].g)*0.5*rs1+(s[jww-1].g-s[jww-2].g)*rd1;
590         rep1.b=(s[jww-1].b+s[jww-2].b)*0.5*rs1+(s[jww-1].b-s[jww-2].b)*rd1;
591         rep2.r=(s[jww-1].r+s[jww-2].r)*0.5*rs2+(s[jww-1].r-s[jww-2].r)*rd2;
592         rep2.g=(s[jww-1].g+s[jww-2].g)*0.5*rs2+(s[jww-1].g-s[jww-2].g)*rd2;
593         rep2.b=(s[jww-1].b+s[jww-2].b)*0.5*rs2+(s[jww-1].b-s[jww-2].b)*rd2;
594 
595         if (ec!=0)
596         {
597             rep1.r=rep1.r+rc1*cr;
598             rep1.g=rep1.g+rc1*cg;
599             rep1.b=rep1.b+rc1*cb;
600             rep2.r=rep2.r+rc2*cr;
601             rep2.g=rep2.g+rc2*cg;
602             rep2.b=rep2.b+rc2*cb;
603         }
604 
605         s[jww-1].r=s[jww-1].r-a1*rep1.r-a2*rep2.r;
606         s[jww-1].g=s[jww-1].g-a1*rep1.g-a2*rep2.g;
607         s[jww-1].b=s[jww-1].b-a1*rep1.b-a2*rep2.b;
608         s[jww-2].r=s[jww-2].r-a1*s[jww-1].r-a2*rep1.r;
609         s[jww-2].g=s[jww-2].g-a1*s[jww-1].g-a2*rep1.g;
610         s[jww-2].b=s[jww-2].b-a1*s[jww-1].b-a2*rep1.b;
611 
612         for (i=w-3;i>=0;i--)		//nazaj
613         {
614             s[jw+i].r=s[jw+i].r-a1*s[jw+i+1].r-a2*s[jw+i+2].r;
615             s[jw+i].g=s[jw+i].g-a1*s[jw+i+1].g-a2*s[jw+i+2].g;
616             s[jw+i].b=s[jw+i].b-a1*s[jw+i+1].b-a2*s[jw+i+2].b;
617         }
618     }	//prvih avg vrstic
619 
620     //edge comp zgoraj za navzdol
621     for (j=0;j<w;j++)	//po stolpcih
622     {
623         cr=0.0;cg=0.0;cb=0.0;
624         if (ec!=0)
625         {	//edge comp (popvprecje prvih)
626             for (i=0;i<avg;i++)
627             {
628                 cr=cr+s[j+w*i].r;
629                 cg=cg+s[j+w*i].g;
630                 cb=cb+s[j+w*i].b;
631             }
632             cr=cr*iavg; cg=cg*iavg; cb=cb*iavg;
633         }
634 
635         //zgornji vrstici
636         s[j].r=s[j].r-(a1+a2)*g*cr;
637         s[j].g=s[j].g-(a1+a2)*g*cg;
638         s[j].b=s[j].b-(a1+a2)*g*cb;
639         s[j+w].r=s[j+w].r-a1*s[j].r-a2*g*cr;
640         s[j+w].g=s[j+w].g-a1*s[j].g-a2*g*cg;
641         s[j+w].b=s[j+w].b-a1*s[j].b-a2*g*cb;
642     }
643 
644     //tretja do avg, samo navzdol (nazaj so ze)
645     for (i=2;i<avg;i++)
646     {
647         iw=i*w; i1w=iw-w;
648         for (j=0;j<w;j++)	//po stolpcih
649         {
650             s[j+iw].r=s[j+iw].r-a1*s[j+i1w].r-a2*s[j+i1w-w].r;
651             s[j+iw].g=s[j+iw].g-a1*s[j+i1w].g-a2*s[j+i1w-w].g;
652             s[j+iw].b=s[j+iw].b-a1*s[j+i1w].b-a2*s[j+i1w-w].b;
653         }
654     }
655 
656     for (j=avg;j<h;j++)	//po vrsticah tja, nazaj in dol
657     {
658         jw=j*w; jww=jw+w;
659         cr=0.0;cg=0.0;cb=0.0;
660         if (ec!=0)
661         {	//edge comp (popvprecje prvih)
662             for (i=0;i<avg;i++)
663             {
664                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
665                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
666                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
667                 cr=cr+s[jw+i].r;
668                 cg=cg+s[jw+i].g;
669                 cb=cb+s[jw+i].b;
670             }
671             cr=cr*gavg; cg=cg*gavg; cb=cb*gavg;
672         }
673         else
674             for (i=0;i<avg;i++)
675             {
676                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
677                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
678                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
679             }
680         s[jw].r=g4*s[jw].r-(a1+a2)*g*cr;
681         s[jw].g=g4*s[jw].g-(a1+a2)*g*cg;
682         s[jw].b=g4*s[jw].b-(a1+a2)*g*cb;
683         s[jw+1].r=g4*s[jw+1].r-a1*s[jw].r-a2*g*cr;
684         s[jw+1].g=g4*s[jw+1].g-a1*s[jw].g-a2*g*cg;
685         s[jw+1].b=g4*s[jw+1].b-a1*s[jw].b-a2*g*cb;
686 
687         if (ec!=0)
688         {	//edge comp za nazaj
689             cr=0.0;cg=0.0;cb=0.0;
690             for (i=w-avg;i<w;i++)
691             {
692                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
693                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
694                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
695                 cr=cr+s[jw+i].r;
696                 cg=cg+s[jw+i].g;
697                 cb=cb+s[jw+i].b;
698             }
699             cr=cr*gavg; cg=cg*gavg; cb=cb*gavg;
700         }
701         else
702             for (i=w-avg;i<w;i++)
703             {
704                 s[jw+i].r=(float)(inframe[jw+i]&0xFF);
705                 s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
706                 s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
707             }
708 
709         for (i=2;i<avg;i++)	//tja  (ze pretv. levo)
710         {
711             s[jw+i].r=g4*s[jw+i].r-a1*s[jw+i-1].r-a2*s[jw+i-2].r;
712             s[jw+i].g=g4*s[jw+i].g-a1*s[jw+i-1].g-a2*s[jw+i-2].g;
713             s[jw+i].b=g4*s[jw+i].b-a1*s[jw+i-1].b-a2*s[jw+i-2].b;
714         }
715 
716         for (i=avg;i<w-avg;i++)	//tja  (s retvorbo)
717         {
718             s[jw+i].r=(float)(inframe[jw+i]&0xFF);
719             s[jw+i].g=(float)((inframe[jw+i]&0xFF00)>>8);
720             s[jw+i].b=(float)((inframe[jw+i]&0xFF0000)>>16);
721             s[jw+i].r=g4*s[jw+i].r-a1*s[jw+i-1].r-a2*s[jw+i-2].r;
722             s[jw+i].g=g4*s[jw+i].g-a1*s[jw+i-1].g-a2*s[jw+i-2].g;
723             s[jw+i].b=g4*s[jw+i].b-a1*s[jw+i-1].b-a2*s[jw+i-2].b;
724         }
725 
726         for (i=w-avg;i<w;i++)	//tja  (ze pretv. desno)
727         {
728             s[jw+i].r=g4*s[jw+i].r-a1*s[jw+i-1].r-a2*s[jw+i-2].r;
729             s[jw+i].g=g4*s[jw+i].g-a1*s[jw+i-1].g-a2*s[jw+i-2].g;
730             s[jw+i].b=g4*s[jw+i].b-a1*s[jw+i-1].b-a2*s[jw+i-2].b;
731         }
732 
733         rep1.r=(s[jww-1].r+s[jww-2].r)*0.5*rs1+(s[jww-1].r-s[jww-2].r)*rd1;
734         rep1.g=(s[jww-1].g+s[jww-2].g)*0.5*rs1+(s[jww-1].g-s[jww-2].g)*rd1;
735         rep1.b=(s[jww-1].b+s[jww-2].b)*0.5*rs1+(s[jww-1].b-s[jww-2].b)*rd1;
736         rep2.r=(s[jww-1].r+s[jww-2].r)*0.5*rs2+(s[jww-1].r-s[jww-2].r)*rd2;
737         rep2.g=(s[jww-1].g+s[jww-2].g)*0.5*rs2+(s[jww-1].g-s[jww-2].g)*rd2;
738         rep2.b=(s[jww-1].b+s[jww-2].b)*0.5*rs2+(s[jww-1].b-s[jww-2].b)*rd2;
739 
740         if (ec!=0)
741         {
742             rep1.r=rep1.r+rc1*cr;
743             rep1.g=rep1.g+rc1*cg;
744             rep1.b=rep1.b+rc1*cb;
745             rep2.r=rep2.r+rc2*cr;
746             rep2.g=rep2.g+rc2*cg;
747             rep2.b=rep2.b+rc2*cb;
748         }
749 
750         s[jww-1].r=s[jww-1].r-a1*rep1.r-a2*rep2.r;
751         s[jww-1].g=s[jww-1].g-a1*rep1.g-a2*rep2.g;
752         s[jww-1].b=s[jww-1].b-a1*rep1.b-a2*rep2.b;
753         s[jww-2].r=s[jww-2].r-a1*s[jww-1].r-a2*rep1.r;
754         s[jww-2].g=s[jww-2].g-a1*s[jww-1].g-a2*rep1.g;
755         s[jww-2].b=s[jww-2].b-a1*s[jww-1].b-a2*rep1.b;
756 
757         for (i=w-3;i>=0;i--)	//po stolpcih
758         { //nazaj
759             s[jw+i].r=s[jw+i].r-a1*s[jw+i+1].r-a2*s[jw+i+2].r;
760             s[jw+i].g=s[jw+i].g-a1*s[jw+i+1].g-a2*s[jw+i+2].g;
761             s[jw+i].b=s[jw+i].b-a1*s[jw+i+1].b-a2*s[jw+i+2].b;
762             //dol
763             s[jw+i+2].r=s[jw+i+2].r-a1*s[jw-w+i+2].r-a2*s[jw-w-w+i+2].r;
764             s[jw+i+2].g=s[jw+i+2].g-a1*s[jw-w+i+2].g-a2*s[jw-w-w+i+2].g;
765             s[jw+i+2].b=s[jw+i+2].b-a1*s[jw-w+i+2].b-a2*s[jw-w-w+i+2].b;
766         }
767 
768         //se leva stolpca dol
769         s[jw+1].r=s[jw+1].r-a1*s[jw-w+1].r-a2*s[jw-w-w+1].r;
770         s[jw+1].g=s[jw+1].g-a1*s[jw-w+1].g-a2*s[jw-w-w+1].g;
771         s[jw+1].b=s[jw+1].b-a1*s[jw-w+1].b-a2*s[jw-w-w+1].b;
772         s[jw].r=s[jw].r-a1*s[jw-w].r-a2*s[jw-w-w].r;
773         s[jw].g=s[jw].g-a1*s[jw-w].g-a2*s[jw-w-w].g;
774         s[jw].b=s[jw].b-a1*s[jw-w].b-a2*s[jw-w-w].b;
775 
776     }	//po vrsticah
777 
778     //pa se navzgor
779     //spodnji dve vrstici
780     h1w=(h-1)*w; h2w=(h-2)*w;
781     for (j=0;j<w;j++)	//po stolpcih
782     {
783         if (ec!=0)
784         {	//edge comp za gor
785             cr=0.0;cg=0.0;cb=0.0;
786             for (i=h-avg;i<h;i++)
787             {
788                 cr=cr+s[j+w*i].r;
789                 cg=cg+s[j+w*i].g;
790                 cb=cb+s[j+w*i].b;
791             }
792             cr=cr*avgg; cg=cg*avgg; cb=cb*avgg;
793         }
794 
795         rep1.r=(s[j+h1w].r+s[j+h2w].r)*0.5*rs1+(s[j+h1w].r-s[j+h2w].r)*rd1;
796         rep1.g=(s[j+h1w].g+s[j+h2w].g)*0.5*rs1+(s[j+h1w].g-s[j+h2w].g)*rd1;
797         rep1.b=(s[j+h1w].b+s[j+h2w].b)*0.5*rs1+(s[j+h1w].b-s[j+h2w].b)*rd1;
798         rep2.r=(s[j+h1w].r+s[j+h2w].r)*0.5*rs2+(s[j+h1w].r-s[j+h2w].r)*rd2;
799         rep2.g=(s[j+h1w].g+s[j+h2w].g)*0.5*rs2+(s[j+h1w].g-s[j+h2w].g)*rd2;
800         rep2.b=(s[j+h1w].b+s[j+h2w].b)*0.5*rs2+(s[j+h1w].b-s[j+h2w].b)*rd2;
801 
802         if (ec!=0)
803         {	//edge comp
804             rep1.r=rep1.r+rc1*cr;
805             rep1.g=rep1.g+rc1*cg;
806             rep1.b=rep1.b+rc1*cb;
807             rep2.r=rep2.r+rc2*cr;
808             rep2.g=rep2.g+rc2*cg;
809             rep2.b=rep2.b+rc2*cb;
810         }
811 
812         s[j+h1w].r=s[j+h1w].r-a1*rep1.r-a2*rep2.r;
813         s[j+h1w].g=s[j+h1w].g-a1*rep1.g-a2*rep2.g;
814         s[j+h1w].b=s[j+h1w].b-a1*rep1.b-a2*rep2.b;
815         if (s[j+h1w].r>255) s[j+h1w].r=255.0;
816         if (s[j+h1w].r<0.0) s[j+h1w].r=0.0;
817         if (s[j+h1w].g>255) s[j+h1w].g=255.0;
818         if (s[j+h1w].g<0.0) s[j+h1w].g=0.0;
819         if (s[j+h1w].b>255) s[j+h1w].b=255.0;
820         if (s[j+h1w].b<0.0) s[j+h1w].b=0.0;
821         outframe[j+h1w]=((uint32_t)s[j+h1w].r&0xFF) + (((uint32_t)s[j+h1w].g&0xFF)<<8) + (((uint32_t)s[j+h1w].b&0xFF)<<16);
822         s[j+h2w].r=s[j+h2w].r-a1*s[j+h1w].r-a2*rep1.r;
823         s[j+h2w].g=s[j+h2w].g-a1*s[j+h1w].g-a2*rep1.g;
824         s[j+h2w].b=s[j+h2w].b-a1*s[j+h1w].b-a2*rep1.b;
825         if (s[j+h2w].r>255) s[j+h2w].r=255.0;
826         if (s[j+h2w].r<0.0) s[j+h2w].r=0.0;
827         if (s[j+h2w].g>255) s[j+h2w].g=255.0;
828         if (s[j+h2w].g<0.0) s[j+h2w].g=0.0;
829         if (s[j+h2w].b>255) s[j+h2w].b=255.0;
830         if (s[j+h2w].b<0.0) s[j+h2w].b=0.0;
831         outframe[j+h2w]=((uint32_t)s[j+h2w].r&0xFF) + (((uint32_t)s[j+h2w].g&0xFF)<<8) + (((uint32_t)s[j+h2w].b&0xFF)<<16);
832     }
833 
834     //ostale vrstice
835     for (i=h-3;i>=0;i--)		//gor
836     {
837         iw=i*w; i1w=iw+w; i2w=i1w+w;
838         for (j=0;j<w;j++)
839         {
840             s[j+iw].r=s[j+iw].r-a1*s[j+i1w].r-a2*s[j+i2w].r;
841             s[j+iw].g=s[j+iw].g-a1*s[j+i1w].g-a2*s[j+i2w].g;
842             s[j+iw].b=s[j+iw].b-a1*s[j+i1w].b-a2*s[j+i2w].b;
843             if (s[j+iw].r>255) s[j+iw].r=255.0;
844             if (s[j+iw].r<0.0) s[j+iw].r=0.0;
845             if (s[j+iw].g>255) s[j+iw].g=255.0;
846             if (s[j+iw].g<0.0) s[j+iw].g=0.0;
847             if (s[j+iw].b>255) s[j+iw].b=255.0;
848             if (s[j+iw].b<0.0) s[j+iw].b=0.0;
849             outframe[j+iw]=((uint32_t)s[j+iw].r&0xFF) + (((uint32_t)s[j+iw].g&0xFF)<<8) + (((uint32_t)s[j+iw].b&0xFF)<<16);
850         }
851     }
852 
853 }
854 
855 //-------------------------------------------------------
856 // 3-tap IIR v stirih smereh
857 //a only verzija, a0=1.0
858 //edge efekt na desni kompenzira tako, da racuna 256 vzorcev
859 //cez rob in in gre potem nazaj
fibe3_8(const uint32_t * inframe,uint32_t * outframe,float_rgba s[],int w,int h,float a1,float a2,float a3,int ec)860 void fibe3_8(const uint32_t* inframe, uint32_t* outframe, float_rgba s[], int w, int h, float a1, float a2, float a3, int ec)
861 {
862     float cr,cg,cb,g,g4;
863     int i,j;
864     const float avg = EDGEAVG; // how many samples for average at edge comp
865     const int cez = 256; // how many samples go right
866     float_rgba *lb = malloc((MAX(w, h) + cez) * sizeof(*lb));
867 
868     g=1.0/(1.0+a1+a2+a3); g4=1.0/g/g/g/g;
869 
870     for (j=0;j<h;j++)	//po vrsticah
871     {
872         cr=0.0;cg=0.0;cb=0.0;
873         if (ec!=0)
874         {	//edge comp (popvprecje prvih)
875             for (i=0;i<avg;i++)
876             {
877                 s[j*w+i].r=(float)(inframe[j*w+i]&0xFF);
878                 s[j*w+i].g=(float)((inframe[j*w+i]&0xFF00)>>8);
879                 s[j*w+i].b=(float)((inframe[j*w+i]&0xFF0000)>>16);
880                 cr=cr+s[j*w+i].r;
881                 cg=cg+s[j*w+i].g;
882                 cb=cb+s[j*w+i].b;
883             }
884             cr=g4*cr/avg; cg=g4*cg/avg; cb=g4*cb/avg;
885         }
886         else
887             for (i=0;i<avg;i++)
888             {
889                 s[j*w+i].r=(float)(inframe[j*w+i]&0xFF);
890                 s[j*w+i].g=(float)((inframe[j*w+i]&0xFF00)>>8);
891                 s[j*w+i].b=(float)((inframe[j*w+i]&0xFF0000)>>16);
892             }
893         lb[0].r=g4*s[j*w].r-(a1+a2+a3)*g*cr;
894         lb[0].g=g4*s[j*w].g-(a1+a2+a3)*g*cg;
895         lb[0].b=g4*s[j*w].b-(a1+a2+a3)*g*cb;
896         lb[1].r=g4*s[j*w+1].r-a1*lb[0].r-(a2+a3)*g*cr;
897         lb[1].g=g4*s[j*w+1].g-a1*lb[0].g-(a2+a3)*g*cg;
898         lb[1].b=g4*s[j*w+1].b-a1*lb[0].b-(a2+a3)*g*cb;
899         lb[2].r=g4*s[j*w+2].r-a1*lb[1].r-a2*lb[0].r-a3*g*cr;
900         lb[2].g=g4*s[j*w+2].g-a1*lb[1].g-a2*lb[0].g-a3*g*cg;
901         lb[2].b=g4*s[j*w+2].b-a1*lb[1].b-a2*lb[0].b-a3*g*cb;
902 
903         for (i=3;i<avg;i++)	//tja  (ze pretvorjeni)
904         {
905             lb[i].r=g4*s[j*w+i].r-a1*lb[i-1].r-a2*lb[i-2].r-a3*lb[i-3].r;
906             lb[i].g=g4*s[j*w+i].g-a1*lb[i-1].g-a2*lb[i-2].g-a3*lb[i-3].g;
907             lb[i].b=g4*s[j*w+i].b-a1*lb[i-1].b-a2*lb[i-2].b-a3*lb[i-3].b;
908         }
909 
910         for (i=avg;i<w;i++)	//tja  (s pretvorbo)
911         {
912             s[j*w+i].r=(float)(inframe[j*w+i]&0xFF);
913             s[j*w+i].g=(float)((inframe[j*w+i]&0xFF00)>>8);
914             s[j*w+i].b=(float)((inframe[j*w+i]&0xFF0000)>>16);
915             lb[i].r=g4*s[j*w+i].r-a1*lb[i-1].r-a2*lb[i-2].r-a3*lb[i-3].r;
916             lb[i].g=g4*s[j*w+i].g-a1*lb[i-1].g-a2*lb[i-2].g-a3*lb[i-3].g;
917             lb[i].b=g4*s[j*w+i].b-a1*lb[i-1].b-a2*lb[i-2].b-a3*lb[i-3].b;
918         }
919 
920         cr=0.0;cg=0.0;cb=0.0;
921         if (ec!=0)
922         {	//edge comp
923             for (i=w-avg;i<w;i++)
924             {
925                 cr=cr+s[j*w+i].r;
926                 cg=cg+s[j*w+i].g;
927                 cb=cb+s[j*w+i].b;
928             }
929             cr=g4*cr/avg; cg=g4*cg/avg; cb=g4*cb/avg;
930         }
931 
932         for (i=w;i<(w+cez);i++)	//naprej cez rob
933         {
934             lb[i].r=cr-a1*lb[i-1].r-a2*lb[i-2].r-a3*lb[i-3].r;
935             lb[i].g=cr-a1*lb[i-1].g-a2*lb[i-2].g-a3*lb[i-3].g;
936             lb[i].b=cr-a1*lb[i-1].b-a2*lb[i-2].b-a3*lb[i-3].b;
937         }
938         //nazaj do roba
939         lb[w+cez-2].r=lb[w+cez-2].r-a1*lb[w+cez-1].r;
940         lb[w+cez-2].g=lb[w+cez-2].g-a1*lb[w+cez-1].g;
941         lb[w+cez-2].b=lb[w+cez-2].b-a1*lb[w+cez-1].b;
942         lb[w+cez-3].r=lb[w+cez-3].r-a1*lb[w+cez-2].r-a2*lb[w+cez-1].r;
943         lb[w+cez-3].g=lb[w+cez-3].g-a1*lb[w+cez-2].g-a2*lb[w+cez-1].g;
944         lb[w+cez-3].b=lb[w+cez-3].b-a1*lb[w+cez-2].b-a2*lb[w+cez-1].b;
945         for (i=(w+cez-4);i>=w;i--)
946         {
947             lb[i].r=lb[i].r-a1*lb[i+1].r-a2*lb[i+2].r-a3*lb[i+3].r;
948             lb[i].g=lb[i].g-a1*lb[i+1].g-a2*lb[i+2].g-a3*lb[i+3].g;
949             lb[i].b=lb[i].b-a1*lb[i+1].b-a2*lb[i+2].b-a3*lb[i+3].b;
950         }
951 
952         s[j*w+w-1].r=lb[w-1].r-a1*lb[w].r-a2*lb[w+1].r-a3*lb[w+2].r;
953         s[j*w+w-1].g=lb[w-1].g-a1*lb[w].g-a2*lb[w+1].g-a3*lb[w+2].g;
954         s[j*w+w-1].b=lb[w-1].b-a1*lb[w].b-a2*lb[w+1].b-a3*lb[w+2].b;
955         s[j*w+w-2].r=lb[w-2].r-a1*s[j*w+w-1].r-a2*lb[w].r-a3*lb[w+1].r;
956         s[j*w+w-2].g=lb[w-2].g-a1*s[j*w+w-1].g-a2*lb[w].g-a3*lb[w+1].g;
957         s[j*w+w-2].b=lb[w-2].b-a1*s[j*w+w-1].b-a2*lb[w].b-a3*lb[w+1].b;
958         s[j*w+w-3].r=lb[w-3].r-a1*s[j*w+w-2].r-a2*s[j*w+w-1].r-a3*lb[w].r;
959         s[j*w+w-3].g=lb[w-3].g-a1*s[j*w+w-2].g-a2*s[j*w+w-1].g-a3*lb[w].g;
960         s[j*w+w-3].b=lb[w-3].b-a1*s[j*w+w-2].b-a2*s[j*w+w-1].b-a3*lb[w].b;
961 
962         for (i=w-4;i>=0;i--)		//nazaj
963         {
964             s[j*w+i].r=lb[i].r-a1*s[j*w+i+1].r-a2*s[j*w+i+2].r-a3*s[j*w+i+3].r;
965             s[j*w+i].g=lb[i].g-a1*s[j*w+i+1].g-a2*s[j*w+i+2].g-a3*s[j*w+i+3].g;
966             s[j*w+i].b=lb[i].b-a1*s[j*w+i+1].b-a2*s[j*w+i+2].b-a3*s[j*w+i+3].b;
967         }
968     }	//po vrsticah
969 
970     for (j=0;j<w;j++)	//po stolpcih
971     {
972         cr=0.0;cg=0.0;cb=0.0;
973         if (ec!=0)
974         {	//edge comp (popvprecje prvih)
975             for (i=0;i<avg;i++)
976             {
977                 cr=cr+s[j+w*i].r;
978                 cg=cg+s[j+w*i].g;
979                 cb=cb+s[j+w*i].b;
980             }
981             cr=cr/avg; cg=cg/avg; cb=cb/avg;
982         }
983         lb[0].r=s[j].r-(a1+a2+a3)*g*cr;
984         lb[0].g=s[j].g-(a1+a2+a3)*g*cg;
985         lb[0].b=s[j].b-(a1+a2+a3)*g*cb;
986         lb[1].r=s[j+w].r-a1*lb[0].r-(a2+a3)*g*cr;
987         lb[1].g=s[j+w].g-a1*lb[0].g-(a2+a3)*g*cg;
988         lb[1].b=s[j+w].b-a1*lb[0].b-(a2+a3)*g*cb;
989         lb[2].r=s[j+2*w].r-a1*lb[1].r-a2*lb[0].r-a3*g*cr;
990         lb[2].g=s[j+2*w].g-a1*lb[1].g-a2*lb[0].g-a3*g*cg;
991         lb[2].b=s[j+2*w].b-a1*lb[1].b-a2*lb[0].b-a3*g*cb;
992 
993         for (i=3;i<h;i++)		//dol
994         {
995             lb[i].r=s[j+w*i].r-a1*lb[i-1].r-a2*lb[i-2].r-a3*lb[i-3].r;
996             lb[i].g=s[j+w*i].g-a1*lb[i-1].g-a2*lb[i-2].g-a3*lb[i-3].g;
997             lb[i].b=s[j+w*i].b-a1*lb[i-1].b-a2*lb[i-2].b-a3*lb[i-3].b;
998         }
999 
1000         cr=0.0;cg=0.0;cb=0.0;
1001         if (ec!=0)
1002         {	//edge comp
1003             for (i=h-avg;i<h;i++)
1004             {
1005                 cr=cr+s[j+w*i].r;
1006                 cg=cg+s[j+w*i].g;
1007                 cb=cb+s[j+w*i].b;
1008             }
1009             cr=cr/avg; cg=cg/avg; cb=cb/avg;
1010         }
1011 
1012         for (i=h;i<(h+cez);i++)	//naprej cez rob
1013         {
1014             lb[i].r=cr-a1*lb[i-1].r-a2*lb[i-2].r-a3*lb[i-3].r;
1015             lb[i].g=cg-a1*lb[i-1].g-a2*lb[i-2].g-a3*lb[i-3].g;
1016             lb[i].b=cb-a1*lb[i-1].b-a2*lb[i-2].b-a3*lb[i-3].b;
1017         }
1018         //nazaj do roba
1019         lb[h+cez-2].r=lb[h+cez-2].r-a1*lb[h+cez-1].r;
1020         lb[h+cez-2].g=lb[h+cez-2].g-a1*lb[h+cez-1].g;
1021         lb[h+cez-2].b=lb[h+cez-2].b-a1*lb[h+cez-1].b;
1022         lb[h+cez-3].r=lb[h+cez-3].r-a1*lb[h+cez-2].r-a2*lb[h+cez-1].r;
1023         lb[h+cez-3].g=lb[h+cez-3].g-a1*lb[h+cez-2].g-a2*lb[h+cez-1].g;
1024         lb[h+cez-3].b=lb[h+cez-3].b-a1*lb[h+cez-2].b-a2*lb[h+cez-1].b;
1025         for (i=(h+cez-4);i>=h;i--)
1026         {
1027             lb[i].r=lb[i].r-a1*lb[i+1].r-a2*lb[i+2].r-a3*lb[i+3].r;
1028             lb[i].g=lb[i].g-a1*lb[i+1].g-a2*lb[i+2].g-a3*lb[i+3].g;
1029             lb[i].b=lb[i].b-a1*lb[i+1].b-a2*lb[i+2].b-a3*lb[i+3].b;
1030         }
1031 
1032         s[j+(h-1)*w].r=lb[h-1].r-a1*lb[h].r-a2*lb[h+1].r-a3*lb[h+2].r;
1033         s[j+(h-1)*w].g=lb[h-1].g-a1*lb[h].g-a2*lb[h+1].g-a3*lb[h+2].g;
1034         s[j+(h-1)*w].b=lb[h-1].b-a1*lb[h].b-a2*lb[h+1].b-a3*lb[h+2].b;
1035         s[j+(h-2)*w].r=lb[h-2].r-a1*s[j+(h-1)*w].r-a2*lb[h].r-a3*lb[h+1].r;
1036         s[j+(h-2)*w].g=lb[h-2].g-a1*s[j+(h-1)*w].g-a2*lb[h].g-a3*lb[h+1].g;
1037         s[j+(h-2)*w].b=lb[h-2].b-a1*s[j+(h-1)*w].b-a2*lb[h].b-a3*lb[h+1].b;
1038         s[j+(h-3)*w].r=lb[h-3].r-a1*s[j+(h-2)*w].r-a2*s[j+(h-1)*w].r-a3*lb[h].r;
1039         s[j+(h-3)*w].g=lb[h-3].g-a1*s[j+(h-2)*w].g-a2*s[j+(h-1)*w].g-a3*lb[h].g;
1040         s[j+(h-3)*w].b=lb[h-3].b-a1*s[j+(h-2)*w].b-a2*s[j+(h-1)*w].b-a3*lb[h].b;
1041 
1042         for (i=h-4;i>=0;i--)		//gor
1043         {
1044             s[j+w*i].r=lb[i].r-a1*s[j+w*(i+1)].r-a2*s[j+w*(i+2)].r-a3*s[j+w*(i+3)].r;
1045             s[j+w*i].g=lb[i].g-a1*s[j+w*(i+1)].g-a2*s[j+w*(i+2)].g-a3*s[j+w*(i+3)].g;
1046             s[j+w*i].b=lb[i].b-a1*s[j+w*(i+1)].b-a2*s[j+w*(i+2)].b-a3*s[j+w*(i+3)].b;
1047             outframe[j+w*i]=((uint32_t)s[j+w*i].r&0xFF) + (((uint32_t)s[j+w*i].g&0xFF)<<8) + (((uint32_t)s[j+w*i].b&0xFF)<<16);
1048         }
1049     }	//po stolpcih
1050 
1051     free(lb);
1052 }
1053