tt.cpp

it will work it is recogntion code only

 }
  //unsigned char Peanoretina[max2Size*max2Size];
  //unsigned char*Peanoretina=new unsigned char[max2Size*max2Size];
{int x,y,l=0,d0y;
    for(y=0;y<ySize;y++)
        {d0y=dilute0[y];l=xSize*y;
        for(x=0;x<xSize;x++)Peanoretina[dilute[x]+d0y]=retina[x+l];
        }
}
#endif  ///////////////END PEANO 
//N2pi=4*Nfi4; Npi=2*Nfi4;
if(Verbose)tic(7);
int MaxInd_p,Enfi1,Ltr;

#ifdef ROTATING_BOX
if ((xSize != xSize0) || (ySize != ySize0))
   {delete[]Ends;Ends=allEndsRotatingBox (xSize,ySize,Nfi4,dp,dt,
                &MaxInd_p,&Enfi1,&Ltr);}
if(Verbose)Toc(7,"allEndsRotatingBox ");
#else 
int Lrulelines,Llines;//these variables are udes in "alllines", 
                      //We do not use "alllines" as it seems does not speed up
if ((xSize != xSize0) || (ySize != ySize0))
   {delete[]Ends;
        Ends=          allEnds(xSize,ySize,Nfi4,dp,dt,
                &MaxInd_p,&Lrulelines,&Llines,&Enfi1,&Ltr);}
if(Verbose)Toc(7,"allEnds ");
#endif 

////////////////ALLLINES IMITATION
#ifdef  USEPEANO
//PSN("                         USEPEANO: Peano indexing")
#else 
//PSN("                         Standard indexing")           
#endif
{using namespace TRACE;  dimp=2*MaxInd_p+1;dlt=dt;//stagetrace();
if(Verbose){NPV(dimp);PVN(MaxInd_p)}

if (Ltr>ML)
      { delete[]tr;tr=new unsigned char[Ltr]; ML=Ltr;}
stagetraceForAnImage();//it is after ML is renewed -- append room for trace matrices
objectsMemoryConstructor();//it is after ML is renewed

unsigned char*trC;//"tr" stands for track

//void alllines(int *Ends,int Nfi4,int xSize,
//              int *rulelines,int *lines,int Enfi1)
if(Verbose)tic(11);
#ifdef ROTATING_BOX
int  NP=2*MaxInd_p+1;// - number of lines for ANY fi in the Box
L=Ltr; lend= (L/2);lbeg= -lend;
    int E=0,rl=0,l=0,xB,yB,Count1,it,xSizeBN=xSize*BN,
        xB_inc,yB_inc;
    for (nfi=0;nfi<=Nfi4;nfi++)
    {E++;//skip nfi 
     npshift= MaxInd_p-((NP-1)>>1);//npshift here should be always=0
     xB_inc=Ends[E++];yB_inc=Ends[E++];
     for (np=0;np<NP;np++)  
         { 
           xB=Ends[E++];
           yB=Ends[E++];
           trC=tr;  
         for (it=L;it--;) // line is traced:
         {     
#else
   int E=0,*EndsTemp=Ends,rl=0,l=0,NP,xB,yB,Count1,it,xSizeBN=xSize*BN,
        xB_inc,yB_inc;
    for (nfi=0;nfi<=Nfi4;nfi++)
    {E++;Ends++;//skip nfi 
     NP=(*Ends++);// - NP=2*Np+1 - number of lines for given fi
     npshift= MaxInd_p-((NP-1)>>1);
     xB_inc=(*Ends++);yB_inc=(*Ends++);
     for (np=0;np<NP;np++)  
         { lbeg=(*Ends++);
           lend=(*Ends++);  
           xB=(*Ends++);
           yB=(*Ends++);
         L=lend-lbeg+1;trC=tr;
         for (int it=L;it;it--) // line is traced:
         {     
#endif 
              
#ifdef USEPEANO
*trC++=Peanoretina[dilute[xB>>dg]+dilute0[yB>>dg]];
#else 
*trC++=retina[(xB>>dg)+xSize*(yB>>dg)];           
#endif
 
                  xB += xB_inc;  
                  yB += yB_inc;
         }
          computeObjects_giveBackAdj();//tr,lbeg,lend,L,dlt; - should be defined
            writeObjects();
         }//for (np
    }//for (nfi=
#ifdef ROTATING_BOX
E=Enfi1;            
Count1=lend-lbeg;
    int Nfi=2*Nfi4;
    for (nfi=Nfi-1;nfi>Nfi4;nfi--)//nfi is in reverse order
    {E++;//skip nfi
     npshift= MaxInd_p-((NP-1)>>1);
     xB_inc=Ends[E++];yB_inc=-Ends[E++];
     for (np=0;np<NP;np++)  
         { xB=xSizeBN-Ends[E++]-xB_inc*Count1;//take out??
           yB=Ends[E++]-yB_inc*Count1;
           trC=tr;
         for (it=L;it--;) // line is traced:
         { 
#else
E=Enfi1; Ends=EndsTemp+Enfi1;
    int Nfi=2*Nfi4;
    for (nfi=Nfi-1;nfi>Nfi4;nfi--)//nfi is in reverse order
    {E++;Ends++;//skip nfi
     NP=(*Ends++);     // should be NP=2*Np+1
     npshift= MaxInd_p-((NP-1)>>1);
     xB_inc=(*Ends++);yB_inc=-(*Ends++);
     for (np=0;np<NP;np++)  
         { lend=-(*Ends++);
           lbeg=-(*Ends++); 
           Count1=lend-lbeg;L=Count1+1;
           xB=xSizeBN-(*Ends++)-xB_inc*Count1;
           yB=(*Ends++)-yB_inc*Count1;
           trC=tr;
         for (it=L;it;it--) // line is traced:
         {     
#endif
    
#ifdef USEPEANO
*trC++=Peanoretina[dilute[xB>>dg]+dilute0[yB>>dg]];
#else 
*trC++=retina[(xB>>dg)+xSize*(yB>>dg)];           
#endif
                  xB += xB_inc;  
                  yB += yB_inc;
         }
        computeObjects_giveBackAdj();//tr,lbeg,lend,L,dlt; - should be defined
        writeObjects();
        }//for (np
    }//for (nfi=
#ifndef ROTATING_BOX
  Ends=EndsTemp;
#endif
////////////////END ALLLINES IMITATION
if(Verbose)Toc(11,"Computations of all trace transforms ");// computeObjects_giveBackAdj();

int nf; //write (ie. write to disk) trace matrice 
if (wTpgm || wTtxt || wT) for (nf=0;nf<LlistU;nf++)
 { int tfnl=listU[nf]; char buf[4];sprintf(buf,"%03i",tfnl);
   char*trMatrfile=concat(concat(DirTrMatrs,nameNE_T),buf);
   if (wTtxt) // write trace matrices in *.txt, for matlab 
    writeMatrixDouble(concat(trMatrfile,".txt"),dimp,N2pi,trMatr[nf]);
   if (wT) // write trace matrices in small format *.sm 
     writeMatrixSmall(concat(trMatrfile, ".sm"),dimp,N2pi,trMatr[nf]);
   if (wTpgm) // write trace matrices in small format *.sm 
    { int LTM=dimp*N2pi,i; unsigned char uc;
      if(LMpgm<LTM){delete[]Mpgm;Mpgm=new unsigned char[LTM];LMpgm=LTM;}
      double *TM=trMatr[nf],t,delta,min=TM[0],max=TM[0];
      for (i=1;i<LTM;i++)
       {t=TM[i];if(t<min)min=t;if(t>max)max=t;}
       delta=max-min; if (delta EQ 0) delta=1;
      for (i=0;i<LTM;i++)
       { uc=round((TM[i]-min)/delta*255);Mpgm[i]=uc;}
      pgmP5Write(concat(trMatrfile, ".pgm"),
                  dimp,N2pi,Mpgm);
      
    }     

 }

}//END:  using namespace TRACE;  dimp=2*MaxInd_p+1;dlt=dt;stagetrace();
rSize0=rSize;xSize0=xSize;ySize0=ySize;
if(Verbose){PN Timer1result("stage TRACE took including reading the image and writing");}
//======== trace matrices are computed for image imNo.
if (LdiList EQ 0) continue; //that is go to consider another image = another imNo

{using namespace DIAM;
if(Verbose)tic(321); 
L=dimp;dlt=dp;lbeg=-MaxInd_p;lend=MaxInd_p;
////the folowing two lines should be done for each image
ML=dimp;//trace length = length of trace matrices along p-coordinate 
objectsMemoryConstructor();


fillCircuses();//before calling this, fill DIAM::dlt=dp;
if(Verbose){PN Toc(321,"plain circuses are filled:");} //exit(0);
//======== Now plain circuses are ready, they are in double**circ 

//______________ writing down (plain) circuses:
double **cW=circ;// - circuses to write down. Change 6 variables in if() only and extensions
if(iwCtxt)  writeAllCircDouble(concat(concat(DirCircuses,nameNE),"C.txt"),
                               LtrList,LdiList,N2pi,cW);
if(iwC)    writeAllCircByParts(concat(concat(DirCircuses,nameNE),"C.sc"),
                               LtrList,LdiList,N2pi,cW);
// circus for trace functional No.No.nf
//        circus functional No.No. cf starts from
//        cW[nf]+cf*N2pi and takes N2pi double elements:
int nf;
if (wCtxt || wC) for (nf=0;nf<TRACE::LlistU;nf++)//writes text for matlab
 { int tfnl=TRACE::listU[nf]; char buf[4];sprintf(buf,"%03i",tfnl);
   char*trMatrfileBare=//"Bare" means it is without extention
                concat(concat(concat(DirCircuses,nameNE_T),buf),"C");
   for (int cf=0;cf<LlistU;cf++)
   { int cfnl=listU[cf]; sprintf(buf,"%03i",cfnl);
     char*circfile=concat(trMatrfileBare,buf);
    if (wCtxt) writeMatrixDouble(concat(circfile,".txt"),N2pi,1,cW[nf]+cf*N2pi);
    if (wC)    writeMatrixSmall (concat(circfile,".sm"), N2pi,1,cW[nf]+cf*N2pi);
   }
 }
//END____________ writing down (plain) circuses.
//_______________ Computing normalisors:
double dfi=2*pi/N2pi;   
if (computeNormalisors) 
for (nf=0;nf<Lt;nf++)
 { int i,tfnl=TRACE::listU[nf]; 
   for (int cf=0;cf<LlistU;cf++)//in this namespace "LlistU" is "Ld";
   { int cfnl=listU[cf];
      if(circAuthorNormalisor[Ld*nf+cf])
       {double *ci=circ[nf]+cf*N2pi;//- pointer to a given circus NoNo=(nf,cf) 
        double *cA=circA[nf]+cf*N2pi;//- pointer to a given assoc.circus NoNo=(nf,cf) 
        //PA(ci,0,N2pi-1);PN; ???
        double uomega = circUomega[Ld*nf+cf];// uomega=lambdaP*kappaT-kappaP; 
     // compute associated circus circA:
        if (uomega)
          for(i=0;i<N2pi;i++)
                 {double absc=fabs(ci[i]);double sign=ci[i]>0?1.0:-1.0; 
                  cA[i]=0;
                  if(absc)cA[i]= sign*exp(uomega*log(absc));//associated circus, simplest case                   
                 }
          else//ie.  uomega EQ 0
            { 
              cA[0]=ci[1]-ci[N2pi-1];cA[N2pi-1]=ci[0]-ci[N2pi-2];
              for(i=1;i<N2pi-1;i++) cA[i]=ci[i+1]-ci[i-1];
              double Co=2*dfi;  
              for(i=0;i<N2pi;i++) //associated circus for uomega EQ 0
                     cA[i]=sqrt(Co*fabs(cA[i]))*(cA[i]>=0?1.0:-1.0);
            }
     // affinely normalising parameters of the associated circus circA:
       double a0=0,a2=0,b2=0,ciA4,fi;
       for(i=0;i<N2pi;i++)//find corresponding normalisor
       {ciA4=cA[i]*cA[i];ciA4=ciA4*ciA4;//cA^4
        fi=i*dfi; 
        a0+=ciA4; a2+=cos(2*fi)*ciA4;b2+=sin(2*fi)*ciA4;        
       }   
     a0=a0/8;a2=a2/8;b2=b2/8; 
     double
     maxAngle = atan2(b2,a2)/2, minAngle=maxAngle+pi/2,
     maxValue = a0+a2*cos(2*maxAngle)+b2*sin(2*maxAngle),//=sigma1
     minValue = a0+a2*cos(2*minAngle)+b2*sin(2*minAngle),//=sigma2
     beta=0;
     if ( maxValue>0 && minValue>0 )beta = sqrt(sqrt(minValue/maxValue));
       //0<=beta<=1 // normalisor is ready
                               /// normalisor is a pair (maxAngle,beta);
                               /// normalisor is found 
     NormMaxAngle[Ld*nf+cf]=maxAngle; NormBeta[Ld*nf+cf]=beta;                               
       }// if(circAuthorNormalisor[Ld*nf+cf])
       else// case when the circus is not the author: no normalisor:
     { double *cA=circA[nf]+cf*N2pi;//- pointer to a given assoc.circus NoNo=(nf,cf) 
       for(i=0;i<N2pi;i++)cA[i]=0;
       NormMaxAngle[Ld*nf+cf]=1; NormBeta[Ld*nf+cf]=0;}//
}}// if (computeNormalisors) Result: if normalor does not exist, beta=0
if (Verbose && computeNormalisors) {PN PA(NormMaxAngle,0,Lc-1);  PA(NormBeta,0,Lc-1);}//???DEBUG
//END____________ Computing normalisors

//_______________ Computing normalised circuses:
if(computeNormalisedCircusesP || computeNormalisedCircusesA)
for (nf=0;nf<Lt;nf++)
 for (int cf=0;cf<LlistU;cf++)//in this namespace "LlistU" is "Ld";
   { double uomegaNative=circUomega[Ld*nf+cf]; double maxAngle,beta;
     double *c=circ[nf]+cf*N2pi,*cA=circA[nf]+cf*N2pi,
            *ciPN=circPN[nf]+cf*N2pi,*ciAN=circAN[nf]+cf*N2pi;
 
      //choose normalisor:  ?????? 
      // if (circAuthorNormalisor[Ld*nf+cf])
         {maxAngle=NormMaxAngle[Ld*nf+cf],beta=NormBeta[Ld*nf+cf];}
      //????????????????????? this works only for invariant circuses       
           
      if(computeNormalisedCircusesP)         
       normalising(maxAngle,beta,N2pi,c,uomegaNative, //normalisation of plain circus
                    ciPN);
      if(computeNormalisedCircusesA)         
       normalising(maxAngle,beta,N2pi,cA,-1.0,//normalisation of assoc. circus
                    ciAN);
   }//for (int cf=0;cf<LlistU;cf++)
//END____________ Computing normalised circuses

//______________ writing down plain normalised circuses:
cW=circPN;// - circuses to write down. Change 6 variables in if() only and extensions
if(iwCPtxt)  writeAllCircDouble(concat(concat(DirCircuses,nameNE),"CP.txt"),
                               LtrList,LdiList,N2pi,cW);
if(iwCP)    writeAllCircByParts(concat(concat(DirCircuses,nameNE),"CP.sc"),
                               LtrList,LdiList,N2pi,cW);

if (wCPtxt || wCP) for (nf=0;nf<TRACE::LlistU;nf++)//writes text for matlab
 { int tfnl=TRACE::listU[nf]; char buf[4];sprintf(buf,"%03i",tfnl);
   char*trMatrfileBare=//"Bare" means it is without extention
                concat(concat(concat(DirCircuses,nameNE_T),buf),"CP");
   for (int cf=0;cf<LlistU;cf++)
   { int cfnl=listU[cf]; sprintf(buf,"%03i",cfnl);
     char*circfile=concat(trMatrfileBare,buf);
    if (wCPtxt) writeMatrixDouble(concat(circfile,".txt"),N2pi,1,cW[nf]+cf*N2pi);
    if (wCP)    writeMatrixSmall (concat(circfile,".sm"), N2pi,1,cW[nf]+cf*N2pi);
   }
 }
//END____________ writing down plain normalised circuses.

//______________ writing down assoc. normalised circuses:
cW=circAN;// - circuses to write down. Change 6 variables in if() only and extensions
if(iwCAtxt)  writeAllCircDouble(concat(concat(DirCircuses,nameNE),"CA.txt"),
                               LtrList,LdiList,N2pi,cW);
if(iwCA)    writeAllCircByParts(concat(concat(DirCircuses,nameNE),"CA.sc"),
                               LtrList,LdiList,N2pi,cW);

if (wCAtxt || wCA) for (nf=0;nf<TRACE::LlistU;nf++)//writes text for matlab
 { int tfnl=TRACE::listU[nf]; char buf[4];sprintf(buf,"%03i",tfnl);
   char*trMatrfileBare=//"Bare" means it is without extention
                concat(concat(concat(DirCircuses,nameNE_T),buf),"CA");
   for (int cf=0;cf<LlistU;cf++)
   { int cfnl=listU[cf]; sprintf(buf,"%03i",cfnl);
     char*circfile=concat(trMatrfileBare,buf);
    if (wCAtxt) writeMatrixDouble(concat(circfile,".txt"),N2pi,1,cW[nf]+cf*N2pi);
    if (wCA)    writeMatrixSmall (concat(circfile,".sm"), N2pi,1,cW[nf]+cf*N2pi);
   }
 }
//END____________ writing down assoc. normalised circuses.

if(Verbose){PN Toc(321,"overall stage DIAM took");}
}//using namespace DIAM;

if (Stages>2)
{using namespace TRIPLE;
if(Verbose)tic(333);
L=N2pi;dlt=2*pi/N2pi;lbeg=0;lend=L-1;

if (iwF+iwFtxt)
  {
   fillTriple(circ,triple); 
   if(iwFtxt)writeMatrixDouble(concat(concat(DirTriple,nameNE),"F.txt"),Ltf,1,triple);
   //if(iwF) ???
  }
if (iwFP+iwFPtxt)
  {
   fillTriple(circPN,triplePN);
   if(iwFPtxt)writeMatrixDouble(concat(concat(DirTriple,nameNE),"FP.txt"),Ltf,1,triplePN);
   //if(iwFP) ???
  }
if (iwFA+iwFAtxt)
  {
    fillTriple(circAN,tripleAN);
   if(iwFAtxt)writeMatrixDouble(concat(concat(DirTriple,nameNE),"FA.txt"),Ltf,1,tripleAN);
   //if(iwFA) ???
  }
if(Verbose)Toc(333,"overall stage TRIPLE took");
}//{using namespace TRIPLE;

if(Verbose){PN Timer2result("working with this image");}



}}//{for (imNo=0;imNo<LimList;imNo++){  //THE BIG IMAGE LOOP

return (0);
}//int main(){