package signal.lib;

import java.util.Vector;

/**
 Finds edge locations in functions.  Searches the entire interval on which the function is given.
 
  @see "Detection of Edges in Spectral Data II. Nonlinear Enhancement. UCLA CAM Report No. 99-23."
*/

 public final class edgeDetect {
    private static chebyshevPolynomial CM;

//       IN
//  u          the function to be postprocessed, the function may be known on the 
//             CGL grid on [-1,1] or on a coord. trans. of this grid
//  JCRIT      critical threshold
//  NE         if an edge is located at x[j], then only one edge may be detected 
//             the interval (x[j]-x[j-NE],x[j]+x[j+NE])
//  Q          enhancement exponent

//       OUT
//  d      vector containing the location of the edges in order -1,..,1
//         with d[0]=-1.0 and d[M]=1.0; thus, there will be (M-1) edges found
//  uE     edge series
//  uNLE   enhanced edge series

/**
 Given a function u, critical threshold JCRIT, and enhancement exponent Q, and a neighborhood parameter NE,
 the method returns the edge series uE, the nonlinear enhancement uNLE, a Vector d containing the edge locations.
*/


 public static void 
  findEdges(double[] u, double JCRIT, int Q, double[] uE, double[] uNLE, Vector d, int NE) {
      int        N = u.length - 1;
      double    Nd = (double)N;
      double  aN[] =new double[N+1];              // Chebyshev coefficients
      double   x[] =new double[N+1];

       for ( int j=0; j<=N; j++) x[j]= Math.cos(Math.PI*j/Nd);   // order 1,..,-1

       
// --- find Chebyshev coefficients ------------------------------

     double dx=(1.0/N), dxP=dx*Math.PI;
     double sum2=0;

     for (int k=0; k<=N; k++ ) {
       sum2=0;
       for (int n=0; n<=N; n++ )   sum2 += CM.CJ(n,N)*u[n]*Math.cos(k*n*dxP);
       aN[k]=2.0*CM.CJ(k,N)*dx*sum2;
     } // k 

// --------------------------------------------------------------
       

       for (int j=1; j<N; j++) {
         for (int k=1; k<=N; k++)   uE[j] +=  aN[k]*CM.Tp(x[j],k);
         uE[j] *= (Math.sqrt(2.0)*Math.PI*Math.sqrt(1.0-x[j]*x[j])/Nd);
       } // j

       uE[0]=0;
       uE[N]=0;
                                                                      // nonlinear enhancement
       double temp=0;
       for (int j=0; j<=N; j++) {
        temp = Math.abs( Math.pow(N,0.5*Q)*Math.pow(uE[j],Q) );
         if (temp<JCRIT) uNLE[j] = 0;                             // the direction of the jump is not important
                                                                  //  so take abs val
         else            uNLE[j] = Math.abs(uE[j]);
        } // j
        
     Vector tempV = new Vector();
     
//  ---  limit one edged per neighborhood of size 2*NE ---

     boolean addEdgeL=false, addEdgeR=false;
     tempV.addElement(new Double(1.0));                                                 // find edges
     for (int j=1; j<N; j++) {
     
       if ( (j-NE)>=0 && (j+NE)<=N ) {
       
         if (uNLE[j]>0) {                                                            // check surronding points
           for (int i=1; i<=NE; i++) {
              if (uNLE[j-i]>0 && (uNLE[j-i]<uNLE[j])) uNLE[j-i]=0;
              if (uNLE[j+i]>0 && (uNLE[j+i]<uNLE[j])) uNLE[j+i]=0;
           }
         }
       }  
        else {  // it is within NE grid points of a boundary
          uNLE[j]=0;
       }  // else
 
       } // j
       
        for (int j=1; j<N; j++) if (uNLE[j]>0) tempV.addElement(new Double(x[j]));
     
     tempV.addElement(new Double(-1.0));
     
// ----------------------
     
     // put edge location in correct order, -1,...,1
     int sizeTempV = tempV.size()-1;
     for (int j=0; j<=sizeTempV; j++) d.addElement(tempV.elementAt(sizeTempV-j));


 } // findEdges

} // edgeDetect