Commit e7310dea authored by David Rohr's avatar David Rohr

add printModel flag, get ride of unnecessary fourier transform for FFTALGO

parent 88ded6c5
......@@ -283,6 +283,13 @@ int bioem::run()
// **************************************************************************************
// **** Main BioEM routine, projects, convolutes and compares with Map using OpenMP ****
// **************************************************************************************
if (param.printModel)
{
//....
return(0);
}
// **** If we want to control the number of threads -> omp_set_num_threads(XX); ******
// ****************** Declarying class of Probability Pointer *************************
......@@ -720,27 +727,45 @@ int bioem::createConvolutedProjectionMap(int iMap, int iConv, mycomplex_t* lproj
// cout << "GG " << i << " " << j << " " << refCTF[i][0] << " " << refCTF[i][1] <<" " <<lproj[i][0] <<" " <<lproj[i][1] << "\n";
}
//FFTW_C2R will destroy the input array, so we have to work on a copy here
memcpy(tmp, localmultFFT, sizeof(mycomplex_t) * param.param_device.NumberPixels * param.param_device.NumberFFTPixels1D);
// **** Bringing convoluted Map to real Space ****
myfftw_execute_dft_c2r(param.fft_plan_c2r_backward, tmp, Mapconv);
// *** Calculating Cross-correlations of cal-convoluted map with its self *****
sumC = 0;
sumC = localmultFFT[0][0];
sumsquareC = 0;
for(int i = 0; i < param.param_device.NumberPixels * param.param_device.NumberPixels; i++)
if (FFTAlgo)
{
sumC += Mapconv[i];
sumsquareC += Mapconv[i] * Mapconv[i];
for(int i = 0; i < param.param_device.NumberPixels; i++)
{
for (int j = 1;j < param.param_device.NumberFFTPixels1D;j++)
{
int k = i * param.param_device.NumberFFTPixels1D + j;
sumsquareC += (localmultFFT[k][0] * localmultFFT[k][0] + localmultFFT[k][1] * localmultFFT[k][1]) * 2;
}
int k = i * param.param_device.NumberFFTPixels1D;
sumsquareC += localmultFFT[k][0] * localmultFFT[k][0] + localmultFFT[k][1] * localmultFFT[k][1];
}
myfloat_t norm2 = (myfloat_t) (param.param_device.NumberPixels * param.param_device.NumberPixels);
sumsquareC = sumsquareC / norm2;
}
// *** The DTF gives an unnormalized value so have to divded by the total number of pixels in Fourier ***
// Normalizing
myfloat_t norm2 = (myfloat_t) (param.param_device.NumberPixels * param.param_device.NumberPixels);
myfloat_t norm4 = norm2 * norm2;
sumC = sumC / norm2;
sumsquareC = sumsquareC / norm4;
else
{
//FFTW_C2R will destroy the input array, so we have to work on a copy here
memcpy(tmp, localmultFFT, sizeof(mycomplex_t) * param.param_device.NumberPixels * param.param_device.NumberFFTPixels1D);
// **** Bringing convoluted Map to real Space ****
myfftw_execute_dft_c2r(param.fft_plan_c2r_backward, tmp, Mapconv);
for(int i = 0; i < param.param_device.NumberPixels * param.param_device.NumberPixels; i++)
{
sumsquareC += Mapconv[i] * Mapconv[i];
sumC += Mapconv[i];
}
myfloat_t norm2 = (myfloat_t) (param.param_device.NumberPixels * param.param_device.NumberPixels);
myfloat_t norm4 = norm2 * norm2;
sumsquareC = sumsquareC / norm4;
}
cuda_custom_timeslot_end;
return(0);
......
......@@ -64,6 +64,10 @@ public:
myfloat3_t* angles;
int nTotGridAngles;
int nTotCTFs;
bool printModel;
int printModelOrientation;
int printModelConvolution;
int fft_plans_created;
myfftw_plan fft_plan_c2c_forward, fft_plan_c2c_backward, fft_plan_r2c_forward, fft_plan_c2r_backward;
......
......@@ -43,6 +43,8 @@ bioem_param::bioem_param()
refCTF = NULL;
CtfParam = NULL;
angles = NULL;
printModel = false;
}
int bioem_param::readParameters(const char* fileinput)
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment