Cross Correlation

6b2dfcef · Pilar Cossio · 8510f0b5 · 6b2dfcef
Commit 6b2dfcef authored Jul 29, 2014 by Pilar Cossio
Hide whitespace changes
Inline Side-by-side

Showing with 61 additions and 5 deletions

bioem.cpp bioem.cpp +61 -5

No files found.
--- a/bioem.cpp
+++ b/bioem.cpp
@@ -362,7 +362,6 @@ int bioem::run()
                         {
                                 for (int cent_y = 0; cent_y < param.param_device.NumberPixels - param.param_device.CCdisplace ; cent_y = cent_y + param.param_device.CCdisplace)
                                 {
-        	                       cout << cc << " " << cent_x << " " << cent_y <<"\n"; 
                                 bioem_Probability_cc& pProbCC = pProb.getProbCC(iRefMap, cc);
                                 pProbCC.forCC = 0.0;
 
@@ -619,16 +618,58 @@ int bioem::run()
 			}
 			if(param.param_device.writeCC)
 			{
-				int  cc=0;
+			int  cc=0;
+                        int halfPix,rx,ry;
+	                myfloat_t localcc[ (param.param_device.NumberPixels+1) * (param.param_device.NumberPixels+1) ];
+
+                        halfPix = param.param_device.NumberPixels / 2 ;
+			// Ordering the centers of the Cross Correlation
+
+		for (int rx = param.param_device.CCdisplace; rx < param.param_device.NumberPixels -param.param_device.CCdisplace ; rx = rx + param.param_device.CCdisplace)
+                         {
+                                 for (int ry = param.param_device.CCdisplace; ry < param.param_device.NumberPixels - param.param_device.CCdisplace ; ry = ry + param.param_device.CCdisplace)
+                                 {
+				localcc[ rx * param.param_device.NumberPixels + ry ] = 0.0;
+                                }
+                        }
+
 			for (int cent_x = 0; cent_x < param.param_device.NumberPixels -param.param_device.CCdisplace ; cent_x = cent_x + param.param_device.CCdisplace)
                         {
                                 for (int cent_y = 0; cent_y < param.param_device.NumberPixels - param.param_device.CCdisplace ; cent_y = cent_y + param.param_device.CCdisplace)
                                 {
+				//localcc[ rx * param.param_device.NumberPixels + ry ] = 0.0;
 						bioem_Probability_cc& pProbCC = pProb.getProbCC(iRefMap, cc);
-						ccProbfile << " " << iRefMap << " " << cent_x << " " << cent_y << " " << pProbCC.forCC <<"\n";
-				}
-
+						if(cent_x < halfPix && cent_y < halfPix){
+					//	ccProbfile << " " << iRefMap << " " << (myfloat_t) halfPix  - cent_x << " " << halfPix - cent_y << " " << pProbCC.forCC <<"\n";
+						rx = halfPix  - cent_x;
+                                                ry = halfPix  - cent_y;}
+						if(cent_x >= halfPix && cent_y < halfPix){
+                                          //      ccProbfile << " " << iRefMap << " " << (myfloat_t) 3 * halfPix  - cent_x << " " << halfPix - cent_y << " " << pProbCC.forCC <<"\n"; 
+						rx = 3 * halfPix  - cent_x;
+                                                ry = halfPix  - cent_y;}
+                                                if(cent_x < halfPix && cent_y >= halfPix){
+                                          //      ccProbfile << " " << iRefMap << " " << (myfloat_t) halfPix  - cent_x << " " << 3 * halfPix - cent_y << " " << pProbCC.forCC <<"\n";
+						rx = halfPix  - cent_x;
+                                                ry = 3 * halfPix  - cent_y;}
+                                                 if(cent_x >= halfPix && cent_y >= halfPix){
+                                        //        ccProbfile << " " << iRefMap << " " << 3* halfPix  - cent_x << " " << 3 * halfPix - cent_y << " " << pProbCC.forCC <<"\n";
+						rx = 3 * halfPix  - cent_x;
+                                                ry = 3 * halfPix  - cent_y;}
+//						cout << " TT " << cent_x << " " << rx << " " << cent_y << " " << ry << " " <<  pProbCC.forCC << "\n";
+						localcc[ rx * param.param_device.NumberPixels + ry ] = pProbCC.forCC;
+                                  cc++;
+	 			}
+                                  //              ccProbfile << "\n";
 			}
+		         for (int rx = param.param_device.CCdisplace; rx < param.param_device.NumberPixels -param.param_device.CCdisplace ; rx = rx + param.param_device.CCdisplace)
+                         {
+                                 for (int ry = param.param_device.CCdisplace; ry < param.param_device.NumberPixels - param.param_device.CCdisplace ; ry = ry + param.param_device.CCdisplace)
+                                 {
+                                    ccProbfile << rx << " " << ry << " " << localcc[ rx * param.param_device.NumberPixels + ry ] << "\n" ; 
+				}
+				ccProbfile << "\n";
+			}			
+  
                       }
 		}

@@ -689,6 +730,21 @@ inline void bioem::calculateCCFFT(int iRefMap, int iOrient, int iConv, myfloat_t
 	myfftw_execute_dft_c2r(param.fft_plan_c2r_backward, localCCT, lCC);

 	doRefMapFFT(iRefMap, iOrient, iConv, lCC, sumC, sumsquareC, pProb, param.param_device, RefMap);
+
+#ifdef PILAR_DEBUG
+         if (param.param_device.writeCC)
+                {      int  cc=0;
+                         for (int cent_x = 0; cent_x < param.param_device.NumberPixels -param.param_device.CCdisplace ; cent_x = cent_x + param.param_device.CCdisplace)
+                         {
+                                 for (int cent_y = 0; cent_y < param.param_device.NumberPixels - param.param_device.CCdisplace ; cent_y = cent_y + param.param_device.CCdisplace)
+                                 {
+                                       cout << "CHECKCC " << " " << cent_x << " " << cent_y <<" " << lCC[cent_x * param.param_device.NumberPixels + cent_y] / (myfloat_t) (param.param_device.NumberPixels * param.param_device.NumberPixels ) << "\n";
+                                 cc++;
+                                 }
+                         }
+                }
+#endif
+
 }

 int bioem::createProjection(int iMap, mycomplex_t* mapFFT)