Source code for the computation of the benchmark metric: Bandwidth



    // Make new tuning curve at 1 deg resolution
    nn = 360;
    xnew = (float *)myalloc(nn*sizeof(float));
    for(i=0;i<nn;i++)
      xnew[i] = (float)i;

    // Period of Hanning window is 54 deg (Ringach et al)
    ynew = interp_mask_circular_farray(gval,sr,ndg->n,360.0,xnew,nn,
		       "hanning",54.0,0.0);

    //append_farray_xy_plot("zzz.pl",gval,sr,ndg->n,"tuning");
    //append_farray_xy_plot("zzz.pl",xnew,ynew,nn,"interp");

    frac = 1.0 / sqrt(2.0);
    width_at_frac_height_circ_farray(ynew,nn,0,nn-1,frac,
			     &b_max,&b_imax,&b_left,&b_right);



/**************************************-**************************************/
/*                                                                           */
/*                           INTERP_MASK_CIRCULAR_W                          */
/*                                                                           */
/*  Get the weight for the x-position of the given mask.                     */
/*                                                                           */
/*****************************************************************************/
float interp_mask_circular_w(x,x0,period,masktype,maskpar1,maskpar2)
     float x;                  // x-coord for weight  [0..period)
     float x0;                 // mask is centered here  [0..period)
     float period;             // Distance around circle
     char *masktype;           // Interpolation mask [mn]
     float maskpar1,maskpar2;  // Parameters for mask
{
  float dx,w;

  dx = x - x0;

  if (dx > period/2.0)
    dx -= period;
  if (dx <= -period/2.0)
    dx += period;

  if (strcmp(masktype,"hanning")==0){
    //  maskpar1 = wavelength of cosine wave
    if ((dx > -maskpar1/2.0) && (dx < maskpar1/2.0))
      w = 0.5 * (1.0 + cos(2.0*M_PI * dx/maskpar1));
    else
      w = 0.0;  // Weights are zero outside of one period
  }else
    exit_error("INTERP_MASK_CIRCULAR_W","Unknown mask type");

  return w;
}
/**************************************-**************************************/
/*                                                                           */
/*                         INTERP_MASK_CIRCULAR_FARRAY                       */
/*                                                                           */
/*  Return a set of y-values that are interpolated from the x,y data,        */
/*  using a mask of the type specified.                                      */
/*                                                                           */
/*****************************************************************************/
float *interp_mask_circular_farray(xdata,ydata,n,period,xnew,nnew,
			   masktype,maskpar1,maskpar2)
     float *xdata,*ydata;      // Original (x,y) data [n]
     int n;                    // Number of points in original data
     float period;             // Distance around circle
     float *xnew;              // X-coords for interpolated data [nnew]
     int nnew;                 // Number of points in interpolated array
     char *masktype;           // Interpolation mask name, e.g., "hanning"
     float maskpar1,maskpar2;  // Parameters for mask
{
  int i,j;
  float *d,totw,maskw;

  d = (float *)myalloc(nnew*sizeof(float));

  for(i=0;i<nnew;i++){  // For each new point

    totw = 0.0;  // Total weight
    for(j=0;j<n;j++){  // Add up weights of each original point

      // Call subroutine to get mask weight
      maskw = interp_mask_circular_w(xdata[j],xnew[i],period,masktype,
		     maskpar1,maskpar2);
      d[i] += maskw * ydata[j];
      totw += maskw;
    }
    d[i] /= totw;
  }
  return d;
}
/**************************************-**************************************/
/*                                                                           */
/*                       WIDTH_AT_FRAC_HEIGHT_CIRC_FARRAY                    */
/*                                                                           */
/*  Like "width_at_frac_height" but assumes the farray wraps around, so the  */
/*  search for left or right rise points might have to wrap around.          */
/*                                                                           */
/*  Find the maximum value between p1 and p2, inclusive.  Search from the    */
/*  maximum to the left for the first value less than or equal to frac*max.  */
/*  Do the same to the right.                                                */
/*                                                                           */
/*  Values returned:                                                         */
/*    "max" the maximum value.                                               */
/*    "imax" the position of the maximum.                                    */
/*    "left" the position of the left fractional height.                     */
/*    "right" the position of the right fractional height.                   */
/*                                                                           */
/*****************************************************************************/
void width_at_frac_height_circ_farray(data,n,p1,p2,frac,max,imax,left,right)
     float *data;
     int n,p1,p2;
     float frac,*max;
     int *imax,*left,*right;
{
  int i;

  *imax = max_coord_farray(data+p1,p2-p1+1);
  *imax += p1;
  *max = data[*imax];

  *left = -1;
  for(i=*imax;i>=0;i--)
    if ((data[i] <= frac*(*max))&&(*left == -1))
      *left = i;
  if (*left == -1) // Wrap around to upper part of array
    for(i=n-1;i>*imax;i--)
      if ((data[i] <= frac*(*max))&&(*left == -1))
      *left = i;

  *right = -1;
  for(i=*imax;i<n;i++)
    if ((data[i] <= frac*(*max))&&(*right == -1))
      *right = i;
  if (*right == -1) // Wrap around to lower part of array
    for(i=0;i<*imax;i++)
      if ((data[i] <= frac*(*max))&&(*right == -1))
      *right = i;
}
/**************************************-**************************************/
/*                                                                           */
/*                              MAX_COORD_FARRAY                             */
/*                                                                           */
/*  Return the index of the *first* occurrence of the maximum value in the   */
/*  floating point array.                                                    */
/*                                                                           */
/*  WYETH - This could be renamed "get_index_max_farray".                    */
/*                                                                           */
/*****************************************************************************/
int max_coord_farray(data,n)
     float *data;
     int n;
{
  int i,k;
  float max;
  
  if (n < 0) {
    printf("(max_coord_farray):  Array length < 0.\n");
    return(-1);
  }
  k = 0;
  max = data[0];
  for (i=1;i<n;i++)
    if (data[i] > max){
      max = data[i];
      k = i;
    }
  return k;
}