/***/
/* findTilt.M: find and apply "optimal" J-resolved tilt.
/*             Expects variable "mid" is Y-Axis center for tilt.
/***/

if (yLoc != 1)
   {
    exit( 0 );
   };

float matrix[xSize*ySize], proj[xSize];

firstTilt = -0.60;
lastTilt  =  0.60;
tiltStep  =  0.05;

ix1       = 583;
ixn       = 800;
iy1       = 94;
iyn       = 165;

bestTilt = firstTilt;
bestHi   = 0.0;

(void) dReadB( inUnit, matrix, wordLen*xSize*ySize );

/***/
/* For each proposed tilt:
/*  Apply the tilt and build a projection over the region of interest.
/*  Quality of tilt is judged by max value of projection.
/***/

for( pass = 1; pass <= 2; pass++ )
   {
    for( tilt = firstTilt; tilt <= lastTilt; tilt += tiltStep )
       {
        (void) vsSet( proj, xSize, 0.0 );

        for( iy = iy1; iy <= iyn; iy++ )
           {
            shift = tilt*(iy - mid);

            (void) vvCopyOff( rdata, matrix, xSize, 0, iy*xSize );
            (void) vsSet( idata, xSize, 0.0 );

            (void) hilbert( rdata, idata, size );
            (void) ift( rdata, idata, size );
            (void) phase( rdata, idata, size, 0.0, -360.0*shift );
            (void) fft( rdata, idata, size );

            (void) vvAdd( proj, rdata, size );
           }

        maxHi = proj[ix1];
    
        for( ix = ix1; ix <= ixn; ix++ )
           {
            if (proj[ix] > maxHi)
               {
                maxHi = proj[ix];
               };
           }
   
        if (maxHi > bestHi)
           {
            bestHi = maxHi;
            bestTilt = tilt;
           };
       }

    firstTilt = bestTilt - 1.5*tiltStep;
    lastTilt  = bestTilt + 1.5*tiltStep;
    tiltStep /= 10.0;
   }

/***/
/* Apply the best tilt:
/***/

(void) printf( "Best tilt was %f\n", bestTilt );

for( iy = 0; iy < ySize; iy++ )
   {
    shift = tilt*(iy - mid);

    (void) vvCopyOff( rdata, matrix, xSize, 0, iy*xSize );
    (void) vsSet( idata, xSize, 0.0 );

    (void) hilbert( rdata, idata, size );
    (void) ift( rdata, idata, size );
    (void) phase( rdata, idata, size, 0.0, -360.0*shift );
    (void) fft( rdata, idata, size );

    (void) vvCopyOff( matrix, rdata, xSize, iy*xSize, 0 );
   }

(void) dWrite( outUnit, matrix, wordLen*xSize*ySize );
