param.cpp 16.6 KB
Newer Older
Pilar Cossio's avatar
License  
Pilar Cossio committed
1 2 3 4
/* ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
        < BioEM software for Bayesian inference of Electron Microscopy images>
            Copyright (C) 2014 Pilar Cossio, David Rohr and Gerhard Hummer.
            Max Planck Institute of Biophysics, Frankfurt, Germany.
5

Pilar Cossio's avatar
License  
Pilar Cossio committed
6 7 8 9
                See license statement for terms of distribution.

   ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/

qon's avatar
qon committed
10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
#include <stdio.h>
#include <stdlib.h>
#include <iostream>
#include <fstream>
#include <cstring>
#include <math.h>
#include <fftw3.h>

#include "param.h"
#include "map.h"

using namespace std;

bioem_param::bioem_param()
{
25
	//Number of Pixels
26
	param_device.NumberPixels = 0;
27
	param_device.NumberFFTPixels1D = 0;
28 29 30 31 32 33 34 35 36
	// Euler angle grid spacing
	angleGridPointsAlpha = 0;
	angleGridPointsBeta = 0;
	//Envelop function paramters
	numberGridPointsEnvelop = 0;
	//Contrast transfer function paramters
	numberGridPointsCTF_amp = 0;
	numberGridPointsCTF_phase = 0;

David Rohr's avatar
David Rohr committed
37
	// ****center displacement paramters Equal in both directions***
38 39
	param_device.maxDisplaceCenter = 0;
	numberGridPointsDisplaceCenter = 0;
40 41

	fft_plans_created = 0;
42 43 44

	refCTF = NULL;
	CtfParam = NULL;
David Rohr's avatar
David Rohr committed
45
	angles = NULL;
qon's avatar
qon committed
46 47
}

48
int bioem_param::readParameters(const char* fileinput)
qon's avatar
qon committed
49
{
David Rohr's avatar
David Rohr committed
50 51 52
	// **************************************************************************************
	// ***************************** Reading Input Parameters ******************************
	// **************************************************************************************
53

David Rohr's avatar
David Rohr committed
54 55 56 57 58 59 60 61 62 63 64 65 66 67 68
	// Control for Parameters
	bool yesPixSi = false;
	bool yesNumPix = false;
	bool yesGPal = false;
	bool yesGPbe = false;
	bool yesGPEnv = false;
	bool yesGPamp = false;
	bool yesGPpha = false;
	bool yesSTEnv = false;
	bool yesSTamp = false;
	bool yesSTpha = false;
	bool yesGSPamp = false ;
	bool yesGSPEnv = false ;
	bool yesGSPpha = false ;
	bool yesMDC = false ;
69
	bool yesGCen = false ;
70

71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86
	ifstream input(fileinput);
	if (!input.good())
	{
		cout << "Failed to open file: " << fileinput << "\n";
		exit(0);
	}

	char line[512] = {0};
	char saveline[512];

	cout << "\n +++++++++++++++++++++++++++++++++++++++++ \n";
	cout << "\n		   READING PARAMETERS  \n\n";

	cout << " +++++++++++++++++++++++++++++++++++++++++ \n";
	while (!input.eof())
	{
87 88 89
		input.getline(line, 512);
		strcpy(saveline, line);
		char *token = strtok(line, " ");
90

91
		if (token == NULL || line[0] == '#' || strlen(token) == 0)
92 93 94
		{
			// comment or blank line
		}
95
		else if (strcmp(token, "PIXEL_SIZE") == 0)
96
		{
97 98
			token = strtok(NULL, " ");
			pixelSize = atof(token);
David Rohr's avatar
David Rohr committed
99 100 101
			if (pixelSize < 0 ) { cout << "*** Error: Negative pixelSize "; exit(1);}
			cout << "Pixel Sixe " << pixelSize << "\n";
			yesPixSi= true;
102
		}
103
		else if (strcmp(token, "NUMBER_PIXELS") == 0)
104
		{
105 106
			token = strtok(NULL, " ");
			param_device.NumberPixels = int(atoi(token));
David Rohr's avatar
David Rohr committed
107
			if (param_device.NumberPixels < 0 ) { cout << "*** Error: Negative Number of Pixels "; exit(1);}
108
			cout << "Number of Pixels " << param_device.NumberPixels << "\n";
David Rohr's avatar
David Rohr committed
109
			yesNumPix= true ;
110
		}
111
		else if (strcmp(token, "GRIDPOINTS_ALPHA") == 0)
112
		{
113 114
			token = strtok(NULL, " ");
			angleGridPointsAlpha = int(atoi(token));
David Rohr's avatar
David Rohr committed
115 116 117
			if (angleGridPointsAlpha < 0 ) { cout << "*** Error: Negative GRIDPOINTS_ALPHA "; exit(1);}
			cout << "Grid points alpha " << angleGridPointsAlpha << "\n";
			yesGPal= true;
118
		}
119
		else if (strcmp(token, "GRIDPOINTS_BETA") == 0)
120
		{
121 122
			token = strtok(NULL, " ");
			angleGridPointsBeta = int(atoi(token));
David Rohr's avatar
David Rohr committed
123
			if (angleGridPointsBeta < 0 ) { cout << "*** Error: Negative GRIDPOINTS_BETA "; exit(1);}
124
			cout << "Grid points in Cosine ( beta ) " << angleGridPointsBeta << "\n";
David Rohr's avatar
David Rohr committed
125
			yesGPbe= true;
126 127
		}

128
		else if (strcmp(token, "GRIDPOINTS_ENVELOPE") == 0)
129
		{
130 131
			token = strtok(NULL, " ");
			numberGridPointsEnvelop = int(atoi(token));
David Rohr's avatar
David Rohr committed
132 133 134
			if (numberGridPointsDisplaceCenter < 0 ) { cout << "*** Error: Negative GRIDPOINTS_ENVELOPE "; exit(1);}
			cout << "Grid points envelope " << numberGridPointsEnvelop << "\n";
			yesGPEnv = true;
135
		}
136
		else if (strcmp(token, "START_ENVELOPE") == 0)
137
		{
138 139
			token = strtok(NULL, " ");
			startGridEnvelop = atof(token);
David Rohr's avatar
David Rohr committed
140 141 142
			if (startGridEnvelop < 0 ) { cout << "*** Error: Negative START_ENVELOPE "; exit(1);}
			cout << "Start Envelope " << startGridEnvelop << "\n";
			yesSTEnv = true ;
143
		}
144
		else if (strcmp(token, "GRIDSPACE_ENVELOPE") == 0)
145
		{
146 147
			token = strtok(NULL, " ");
			gridEnvelop = atof(token);
David Rohr's avatar
David Rohr committed
148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193
			if (gridEnvelop < 0 ) { cout << "*** Error: Negative GRIDSPACE_ENVELOPE "; exit(1);}
			cout << "Grid spacing Envelope " << gridEnvelop << "\n";
			yesGSPEnv = true ;
		}
		else if (strcmp(token,"GRIDPOINTS_PSF_PHASE")==0)
		{
			token = strtok(NULL," ");
			numberGridPointsCTF_phase=int(atoi(token));
			cout << "Grid points PSF " << numberGridPointsCTF_phase << "\n";
			yesGPpha = true;
		}
		else if (strcmp(token,"START_PSF_PHASE")==0)
		{
			token = strtok(NULL," ");
			startGridCTF_phase=atof(token);
			cout << "Start PSF " << startGridCTF_phase << "\n";
			yesSTpha = true ;
		}
		else if (strcmp(token,"GRIDSPACE_PSF_PHASE")==0)
		{
			token = strtok(NULL," ");
			gridCTF_phase=atof(token);
			cout << "Grid Space PSF " << gridCTF_phase << "\n";
			yesGSPpha = true ;
		}
		else if (strcmp(token,"GRIDPOINTS_PSF_AMP")==0)
		{
			token = strtok(NULL," ");
			numberGridPointsCTF_amp=int(atoi(token));
			cout << "Grid points PSF " << numberGridPointsCTF_amp << "\n";
			yesGPamp = true ;
		}
		else if (strcmp(token,"START_PSF_AMP")==0)
		{
			token = strtok(NULL," ");
			startGridCTF_amp=atof(token);
			cout << "Start PSF " << startGridCTF_amp << "\n";
			yesSTamp = true ;
		}
		else if (strcmp(token,"GRIDSPACE_PSF_AMP")==0)
		{
			token = strtok(NULL," ");
			gridCTF_amp=atof(token);
			cout << "Grid Space PSF " << gridCTF_amp << "\n";
			yesGSPamp = true ;
		}
194
		else if (strcmp(token, "MAX_D_CENTER") == 0)
195
		{
196 197
			token = strtok(NULL, " ");
			param_device.maxDisplaceCenter = int(atoi(token));
David Rohr's avatar
David Rohr committed
198 199 200
			if (param_device.maxDisplaceCenter < 0 ) { cout << "*** Error: Negative MAX_D_CENTER "; exit(1);}
			cout << "Maximum displacement Center " <<  param_device.maxDisplaceCenter << "\n";
			yesMDC = true;
201
		}
202
		else if (strcmp(token, "PIXEL_GRID_CENTER") == 0)
203
		{
204 205
			token = strtok(NULL, " ");
			param_device.GridSpaceCenter = int(atoi(token));
David Rohr's avatar
David Rohr committed
206 207 208
			if (param_device.GridSpaceCenter < 0 ) { cout << "*** Error: Negative PIXEL_GRID_CENTER "; exit(1);}
			cout << "Grid space displacement center " <<   param_device.GridSpaceCenter << "\n";
			yesGCen = true;
209
		}
210
		else if (strcmp(token, "WRITE_PROB_ANGLES") == 0)
211
		{
212
			param_device.writeAngles = true;
213 214 215 216
			cout << "Writing Probabilies of each angle \n";
		}
	}
	input.close();
217

218
	// Checks for ALL INPUT
David Rohr's avatar
David Rohr committed
219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235
	if( not ( yesPixSi ) ){ cout << "**** INPUT MISSING: Please provide PIXEL_SIZE\n" ; exit (1);};
	if( not ( yesNumPix ) ){ cout << "**** INPUT MISSING: Please provide NUMBER_PIXELS \n" ; exit (1);};
	if( not ( yesGPal ) ) { cout << "**** INPUT MISSING: Please provide GRIDPOINTS_ALPHA \n" ; exit (1);};
	if( not ( yesGPbe ) ) { cout << "**** INPUT MISSING: Please provide GRIDPOINTS_BETA \n" ; exit (1);};
	if( not (  yesGPEnv ) ) { cout << "**** INPUT MISSING: Please provide GRIDPOINTS_ENVELOPE \n" ; exit (1);};
	if( not (  yesGPamp ) ) { cout << "**** INPUT MISSING: Please provide GRIDPOINTS_PSF_AMP \n" ; exit (1);};
	if( not (  yesGPpha ) ) { cout << "**** INPUT MISSING: Please provide GRIDPOINTS_PSF_PHASE \n" ; exit (1);};
	if( not (  yesSTEnv  ) ) { cout << "**** INPUT MISSING: Please provide START_ENVELOPE \n" ; exit (1);};
	if( not (  yesSTamp  ) ) { cout << "**** INPUT MISSING: Please provide START_PSF_AMP \n" ; exit (1);};
	if( not (  yesSTpha  ) ) { cout << "**** INPUT MISSING: Please provide START_PSF_PHASE \n" ; exit (1);};
	if( not (  yesGSPamp  ) ) { cout << "**** INPUT MISSING: Please provide GRIDSPACE_PSF_AMP \n" ; exit (1);};
	if( not ( yesGSPEnv  ) ) { cout << "**** INPUT MISSING: Please provide GRIDSPACE_ENVELOPE \n" ; exit (1);};
	if( not ( yesGSPpha  ) ) { cout << "**** INPUT MISSING: Please provide GRIDSPACE_PSF_PHASE \n" ; exit (1);};
	if( not (  yesMDC  ) ) { cout << "**** INPUT MISSING: Please provide MAX_D_CENTER \n" ; exit (1);};
	if( not (  yesGCen  ) ) { cout << "**** INPUT MISSING: Please provide PIXEL_GRID_CENTER \n" ; exit (1);};


236 237 238 239
    //if only one grid point for PSF kernel:
    if( (myfloat_t) numberGridPointsCTF_amp == 1 ) gridCTF_amp = startGridCTF_amp;
    if( (myfloat_t) numberGridPointsCTF_phase == 1 ) gridCTF_phase = startGridCTF_phase;
    if( (myfloat_t) numberGridPointsEnvelop == 1 ) gridEnvelop = startGridEnvelop;
240 241

	//More checks with input parameters
242
	// Envelope should not have a standard deviation greater than Npix/2
David Rohr's avatar
David Rohr committed
243
	if(sqrt(1./( (myfloat_t) numberGridPointsDisplaceCenter  * gridEnvelop + startGridEnvelop))> float(param_device.NumberPixels)/2.0) {
244
		cout << "MAX Standard deviation of envelope is larger than Allowed KERNEL Length\n";
David Rohr's avatar
David Rohr committed
245 246 247 248
		exit(1);
	}
	// Envelop param should be positive
	if( startGridCTF_amp < 0 || startGridCTF_amp > 1){
David Rohr's avatar
David Rohr committed
249
		cout << "PSF Amplitud should be between 0 and 1\n";
David Rohr's avatar
David Rohr committed
250 251 252
		exit(1);
	}

253
	if( (myfloat_t) numberGridPointsCTF_amp * gridCTF_amp + startGridCTF_amp < 0 || (myfloat_t) (numberGridPointsCTF_amp - 1) * gridCTF_amp + startGridCTF_amp > 1){
David Rohr's avatar
David Rohr committed
254
		cout << "Value should be between 0 and 1\n" ;
David Rohr's avatar
David Rohr committed
255 256 257
		exit(1);
	}

258

259 260 261 262 263 264 265 266 267 268
	param_device.NumberFFTPixels1D = param_device.NumberPixels / 2 + 1;
	FFTMapSize = param_device.NumberPixels * param_device.NumberFFTPixels1D;

	//Does currently not work with custom alignment on GPU
	/*if (FFTMapSize % Alignment)
	{
		FFTMapSize += Alignment - FFTMapSize % Alignment;
	}
	cout << "Using MAP Size " << FFTMapSize << " (Alignment " << Alignment << ", Unaligned Size " << param_device.NumberPixels * param_device.NumberFFTPixels1D << ")\n";*/

269
	cout << " +++++++++++++++++++++++++++++++++++++++++ \n";
270

271 272 273 274 275
	return(0);
}

void bioem_param::PrepareFFTs()
{
276 277 278 279 280
	cout << "Preparing FFTs\n";
	releaseFFTPlans();
	mycomplex_t *tmp_map, *tmp_map2;
	tmp_map = (mycomplex_t *) myfftw_malloc(sizeof(mycomplex_t) * param_device.NumberPixels * param_device.NumberPixels);
	tmp_map2 = (mycomplex_t *) myfftw_malloc(sizeof(mycomplex_t) * param_device.NumberPixels * param_device.NumberPixels);
281
	Alignment = 64;
282

283 284 285 286
	fft_plan_c2c_forward = myfftw_plan_dft_2d(param_device.NumberPixels, param_device.NumberPixels, tmp_map, tmp_map2, FFTW_FORWARD, FFTW_MEASURE | FFTW_DESTROY_INPUT);
	fft_plan_c2c_backward = myfftw_plan_dft_2d(param_device.NumberPixels, param_device.NumberPixels, tmp_map, tmp_map2, FFTW_BACKWARD, FFTW_MEASURE | FFTW_DESTROY_INPUT);
	fft_plan_r2c_forward = myfftw_plan_dft_r2c_2d(param_device.NumberPixels, param_device.NumberPixels, (myfloat_t*) tmp_map, tmp_map2, FFTW_MEASURE | FFTW_DESTROY_INPUT);
	fft_plan_c2r_backward = myfftw_plan_dft_c2r_2d(param_device.NumberPixels, param_device.NumberPixels, tmp_map, (myfloat_t*) tmp_map2, FFTW_MEASURE | FFTW_DESTROY_INPUT);
287

288
	if (fft_plan_c2c_forward == 0 || fft_plan_c2c_backward == 0 || fft_plan_r2c_forward == 0 || fft_plan_c2r_backward == 0)
289 290 291 292 293 294 295
	{
		cout << "Error planing FFTs\n";
		exit(1);
	}

	myfftw_free(tmp_map);
	myfftw_free(tmp_map2);
296
	fft_plans_created = 1;	
qon's avatar
qon committed
297 298
}

299 300 301 302 303 304
void bioem_param::releaseFFTPlans()
{
	if (fft_plans_created)
	{
		myfftw_destroy_plan(fft_plan_c2c_forward);
		myfftw_destroy_plan(fft_plan_c2c_backward);
305 306
		myfftw_destroy_plan(fft_plan_r2c_forward);
		myfftw_destroy_plan(fft_plan_c2r_backward);
307
		myfftw_cleanup();
308 309 310 311
	}
	fft_plans_created = 0;
}

qon's avatar
qon committed
312 313
int bioem_param::CalculateGridsParam() //TO DO FOR QUATERNIONS
{
David Rohr's avatar
David Rohr committed
314 315 316
	// **************************************************************************************
	// **************** Routine that pre-calculates Euler angle grids **********************
	// ************************************************************************************
317 318
	myfloat_t grid_alpha, cos_grid_beta;
	int n = 0;
319 320

	//alpha and gamma are uniform in -PI,PI
321
	grid_alpha = 2.f * M_PI / (myfloat_t) angleGridPointsAlpha;
322 323

	//cosine beta is uniform in -1,1
324
	cos_grid_beta = 2.f / (myfloat_t) angleGridPointsBeta;
325 326

	// Euler Angle Array
327
	delete[] angles;
David Rohr's avatar
David Rohr committed
328
	angles = new myfloat3_t[angleGridPointsAlpha * angleGridPointsBeta * angleGridPointsAlpha];
329 330 331 332 333 334
	for (int ialpha = 0; ialpha < angleGridPointsAlpha; ialpha ++)
	{
		for (int ibeta = 0; ibeta < angleGridPointsBeta; ibeta ++)
		{
			for (int igamma = 0; igamma < angleGridPointsAlpha; igamma ++)
			{
335 336 337
				angles[n].pos[0] = (myfloat_t) ialpha * grid_alpha - M_PI + grid_alpha * 0.5f; //ALPHA centered in the middle
				angles[n].pos[1] = acos((myfloat_t) ibeta * cos_grid_beta - 1 + cos_grid_beta * 0.5f); //BETA centered in the middle
				angles[n].pos[2] = (myfloat_t) igamma * grid_alpha - M_PI + grid_alpha * 0.5f; //GAMMA centered in the middle
338 339 340 341
				n++;
			}
		}
	}
342
	nTotGridAngles = n;
343

344

David Rohr's avatar
David Rohr committed
345
	// ********** Calculating normalized volumen element *********
346

347
	param_device.volu = grid_alpha * grid_alpha * cos_grid_beta * (myfloat_t) param_device.GridSpaceCenter * pixelSize * (myfloat_t) param_device.GridSpaceCenter * pixelSize
348 349
						* gridCTF_phase * gridCTF_amp * gridEnvelop / (2.f * M_PI) / (2.f * M_PI) / 2.f / (2.f * (myfloat_t) param_device.maxDisplaceCenter) / (2.f * (myfloat_t) param_device.maxDisplaceCenter) / ((myfloat_t) numberGridPointsCTF_amp * gridCTF_amp )
						/ ((myfloat_t) numberGridPointsCTF_phase * gridCTF_phase) / ((myfloat_t) numberGridPointsEnvelop * gridEnvelop );
350

351
      //  cout << "VOLU " << param_device.volu  << " " << gridCTF_amp << "\n";
David Rohr's avatar
David Rohr committed
352
	// *** Number of total pixels***
353

354 355
	param_device.Ntotpi = (myfloat_t) (param_device.NumberPixels * param_device.NumberPixels);
	param_device.NtotDist = (2 * (int) (param_device.maxDisplaceCenter / param_device.GridSpaceCenter) + 1 ) * (2 * (int) (param_device.maxDisplaceCenter / param_device.GridSpaceCenter) + 1);
356 357

	return(0);
qon's avatar
qon committed
358 359 360 361 362

}

int bioem_param::CalculateRefCTF()
{
David Rohr's avatar
David Rohr committed
363 364 365
	// **************************************************************************************
	// ********** Routine that pre-calculates Kernels for Convolution **********************
	// ************************************************************************************
366

367
	myfloat_t amp, env, phase, ctf, radsq;
368
	myfloat_t* localCTF;
David Rohr's avatar
David Rohr committed
369
	mycomplex_t* localout;
370 371
	int nctfmax = param_device.NumberPixels / 2;
	int n = 0;
372

373
	localCTF = (myfloat_t *) myfftw_malloc(sizeof(myfloat_t) * param_device.NumberPixels * param_device.NumberPixels);
David Rohr's avatar
David Rohr committed
374
	localout = (mycomplex_t *) myfftw_malloc(sizeof(mycomplex_t) * param_device.NumberPixels * param_device.NumberFFTPixels1D);
375 376 377

	nTotCTFs = numberGridPointsCTF_amp * numberGridPointsCTF_phase * numberGridPointsEnvelop;
	delete[] refCTF;
378
	refCTF = new mycomplex_t[nTotCTFs * FFTMapSize];
379 380
	delete[] CtfParam;
	CtfParam = new myfloat3_t[nTotCTFs];
381 382 383

	for (int iamp = 0; iamp <  numberGridPointsCTF_amp ; iamp++) //Loop over amplitud
	{
384
		amp = (myfloat_t) iamp * gridCTF_amp + startGridCTF_amp;
385 386 387

		for (int iphase = 0; iphase <  numberGridPointsCTF_phase ; iphase++)//Loop over phase
		{
388
			phase = (myfloat_t) iphase * gridCTF_phase + startGridCTF_phase;
389

390
			for (int ienv = 0; ienv <  numberGridPointsEnvelop ; ienv++)//Loop over envelope
391
			{
392
				env = (myfloat_t) ienv * gridEnvelop + startGridEnvelop;
393

394
				memset(localCTF, 0, param_device.NumberPixels * param_device.NumberPixels * sizeof(myfloat_t));
395 396

				//Assigning CTF values
397
				for(int i = 0; i < nctfmax; i++)
398
				{
399
					for(int j = 0; j < nctfmax; j++)
400
					{
401 402
						radsq = (myfloat_t) (i * i + j * j) * pixelSize * pixelSize;
						ctf = exp(-radsq * env / 2.0) * (amp * cos(radsq * phase / 2.0) - sqrt((1 - amp * amp)) * sin(radsq * phase / 2.0));
403

404 405 406 407
						localCTF[i * param_device.NumberPixels + j] = (myfloat_t) ctf;
						localCTF[i * param_device.NumberPixels + param_device.NumberPixels - j - 1] = (myfloat_t) ctf;
						localCTF[(param_device.NumberPixels - i - 1)*param_device.NumberPixels + j] = (myfloat_t) ctf;
						localCTF[(param_device.NumberPixels - i - 1)*param_device.NumberPixels + param_device.NumberPixels - j - 1] = (myfloat_t) ctf;
408 409 410
					}
				}
				//Calling FFT_Forward
411
				myfftw_execute_dft_r2c(fft_plan_r2c_forward, localCTF, localout);
412 413

				// Normalizing and saving the Reference CTFs
414
				mycomplex_t* curRef = &refCTF[n * FFTMapSize];
415
				for(int i = 0; i < param_device.NumberPixels * param_device.NumberFFTPixels1D; i++ )
416
				{
417 418
					curRef[i][0] = localout[i][0];
					curRef[i][1] = localout[i][1];
419
				}
420 421 422
				CtfParam[n].pos[0] = amp;
				CtfParam[n].pos[1] = phase;
				CtfParam[n].pos[2] = env;
423 424 425 426 427 428
				n++;
			}
		}
	}

	myfftw_free(localCTF);
David Rohr's avatar
David Rohr committed
429
	myfftw_free(localout);
430 431 432 433 434
	if (nTotCTFs != n)
	{
		cout << "Internal error during CTF preparation\n";
		exit(1);
	}
435 436

	return(0);
qon's avatar
qon committed
437 438 439 440
}

bioem_param::~bioem_param()
{
441
	releaseFFTPlans();
442
	param_device.NumberPixels = 0;
443 444 445 446 447 448 449
	angleGridPointsAlpha = 0;
	angleGridPointsBeta = 0;
	numberGridPointsEnvelop = 0;
	numberGridPointsCTF_amp = 0;
	numberGridPointsCTF_phase = 0;
	param_device.maxDisplaceCenter = 0;
	numberGridPointsDisplaceCenter = 0;
David Rohr's avatar
David Rohr committed
450 451 452
	if (refCTF) delete[] refCTF;
	if (CtfParam) delete[] CtfParam;
	if (angles) delete[] angles;
qon's avatar
qon committed
453
}