Commit 1d1f08c6 authored by David Rohr's avatar David Rohr
Browse files

cleanups, replace float by myfloat_t where applicable

parent e6d44f5b
......@@ -418,22 +418,22 @@ inline int bioem::calculateCCFFT(int iRefMap, int iOrient, int iConv, myfloat_t
{
for (int cent_y = 0; cent_y <= param.param_device.maxDisplaceCenter; cent_y=cent_y+param.param_device.GridSpaceCenter)
{
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, float(lCC[cent_x*param.param_device.NumberPixels+cent_y][0])/float(param.param_device.NumberPixels * param.param_device.NumberPixels), cent_x, cent_y);
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, (myfloat_t) lCC[cent_x*param.param_device.NumberPixels+cent_y][0]/ (myfloat_t) (param.param_device.NumberPixels * param.param_device.NumberPixels), cent_x, cent_y);
}
for (int cent_y = param.param_device.NumberPixels-param.param_device.maxDisplaceCenter; cent_y < param.param_device.NumberPixels; cent_y=cent_y+param.param_device.GridSpaceCenter)
{
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, float(lCC[cent_x*param.param_device.NumberPixels+cent_y][0])/float(param.param_device.NumberPixels*param.param_device.NumberPixels), cent_x, param.param_device.NumberPixels-cent_y);
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, (myfloat_t) lCC[cent_x*param.param_device.NumberPixels+cent_y][0]/ (myfloat_t) (param.param_device.NumberPixels*param.param_device.NumberPixels), cent_x, param.param_device.NumberPixels-cent_y);
}
}
for (int cent_x = param.param_device.NumberPixels-param.param_device.maxDisplaceCenter; cent_x < param.param_device.NumberPixels; cent_x=cent_x+param.param_device.GridSpaceCenter)
{
for (int cent_y = 0; cent_y < param.param_device.maxDisplaceCenter; cent_y=cent_y+param.param_device.GridSpaceCenter)
{
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, float(lCC[cent_x*param.param_device.NumberPixels+cent_y][0])/float(param.param_device.NumberPixels*param.param_device.NumberPixels), param.param_device.NumberPixels-cent_x, cent_y);
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, (myfloat_t) lCC[cent_x*param.param_device.NumberPixels+cent_y][0]/ (myfloat_t) (param.param_device.NumberPixels*param.param_device.NumberPixels), param.param_device.NumberPixels-cent_x, cent_y);
}
for (int cent_y = param.param_device.NumberPixels-param.param_device.maxDisplaceCenter; cent_y <= param.param_device.NumberPixels; cent_y=cent_y+param.param_device.GridSpaceCenter)
{
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, float(lCC[cent_x*param.param_device.NumberPixels+cent_y][0])/float(param.param_device.NumberPixels*param.param_device.NumberPixels), param.param_device.NumberPixels-cent_x, param.param_device.NumberPixels-cent_y);
calProb(iRefMap, iOrient, iConv, sumC, sumsquareC, (myfloat_t) lCC[cent_x*param.param_device.NumberPixels+cent_y][0]/ (myfloat_t) (param.param_device.NumberPixels*param.param_device.NumberPixels), param.param_device.NumberPixels-cent_x, param.param_device.NumberPixels-cent_y);
}
}
......@@ -491,7 +491,6 @@ inline int bioem::calProb(int iRefMap,int iOrient, int iConv,myfloat_t sumC,myfl
int bioem::createProjection(int iMap,mycomplex_t* mapFFT)
{
/**************************************************************************************/
/**** BioEM Create Projection routine in Euler angle predefined grid****************
********************* and turns projection into Fourier space **********************/
......@@ -547,8 +546,8 @@ int bioem::createProjection(int iMap,mycomplex_t* mapFFT)
for(int n=0; n< Model.nPointsModel; n++)
{
//Getting pixel that represents coordinates & shifting the start at to Numpix/2,Numpix/2 )
i=floor(RotatedPointsModel[n].pos[0]/param.pixelSize+float(param.param_device.NumberPixels)/2.0+0.5);
j=floor(RotatedPointsModel[n].pos[1]/param.pixelSize+float(param.param_device.NumberPixels)/2.0+0.5);
i=floor(RotatedPointsModel[n].pos[0]/param.pixelSize+ (myfloat_t) param.param_device.NumberPixels / 2.0f + 0.5f);
j=floor(RotatedPointsModel[n].pos[1]/param.pixelSize+ (myfloat_t) param.param_device.NumberPixels / 2.0f + 0.5f);
localproj[i*param.param_device.NumberPixels+j][0]+=Model.densityPointsModel[n]/Model.NormDen;
}
......@@ -627,8 +626,8 @@ int bioem::createConvolutedProjectionMap(int iMap,int iConv,mycomplex_t* lproj,b
}
/*** The DTF gives an unnormalized value so have to divded by the total number of pixels in Fourier ***/
// Normalizing
sumC=sumC/float(param.param_device.NumberPixels*param.param_device.NumberPixels);
sumsquareC=sumsquareC/pow(float(param.param_device.NumberPixels),4);
sumC=sumC/ (myfloat_t) (param.param_device.NumberPixels*param.param_device.NumberPixels);
sumsquareC=sumsquareC / pow((myfloat_t) param.param_device.NumberPixels,4);
/**** Freeing fftw_complex created (dont know if omp critical is necessary) ****/
myfftw_free(localconvFFT);
......
......@@ -3,7 +3,7 @@
//#include <boost/iterator/iterator_concepts.hpp>
#ifndef BIOEM_GPUCODE
//#define SSECODE //Explicit SSE code, not correct yet since loop counter is assumed multiple of 4, anyway not faster than autovectorized code.
//#define SSECODE //Explicit SSE code, not correct yet since loop counter is assumed multiple of 4, anyway not faster than autovectorized code, only implemented for float, not for double.
#endif
#ifdef SSECODE
......@@ -156,7 +156,7 @@ __device__ static inline void compareRefMap(const int iRefMap, const int iOrient
if (myShift < 32) //Warp Size is 32, threads are synched automatically
{
volatile float* vbuf = buf; //Mem must be volatile such that memory access is not reordered
volatile myfloat_t* vbuf = buf; //Mem must be volatile such that memory access is not reordered
if (nShifts2 >= 64 && vbuf[myThreadIdxX + 32] > vbuf[myThreadIdxX]) vbuf[myThreadIdxX] = vbuf[myThreadIdxX + 32];
if (nShifts2 >= 32 && vbuf[myThreadIdxX + 16] > vbuf[myThreadIdxX]) vbuf[myThreadIdxX] = vbuf[myThreadIdxX + 16];
if (nShifts2 >= 16 && vbuf[myThreadIdxX + 8] > vbuf[myThreadIdxX]) vbuf[myThreadIdxX] = vbuf[myThreadIdxX + 8];
......@@ -245,8 +245,8 @@ __device__ static inline void compareRefMap(const int iRefMap, const int iOrient
if (myShift < 32) //Warp Size is 32, threads are synched automatically
{
volatile float* vbuf = buf; //Mem must be volatile such that memory access is not reordered
volatile float* vbuf2 = buf2;
volatile myfloat_t* vbuf = buf; //Mem must be volatile such that memory access is not reordered
volatile myfloat_t* vbuf2 = buf2;
if (nShifts2 >= 64)
{
vbuf[myThreadIdxX] += vbuf[myThreadIdxX + 32];
......
#ifndef BIOEM_DEFS_H
#define BIOEM_DEFS_H
//#define BIOEM_USE_DOUBLE
#ifndef BIOEM_USE_DOUBLE
typedef float myfloat_t;
#define myfftw_malloc fftwf_malloc
......
......@@ -80,7 +80,7 @@ int bioem_RefMap::readRefMaps(bioem_param& param)
else
{
int i,j;
myfloat_t z;
float z;
char tmpVals[36] = {' '};
......@@ -95,7 +95,7 @@ int bioem_RefMap::readRefMaps(bioem_param& param)
//checking for Map limits
if(i>0 && i-1 <BIOEM_MAP_SIZE_X&& j>0 && j-1<BIOEM_MAP_SIZE_Y && nummap < BIOEM_MAX_MAPS)
{
Ref[nummap].points[i-1][j-1]=z;
Ref[nummap].points[i-1][j-1] = (myfloat_t) z;
lasti=i;
lastj=j;
}
......
......@@ -60,9 +60,9 @@ int bioem_model::readModel()
char name[5] = {"PPP"};
char resName[3] = {' '};
float x=0.0;
float y=0.0;
float z=0.0;
float x = 0.0;
float y = 0.0;
float z = 0.0;
char tmp[6] = {' '};
......@@ -99,9 +99,9 @@ int bioem_model::readModel()
NormDen+=densityPointsModel[numres];
//Getting the coordinates
PointsModel[numres].pos[0]=x;
PointsModel[numres].pos[1]=y;
PointsModel[numres].pos[2]=z;
PointsModel[numres].pos[0] = (myfloat_t) x;
PointsModel[numres].pos[1] = (myfloat_t) y;
PointsModel[numres].pos[2] = (myfloat_t) z;
exampleReadCoor << "RESIDUE " << numres << " " << PointsModel[numres].pos[0] << " " << PointsModel[numres].pos[1] << " " << PointsModel[numres].pos[2] << " " <<densityPointsModel[numres] << "\n";
numres++;
}
......@@ -128,7 +128,14 @@ int bioem_model::readModel()
cout << "BIOEM_MODEL_SIZE too small\n";
exit(1);
}
sscanf(line, "%f %f %f %f %f",&PointsModel[numres].pos[0],&PointsModel[numres].pos[1],&PointsModel[numres].pos[2],&radiusPointsModel[numres],&densityPointsModel[numres]);
float tmpval[5];
sscanf(line, "%f %f %f %f %f", &tmpval[0], &tmpval[1], &tmpval[2], &tmpval[3], &tmpval[4]);
PointsModel[numres].pos[0] = (myfloat_t) tmpval[0];
PointsModel[numres].pos[1] = (myfloat_t) tmpval[1];
PointsModel[numres].pos[2] = (myfloat_t) tmpval[2];
radiusPointsModel[numres] = (myfloat_t) tmpval[3];
densityPointsModel[numres] = (myfloat_t) tmpval[4];
exampleReadCoor << "RESIDUE " << numres << " " << PointsModel[numres].pos[0] << " " << PointsModel[numres].pos[1] << " " << PointsModel[numres].pos[2] << " " <<densityPointsModel[numres]<< "\n";
NormDen+=densityPointsModel[numres];
numres++;
......@@ -151,9 +158,9 @@ int bioem_model::readModel()
r_cm.pos[1] += PointsModel[n].pos[1];
r_cm.pos[2] += PointsModel[n].pos[2];
}
r_cm.pos[0]=r_cm.pos[0]/float(nPointsModel);
r_cm.pos[1]=r_cm.pos[1]/float(nPointsModel);
r_cm.pos[2]=r_cm.pos[2]/float(nPointsModel);
r_cm.pos[0]=r_cm.pos[0] / (myfloat_t) nPointsModel;
r_cm.pos[1]=r_cm.pos[1] / (myfloat_t) nPointsModel;
r_cm.pos[2]=r_cm.pos[2] / (myfloat_t) nPointsModel;
for(int n=0; n < nPointsModel; n++)
{
PointsModel[n].pos[0]-= r_cm.pos[0];
......
......@@ -219,22 +219,21 @@ int bioem_param::CalculateGridsParam() //TO DO FOR QUATERNIONS
int n=0;
//alpha and gamma are uniform in -PI,PI
grid_alpha=2*M_PI/float(angleGridPointsAlpha);
grid_alpha=2.f * M_PI / (myfloat_t) angleGridPointsAlpha;
//cosine beta is uniform in -1,1
cos_grid_beta=2.0/float(angleGridPointsBeta);
cos_grid_beta=2.f / (myfloat_t) angleGridPointsBeta;
// Euler Angle Array
for (int ialpha = 0; ialpha < angleGridPointsAlpha; ialpha ++)
{
for (int ibeta = 0; ibeta < angleGridPointsBeta; ibeta ++)
{
for (int igamma = 0; igamma < angleGridPointsAlpha; igamma ++)
{
angles[n].pos[0]= float(ialpha)*grid_alpha-M_PI+grid_alpha*0.5; //ALPHA centered in the middle
angles[n].pos[1]= acos(float(ibeta)*cos_grid_beta-1+cos_grid_beta*0.5); //BETA centered in the middle
angles[n].pos[2]= float(igamma)*grid_alpha-M_PI+grid_alpha*0.5; //GAMMA centered in the middle
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
n++;
}
}
......@@ -243,14 +242,14 @@ int bioem_param::CalculateGridsParam() //TO DO FOR QUATERNIONS
/********** Calculating normalized volumen element *********/
param_device.volu=grid_alpha*grid_alpha*cos_grid_beta*float(param_device.GridSpaceCenter)*pixelSize*float(param_device.GridSpaceCenter)*pixelSize
*gridCTF_phase*gridCTF_amp*gridEnvelop/(2*M_PI)/(2*M_PI)/2/(2*float(param_device.maxDisplaceCenter))/(2*float(param_device.maxDisplaceCenter))/(float(numberGridPointsCTF_amp)*gridCTF_amp+startGridCTF_amp)
/(float(numberGridPointsCTF_phase)*gridCTF_phase+startGridCTF_phase)/(float(numberGridPointsEnvelop)*gridEnvelop+startGridEnvelop);
param_device.volu = grid_alpha * grid_alpha * cos_grid_beta * (myfloat_t) param_device.GridSpaceCenter * pixelSize * (myfloat_t) param_device.GridSpaceCenter * pixelSize
* 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 + startGridCTF_amp)
/ ((myfloat_t) numberGridPointsCTF_phase * gridCTF_phase + startGridCTF_phase) / ((myfloat_t) numberGridPointsEnvelop * gridEnvelop + startGridEnvelop);
/*** Number of total pixels***/
param_device.Ntotpi=float(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);
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);
return(0);
......@@ -265,7 +264,7 @@ int bioem_param::CalculateRefCTF()
myfloat_t amp,env,phase,ctf,radsq;
mycomplex_t* localCTF;
int nctfmax=int(float(param_device.NumberPixels)/2.0);
int nctfmax= param_device.NumberPixels / 2;
int n=0;
localCTF= (mycomplex_t *) myfftw_malloc(sizeof(mycomplex_t) * param_device.NumberPixels*param_device.NumberPixels);
......@@ -274,15 +273,15 @@ int bioem_param::CalculateRefCTF()
for (int iamp = 0; iamp < numberGridPointsCTF_amp ; iamp++) //Loop over amplitud
{
amp=float(iamp)*gridCTF_amp + startGridCTF_amp;
amp = (myfloat_t) iamp * gridCTF_amp + startGridCTF_amp;
for (int iphase = 0; iphase < numberGridPointsCTF_phase ; iphase++)//Loop over phase
{
phase=float(iphase)*gridCTF_phase + startGridCTF_phase;
phase = (myfloat_t) iphase * gridCTF_phase + startGridCTF_phase;
for ( int ienv = 0; ienv < numberGridPointsEnvelop ; ienv++)//Loop over envelope
for (int ienv = 0; ienv < numberGridPointsEnvelop ; ienv++)//Loop over envelope
{
env=float(ienv)*gridEnvelop + startGridEnvelop;
env= (myfloat_t) ienv * gridEnvelop + startGridEnvelop;
memset(localCTF,0,param_device.NumberPixels*param_device.NumberPixels*sizeof(mycomplex_t));
......@@ -291,7 +290,7 @@ int bioem_param::CalculateRefCTF()
{
for(int j=0; j< nctfmax; j++)
{
radsq=float(i*i+j*j)*pixelSize*pixelSize;
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));
localCTF[i*param_device.NumberPixels+j][0]=(myfloat_t) ctf;
......
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