1 /*
2  *  transformtype.c
3  *
4  *  Copyright (C) Georg Martius - June 2007
5  *
6  *  This file is part of transcode, a video stream processing tool
7  *
8  *  transcode is free software; you can redistribute it and/or modify
9  *  it under the terms of the GNU General Public License as published by
10  *  the Free Software Foundation; either version 2, or (at your option)
11  *  any later version.
12  *
13  *  transcode is distributed in the hope that it will be useful,
14  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16  *  GNU General Public License for more details.
17  *
18  *  You should have received a copy of the GNU General Public License
19  *  along with GNU Make; see the file COPYING.  If not, write to
20  *  the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
21  *
22  */
23 #include <stdlib.h>
24 #include <math.h>
25 #include <string.h>
26 #include "transformtype.h"
27 #include "transformtype_operations.h"
28 #include "vidstabdefines.h"
29 
30 /***********************************************************************
31  * helper functions to create and operate with transforms.
32  * all functions are non-destructive
33  */
34 
35 /* create an initialized transform*/
new_transform(double x,double y,double alpha,double zoom,double barrel,double rshutter,int extra)36 VSTransform new_transform(double x, double y, double alpha,
37                           double zoom, double barrel, double rshutter, int extra)
38 {
39   VSTransform t;
40   t.x        = x;
41   t.y        = y;
42   t.alpha    = alpha;
43   t.zoom     = zoom;
44   t.barrel   = barrel;
45   t.rshutter = rshutter;
46   t.extra    = extra;
47   return t;
48 }
49 
50 /* create a zero initialized transform*/
null_transform(void)51 VSTransform null_transform(void)
52 {
53   return new_transform(0, 0, 0, 0, 0, 0, 0);
54 }
55 
56 /* adds two transforms */
add_transforms(const VSTransform * t1,const VSTransform * t2)57 VSTransform add_transforms(const VSTransform* t1, const VSTransform* t2)
58 {
59   VSTransform t;
60   t.x        = t1->x + t2->x;
61   t.y        = t1->y + t2->y;
62   t.alpha    = t1->alpha + t2->alpha;
63   t.zoom     = t1->zoom + t2->zoom;
64   t.barrel   = t1->barrel + t2->barrel;
65   t.rshutter = t1->rshutter + t2->rshutter;
66   t.extra    = t1->extra || t2->extra;
67   return t;
68 }
69 
70 /* like add_transform but with non-pointer signature */
add_transforms_(const VSTransform t1,const VSTransform t2)71 VSTransform add_transforms_(const VSTransform t1, const VSTransform t2)
72 {
73   return add_transforms(&t1, &t2);
74 }
75 
76 /* subtracts two transforms */
sub_transforms(const VSTransform * t1,const VSTransform * t2)77 VSTransform sub_transforms(const VSTransform* t1, const VSTransform* t2)
78 {
79   VSTransform t;
80   t.x        = t1->x - t2->x;
81   t.y        = t1->y - t2->y;
82   t.alpha    = t1->alpha - t2->alpha;
83   t.zoom     = t1->zoom - t2->zoom;
84   t.barrel   = t1->barrel - t2->barrel;
85   t.rshutter = t1->rshutter - t2->rshutter;
86   t.extra    = t1->extra || t2->extra;
87   return t;
88 }
89 
90 /* multiplies a transforms with a scalar */
mult_transform(const VSTransform * t1,double f)91 VSTransform mult_transform(const VSTransform* t1, double f)
92 {
93   VSTransform t;
94   t.x        = t1->x        * f;
95   t.y        = t1->y        * f;
96   t.alpha    = t1->alpha    * f;
97   t.zoom     = t1->zoom     * f;
98   t.barrel   = t1->barrel   * f;
99   t.rshutter = t1->rshutter * f;
100   t.extra    = t1->extra;
101   return t;
102 }
103 
104 /* like mult_transform but with non-pointer signature */
mult_transform_(const VSTransform t1,double f)105 VSTransform mult_transform_(const VSTransform t1, double f)
106 {
107   return mult_transform(&t1,f);
108 }
109 
storeVSTransform(FILE * f,const VSTransform * t)110 void storeVSTransform(FILE* f, const VSTransform* t){
111   fprintf(f,"Trans %lf %lf %lf %lf %i\n", t->x, t->y, t->alpha, t->zoom, t->extra);
112 }
113 
114 
prepare_transform(const VSTransform * t,const VSFrameInfo * fi)115 PreparedTransform prepare_transform(const VSTransform* t, const VSFrameInfo* fi){
116   PreparedTransform pt;
117   pt.t = t;
118   double z = 1.0+t->zoom/100.0;
119   pt.zcos_a = z*cos(t->alpha); // scaled cos
120   pt.zsin_a = z*sin(t->alpha); // scaled sin
121   pt.c_x    = fi->width / 2;
122   pt.c_y    = fi->height / 2;
123   return pt;
124 }
125 
transform_vec(const PreparedTransform * pt,const Vec * v)126 Vec transform_vec(const PreparedTransform* pt, const Vec* v){
127   double x,y;
128   transform_vec_double(&x, &y, pt, v);
129   Vec res = {x,y};
130   return res;
131 }
132 
transform_vec_double(double * x,double * y,const PreparedTransform * pt,const Vec * v)133 void transform_vec_double(double* x, double* y, const PreparedTransform* pt, const Vec* v){
134   double rx = v->x - pt->c_x;
135   double ry = v->y - pt->c_y;
136   *x =  pt->zcos_a * rx + pt->zsin_a * ry + pt->t->x + pt->c_x;
137   *y = -pt->zsin_a * rx + pt->zcos_a * ry + pt->t->y + pt->c_y;
138 }
139 
sub_vec(Vec v1,Vec v2)140 Vec sub_vec(Vec v1, Vec v2){
141   Vec r = {v1.x - v2.x, v1.y - v2.y};
142   return r;
143 }
add_vec(Vec v1,Vec v2)144 Vec add_vec(Vec v1, Vec v2){
145   Vec r = {v1.x + v2.x, v1.y + v2.y};
146   return r;
147 }
field_to_vec(Field f)148 Vec field_to_vec(Field f){
149   Vec r = {f.x , f.y};
150   return r;
151 }
152 
153 /* compares a transform with respect to x (for sort function) */
cmp_trans_x(const void * t1,const void * t2)154 int cmp_trans_x(const void *t1, const void* t2)
155 {
156   double a = ((VSTransform*)t1)->x;
157   double b = ((VSTransform*)t2)->x;
158   return a < b ? -1 : ( a > b ? 1 : 0 );
159 }
160 
161 /* compares a transform with respect to y (for sort function) */
cmp_trans_y(const void * t1,const void * t2)162 int cmp_trans_y(const void *t1, const void* t2)
163 {
164   double a = ((VSTransform*)t1)->y;
165   double b = ((VSTransform*)t2)->y;
166   return a < b ? -1 : ( a > b ? 1: 0 );
167 }
168 
169 /* static int cmp_trans_alpha(const void *t1, const void* t2){ */
170 /*   double a = ((VSTransform*)t1)->alpha; */
171 /*   double b = ((VSTransform*)t2)->alpha; */
172 /*   return a < b ? -1 : ( a > b ? 1 : 0 ); */
173 /* } */
174 
175 
176 /* compares two double values (for sort function)*/
cmp_double(const void * t1,const void * t2)177 int cmp_double(const void *t1, const void* t2)
178 {
179   double a = *((double*)t1);
180   double b = *((double*)t2);
181   return a < b ? -1 : ( a > b ? 1 : 0 );
182 }
183 
184 /* compares two int values (for sort function)*/
cmp_int(const void * t1,const void * t2)185 int cmp_int(const void *t1, const void* t2)
186 {
187   int a = *((int*)t1);
188   int b = *((int*)t2);
189   return a < b ? -1 : ( a > b ? 1 : 0 );
190 }
191 
192 /**
193  * median_xy_transform: calulcates the median of an array
194  * of transforms, considering only x and y
195  *
196  * Parameters:
197  *    transforms: array of transforms.
198  *           len: length  of array
199  * Return value:
200  *     A new transform with x and y beeing the median of
201  *     all transforms. alpha and other fields are 0.
202  * Preconditions:
203  *     len>0
204  * Side effects:
205  *     None
206  */
median_xy_transform(const VSTransform * transforms,int len)207 VSTransform median_xy_transform(const VSTransform* transforms, int len)
208 {
209   VSTransform* ts = vs_malloc(sizeof(VSTransform) * len);
210   VSTransform t   = null_transform();
211   memcpy(ts,transforms, sizeof(VSTransform)*len );
212   int half = len/2;
213   qsort(ts, len, sizeof(VSTransform), cmp_trans_x);
214   t.x = len % 2 == 0 ? ts[half].x : (ts[half].x + ts[half+1].x)/2;
215   qsort(ts, len, sizeof(VSTransform), cmp_trans_y);
216   t.y = len % 2 == 0 ? ts[half].y : (ts[half].y + ts[half+1].y)/2;
217   vs_free(ts);
218   return t;
219 }
220 
221 /**
222  * cleanmean_xy_transform: calulcates the cleaned mean of an array
223  * of transforms, considering only x and y
224  *
225  * Parameters:
226  *    transforms: array of transforms.
227  *           len: length  of array
228  * Return value:
229  *     A new transform with x and y beeing the cleaned mean
230  *     (meaning upper and lower pentile are removed) of
231  *     all transforms. alpha and other fields are 0.
232  * Preconditions:
233  *     len>0
234  * Side effects:
235  *     None
236  */
cleanmean_xy_transform(const VSTransform * transforms,int len)237 VSTransform cleanmean_xy_transform(const VSTransform* transforms, int len)
238 {
239   VSTransform* ts = vs_malloc(sizeof(VSTransform) * len);
240   VSTransform t = null_transform();
241   int i, cut = len / 5;
242   memcpy(ts, transforms, sizeof(VSTransform) * len);
243   qsort(ts,len, sizeof(VSTransform), cmp_trans_x);
244   for (i = cut; i < len - cut; i++){ // all but cutted
245     t.x += ts[i].x;
246   }
247   qsort(ts, len, sizeof(VSTransform), cmp_trans_y);
248   for (i = cut; i < len - cut; i++){ // all but cutted
249     t.y += ts[i].y;
250   }
251   vs_free(ts);
252   return mult_transform(&t, 1.0 / (len - (2.0 * cut)));
253 }
254 
255 
256 /**
257  * calulcates the cleaned maximum and minimum of an array of transforms,
258  * considerung only x and y
259  * It cuts off the upper and lower x-th percentil
260  *
261  * Parameters:
262  *    transforms: array of transforms.
263  *           len: length  of array
264  *     percentil: the x-th percentil to cut off
265  *           min: pointer to min (return value)
266  *           max: pointer to max (return value)
267  * Return value:
268  *     call by reference in min and max
269  * Preconditions:
270  *     len>0, 0<=percentil<50
271  * Side effects:
272  *     only on min and max
273  */
cleanmaxmin_xy_transform(const VSTransform * transforms,int len,int percentil,VSTransform * min,VSTransform * max)274 void cleanmaxmin_xy_transform(const VSTransform* transforms, int len,
275                               int percentil,
276                               VSTransform* min, VSTransform* max){
277   VSTransform* ts = vs_malloc(sizeof(VSTransform) * len);
278   int cut = len * percentil / 100;
279   memcpy(ts, transforms, sizeof(VSTransform) * len);
280   qsort(ts,len, sizeof(VSTransform), cmp_trans_x);
281   min->x = ts[cut].x;
282   max->x = ts[len-cut-1].x;
283   qsort(ts, len, sizeof(VSTransform), cmp_trans_y);
284   min->y = ts[cut].y;
285   max->y = ts[len-cut-1].y;
286   vs_free(ts);
287 }
288 
289 /* calculates the required zoom value to have no borders visible
290  */
transform_get_required_zoom(const VSTransform * transform,int width,int height)291 double transform_get_required_zoom(const VSTransform* transform, int width, int height){
292   return 100.0*(2.0*VS_MAX(fabs(transform->x)/width,fabs(transform->y)/height)  // translation part
293                 + fabs(sin(transform->alpha)));          // rotation part
294 
295 }
296 
297 
298 /**
299  * media: median of a double array
300  *
301  * Parameters:
302  *            ds: array of values
303  *           len: length  of array
304  * Return value:
305  *     the median value of the array
306  * Preconditions: len>0
307  * Side effects:  ds will be sorted!
308  */
median(double * ds,int len)309 double median(double* ds, int len)
310 {
311   int half=len/2;
312   qsort(ds,len, sizeof(double), cmp_double);
313   return len % 2 == 0 ? ds[half] : (ds[half] + ds[half+1])/2;
314 }
315 
316 
317 /** square of a number */
sqr(double x)318 double sqr(double x){ return x*x; }
319 
320 /**
321  * mean: mean of a double array
322  *
323  * Parameters:
324  *            ds: array of values
325  *           len: length  of array
326  * Return value: the mean value of the array
327  * Preconditions: len>0
328  * Side effects:  None
329  */
mean(const double * ds,int len)330 double mean(const double* ds, int len)
331 {
332   double sum=0;
333   int i = 0;
334   for (i = 0; i < len; i++)
335     sum += ds[i];
336   return sum / len;
337 }
338 
339 /**
340  * stddev: standard deviation of a double array
341  *
342  * Parameters:
343  *            ds: array of values
344  *           len: length  of array
345  *          mean: mean of the array (@see mean())
346  * Return value: the standard deviation value of the array
347  * Preconditions: len>0
348  * Side effects:  None
349  */
stddev(const double * ds,int len,double mean)350 double stddev(const double* ds, int len, double mean)
351 {
352   double sum=0;
353   int i = 0;
354   for (i = 0; i < len; i++)
355     sum += sqr(ds[i]-mean);
356   return sqrt(sum / len);
357 }
358 
359 /**
360  * cleanmean: mean with cutted upper and lower pentile
361  *
362  * Parameters:
363  *            ds: array of values
364  *           len: length  of array
365  *       minimum: minimal value (after cleaning) if not NULL
366  *       maximum: maximal value (after cleaning) if not NULL
367  * Return value:
368  *     the mean value of the array without the upper
369  *     and lower pentile (20% each)
370  *     and minimum and maximum without the pentiles
371  * Preconditions: len>0
372  * Side effects:  ds will be sorted!
373  */
cleanmean(double * ds,int len,double * minimum,double * maximum)374 double cleanmean(double* ds, int len, double* minimum, double* maximum)
375 {
376   int cut    = len / 5;
377   double sum = 0;
378   int i      = 0;
379   qsort(ds, len, sizeof(double), cmp_double);
380   for (i = cut; i < len - cut; i++) { // all but first and last
381     sum += ds[i];
382   }
383   if (minimum)
384     *minimum = ds[cut];
385   if (maximum)
386     *maximum = ds[len-cut-1];
387   return sum / (len - (2.0 * cut));
388 }
389 
390 /************************************************/
391 /***************LOCALMOTION**********************/
392 
null_localmotion()393 LocalMotion null_localmotion(){
394   LocalMotion lm;
395   memset(&lm,0,sizeof(lm));
396   return lm;
397 }
398 
localmotions_getx(const LocalMotions * localmotions)399 int* localmotions_getx(const LocalMotions* localmotions){
400   int len = vs_vector_size(localmotions);
401   int* xs = vs_malloc(sizeof(int) * len);
402   int i;
403   for (i=0; i<len; i++){
404     xs[i]=LMGet(localmotions,i)->v.x;
405   }
406   return xs;
407 }
408 
localmotions_gety(const LocalMotions * localmotions)409 int* localmotions_gety(const LocalMotions* localmotions){
410   int len = vs_vector_size(localmotions);
411   int* ys = vs_malloc(sizeof(int) * len);
412   int i;
413   for (i=0; i<len; i++){
414     ys[i]=LMGet(localmotions,i)->v.y;
415   }
416   return ys;
417 }
418 
sub_localmotion(const LocalMotion * lm1,const LocalMotion * lm2)419 LocalMotion sub_localmotion(const LocalMotion* lm1, const LocalMotion* lm2){
420   LocalMotion res = *lm1;
421   res.v.x -= lm2->v.x;
422   res.v.y -= lm2->v.y;
423   return res;
424 }
425 
426 
427 /**
428  * cleanmean_localmotions: calulcates the cleaned mean of a vector
429  * of local motions considering
430  *
431  * Parameters:
432  *    localmotions : vs_vector of local motions
433  * Return value:
434  *     A localmotion with vec with x and y being the cleaned mean
435  *     (meaning upper and lower pentile are removed) of
436  *     all local motions. all other fields are 0.
437  * Preconditions:
438  *     size of vector >0
439  * Side effects:
440  *     None
441  */
cleanmean_localmotions(const LocalMotions * localmotions)442 LocalMotion cleanmean_localmotions(const LocalMotions* localmotions)
443 {
444   int len = vs_vector_size(localmotions);
445   int i, cut = len / 5;
446   int* xs = localmotions_getx(localmotions);
447   int* ys = localmotions_gety(localmotions);
448   LocalMotion m = null_localmotion();
449   m.v.x=0; m.v.y=0;
450   qsort(xs,len, sizeof(int), cmp_int);
451   for (i = cut; i < len - cut; i++){ // all but cutted
452     m.v.x += xs[i];
453   }
454   qsort(ys, len, sizeof(int), cmp_int);
455   for (i = cut; i < len - cut; i++){ // all but cutted
456     m.v.y += ys[i];
457   }
458   vs_free(xs);
459   vs_free(ys);
460   m.v.x/=(len - (2.0 * cut));
461   m.v.y/=(len - (2.0 * cut));
462   return m;
463 }
464 
localmotionsGetMatch(const LocalMotions * localmotions)465 VSArray localmotionsGetMatch(const LocalMotions* localmotions){
466   VSArray m = vs_array_new(vs_vector_size(localmotions));
467   for (int i=0; i<m.len; i++){
468     m.dat[i]=LMGet(localmotions,i)->match;
469   }
470   return m;
471 }
472 
473 
474 /*
475  * Local variables:
476  *   c-file-style: "stroustrup"
477  *   c-file-offsets: ((case-label . *) (statement-case-intro . *))
478  *   indent-tabs-mode: nil
479  *   c-basic-offset: 2 t
480  * End:
481  *
482  * vim: expandtab shiftwidth=2:
483  */
484