(*ker).k[0] = 1.0 ;    (*ker).f = 1.0 ;  }  (*ker).m = h ; }/*            NAME : init_beaudet(Ix,Ixx,Ixy)   PARAMETER(S) :  Ix : operator for 1st order x derivative;                  Ixx : operator for 2nd order x derivative;                  Ixy : operator for partial x and y derivatives.         PURPOSE : Initialize the image operators defined by Beaudet                  to obtain numerical derivatives.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : September 13 1990*/init_beaudet(Ix,Ixx,Ixy)beaudet_t *Ix, *Ixx, *Ixy ;{  (*Ix).g[0][0] = -2 ; (*Ix).g[1][0] = -1 ; (*Ix).g[2][0] = 0 ;   (*Ix).g[3][0] = 1 ;  (*Ix).g[4][0] = 2 ;  (*Ix).g[0][1] = -2 ; (*Ix).g[1][1] = -1 ; (*Ix).g[2][1] = 0 ;   (*Ix).g[3][1] = 1 ;  (*Ix).g[4][1] = 2 ;  (*Ix).g[0][2] = -2 ; (*Ix).g[1][2] = -1 ; (*Ix).g[2][2] = 0 ;   (*Ix).g[3][2] = 1 ;  (*Ix).g[4][2] = 2 ;  (*Ix).g[0][3] = -2 ; (*Ix).g[1][3] = -1 ; (*Ix).g[2][3] = 0 ;   (*Ix).g[3][3] = 1 ;  (*Ix).g[4][3] = 2 ;  (*Ix).g[0][4] = -2 ; (*Ix).g[1][4] = -1 ; (*Ix).g[2][4] = 0 ;   (*Ix).g[3][4] = 1 ;  (*Ix).g[4][4] = 2 ;  (*Ix).f = 50.0 ;  (*Ix).m = 5 ;        (*Ixx).g[0][0] = 2 ;  (*Ixx).g[1][0] = -1 ;  (*Ixx).g[2][0] = -2 ;   (*Ixx).g[3][0] = -1 ;  (*Ixx).g[4][0] = 2 ;  (*Ixx).g[0][1] = 2 ; (*Ixx).g[1][1] = -1 ; (*Ixx).g[2][1] = -2;   (*Ixx).g[3][1] = -1 ; (*Ixx).g[4][1] = 2 ;  (*Ixx).g[0][2] = 2 ; (*Ixx).g[1][2] = -1 ; (*Ixx).g[2][2] = -2 ;   (*Ixx).g[3][2] = -1 ; (*Ixx).g[4][2] = 2 ;  (*Ixx).g[0][3] = 2 ; (*Ixx).g[1][3] = -1 ; (*Ixx).g[2][3] = -2 ;   (*Ixx).g[3][3] = -1 ; (*Ixx).g[4][3] = 2 ;  (*Ixx).g[0][4] = 2 ;  (*Ixx).g[1][4] = -1 ;  (*Ixx).g[2][4] = -2 ;   (*Ixx).g[3][4] = -1 ;  (*Ixx).g[4][4] = 2 ;  (*Ixx).f = 35.0 ;  (*Ixx).m = 5 ;  (*Ixy).g[0][0] = 4 ;  (*Ixy).g[1][0] = 2 ;  (*Ixy).g[2][0] = 0 ;   (*Ixy).g[3][0] = -2 ; (*Ixy).g[4][0] = -4 ;  (*Ixy).g[0][1] = 2 ;  (*Ixy).g[1][1] = 1 ;  (*Ixy).g[2][1] = 0 ;   (*Ixy).g[3][1] = -1 ; (*Ixy).g[4][1] = -2 ;  (*Ixy).g[0][2] = 0 ;  (*Ixy).g[1][2] = 0 ;  (*Ixy).g[2][2] = 0 ;   (*Ixy).g[3][2] = 0 ;  (*Ixy).g[4][2] = 0 ;  (*Ixy).g[0][3] = -2 ; (*Ixy).g[1][3] = -1 ; (*Ixy).g[2][3] = 0 ;   (*Ixy).g[3][3] = 1 ;  (*Ixy).g[4][3] = 2 ;  (*Ixy).g[0][4] = -4 ; (*Ixy).g[1][4] = -2 ; (*Ixy).g[2][4] = 0 ;   (*Ixy).g[3][4] = 2 ;  (*Ixy).g[4][4] = 4 ;  (*Ixy).f = 100.0 ;  (*Ixy).m = 5 ;}/*            NAME : init_central(C)   PARAMETER(S) : C : 1-D mask contraining 4-point central difference                      factors.         PURPOSE : Inits a 4-point central difference mask.          AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : April 25 1992*/init_central(C)kernel_t *C ;{  (*C).k[0] = -1.0 ; (*C).k[1] = 8.0 ;   (*C).k[2] = 0.0 ; (*C).k[3] = -8.0 ; (*C).k[4] = 1.0 ;  (*C).m = 5 ;  (*C).f = 12.0 ;}/*            NAME : init_cube(p,n) ;   PARAMETER(S) : p : pointer on an image cube node;                  n : number of image frames to allocate.         PURPOSE : allocates n images to p->gauss_ptr[i].         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : September 15 1990*/init_cube(p,n)qnode_ptr_t p ;int n ;{ int i ;  for (i = 0 ; i < n ; i++) {     p->gauss_ptr[i] = (image512_t *)malloc(sizeof(image512_t)) ;  }}/*            NAME : init_flow(p,n) ;   PARAMETER(S) : p : pointer on an parameter cube node;                  n : number of parameter frames to allocate.         PURPOSE : allocates n parameter frames to p->param_ptr[i]                  and a flow field to p->flow_ptr.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : September 15 1990*/init_flow(p,n)qnode_ptr_t p ;int n ;{ int i ;  p->flow_ptr = (disp_field512_t *)malloc(sizeof(disp_field512_t)) ;  for (i = 0 ; i < n ; i++) {     p->param_ptr[i] = (param512_t *)malloc(sizeof(param512_t)) ;  }}/*            NAME : init_list(h,q)    PARAMETER(S) : h : head list pointer;                  q : tail list pointer.         PURPOSE : Initialization of the head                  and tail pointers to NULL.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : March 17 1990*/init_list(h,q)qnode_ptr_t *h, *q ;{   *h = (qnode_ptr_t)NULL ;  *q = (qnode_ptr_t)NULL ;}/*            NAME : insert_node(p,h,q)    PARAMETER(S) : p : pointer on the node to insert;                  h : head pointer of the list;                  q : tail pointer of the list.         PURPOSE : Insertion of the node p in a doubly                   linked list in order of the level                  numbers (p->level).         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : March 17 1990*/insert_node(p,h,q) qnode_ptr_t p, *h, *q ;{ qnode_ptr_t s, t ;  t = *h ;  s = (qnode_ptr_t)NULL ;  while ((t != (qnode_ptr_t)NULL) && (p->level < t->level)) {    s = t ;    t = t->forth ;  }  if (s != (qnode_ptr_t)NULL) {    s->forth = p ;  }  else {    *h = p ;  }  if (t != (qnode_ptr_t)NULL) {    t->back = p ;  }  else {    *q = p ;  }  p->back = s ;  p->forth = t ;}/*            NAME : itoa(n,s)   PARAMETER(S) : n : integer value ;                  s : output string.         PURPOSE : Converts integer into string.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : June 25 1990*/itoa(n,s)int n ;string s ; { int i, j, digit, neg ;  string t ;  if (n < 0) {     neg = TRUE ;  }  else {    neg = FALSE ;  }  if (n == 0) {    s[0] = '0' ;    s[1] = '\0' ;  }  else {    i = 0 ;    while (abs(n) > 0) {      digit = abs(n) % 10 ;      t[i] = (char)((int)'0' + digit) ;      n = n / 10 ;      i++ ;    }    if (neg) {      t[i] = '-' ;       i++ ;    }    for (j = 0 ; j < i ; j++) {      s[j] = t[i - j - 1] ;    }    s[i] = '\0' ;  }}/*            NAME : mag(u)   PARAMETER(S) : u : 2D flow vector.         PURPOSE : returns the magnitude of u.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : October 1 1990*/float mag(u)disp_vect_t u ;{   return(sqrt(pow(u.x,2.0) + pow(u.y,2.0))) ;}/*            NAME : min_cond(s,f,n,m,trsh,sf)   PARAMETER(S) : s : 1D array containing at most n image flow locations;                  f : image parameter node;                  n : maximum number of image locations in s;                  m : actual number of image locations in s;               trsh : threshold for condition number;                 sf : thresholding boolean value.         PURPOSE : determines among the locations in s, which one has a                  minimal condition number.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : April 4 1992*/pos_t min_cond(s,f,n,m,trsh,sf)pos_t s[X] ;qnode_ptr_t f ;int n, m ;float trsh ;int sf ;{ pos_t loc, err_loc ;  float min_c, gauss ;  int l ;  min_c = MAX_COND + 1 ;  loc.i = -1 ;  loc.j = -1 ;  err_loc.i = -1 ;  err_loc.j = -1 ;  if (m < n) {    n = m ;  }  for (l = 0 ; l < n ; l++) {    if ((*f->param_ptr[f->sizz/2])[f->res*s[l].i + s[l].j].cond < min_c) {      min_c = (*f->param_ptr[f->sizz/2])[f->res*s[l].i + s[l].j].cond ;      gauss = (*f->param_ptr[f->sizz/2])[f->res*s[l].i + s[l].j].gauss ;      loc = s[l] ;    }  }  if (sf) {    if (gauss > trsh) {      return(loc) ;    }    else {      return(err_loc) ;    }  }  else {    return(loc) ;  }}/*            NAME : sort(discr_val,sample,num,n)   PARAMETER(S) : discr_val : discriminant values for image positions                              contained in sample;                  sample    : image positions to be sorted with respect                              to their discriminant values;                  num       : number of elements in sample;                  n         : image region size.         PURPOSE : sorts the n first elements of sample[1..m] with respect                  to the discriminant value.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : September 15 1990*/sort(discr_val,sample,num,n) float discr_val[X] ;pos_t sample[X] ;int num, n ;{ pos_t ts ;   float td ;  int i, j, index ;  if (num < n) {    n = num ;  }  for (i = 0 ; i < n - 1 ; i++) {    index = i ;    for (j = index + 1 ; j < num ; j++) {      if (discr_val[index] > discr_val[j]) {        index = j ;      }    }    td = discr_val[index] ;    discr_val[index] = discr_val[i] ;    discr_val[i] = td ;       ts = sample[index] ;    sample[index] = sample[i] ;    sample[i] = ts ;  }}/*            NAME : min_discr(f,i,j,n,sample,num)   PARAMETER(S) :   f      : image parameter node;                  i,j      : upper left corner location of a n*n image area;                    n      : image area size;                    sample : valid image locations of image area;                    num    : number of image locations in sample.         PURPOSE : extracts the n image locations from the n*n image area                  for which the discriminant value is minimal.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : November 7 1990*/min_discr(f,i,j,n,sample,num)qnode_ptr_t f ;int i, j, n ;pos_t sample[X] ;int *num ;{ float discr_val[X] ;  int k, l ;  *num = 0 ;   for (k = 0 ; k < n ; k++) {    for (l = 0 ; l < n ; l++) {      if ((*f->param_ptr[f->sizz/2])[f->res*(i+k) + j+l].err == NO_ERROR) {        sample[*num].i = i + k ;        sample[*num].j = j + l ;        discr_val[*num] = (*f->param_ptr[f->sizz/2])[f->res*(i+k) + j+l].discr ;        (*num)++ ;      }    }  }  sort(discr_val,sample,*num,n) ;}/*            NAME : norm(V,n)   PARAMETER(S) : V : vector;                  n : length of V.         PURPOSE : returns L2 norm of V         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : November 7 1990*/float norm(V,n) float V[] ;int n ;{ int i ;  float sum ;   sum = 0.0 ;  for (i = 0 ; i < n ; i++) {    sum += (V[i]*V[i]) ;  }  return(sqrt(sum)) ;}/*            NAME : pgetrast(fn,hd,bf,sx,sy,r)    PARAMETER(S) : fn    : filename;                  hd    : image header (raster file);                  bf    : Pointer on a 2D array containing the image;                  sx,sy : image size to read;                   r    : resolution of the array.        PURPOSE : Reads rasterfile specified                  by fn into bf.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : October 27 1989, modified March 20 1990.*/pgetrast(fn,hd,bf,sx,sy,r) char *fn ;unsigned char hd[H] ;float bf[X*Y] ;int sx, sy, r ;{ int fd, i, j ;  unsigned char c ;  if ((fd = open(fn,O_RDONLY)) > 0) { ;    if (read(fd,hd,H) == H) {       for (i = 0 ; i < sy ; i++) {        for (j = 0 ; j < sx ; j++) {          if (read(fd,&c,sizeof(c)) != 1) {            error(4) ;          }          bf[r*i + j] = (float)c ;        }      }    }    else {      error(5) ;    }  }  else {    error(6) ;  }  close(fd) ;}/*            NAME : Bpgetrast(fn,bf,sx,sy,r)    PARAMETER(S) : fn    : filename;                  bf    : Pointer on a 2D array containing the image;                  sx,sy : image size to read;                   r    : resolution of the array.        PURPOSE : Reads rasterfile without header specified                  by fn into bf.         AUTHOR : Steven Beauchemin             AT : University of Western Ontario           DATE : October 27 1989, modified March 20 1990.*/Bpgetrast(fn,bf,sx,sy,r) char *fn ;float bf[X*Y] ;int sx, sy, r ;{ int fd, i, j ;  unsigned char c ;  if ((fd = open(fn,O_RDONLY)) > 0) { ;    for (i = 0 ; i < sy ; i++) {      for (j = 0 ; j < sx ; j++) {        if (read(fd,&c,sizeof(c)) != 1) {          error(4) ;        }        bf[r*i + j] = (float)c ;      }    }  }