ProSHADE  0.7.5.4 (MAR 2021)
Protein Shape Detection
ProSHADE_internal_overlay Namespace Reference

This namespace contains all the functions required for map overlays. More...

Functions

void getOptimalRotation (ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *staticStructure, ProSHADE_internal_data::ProSHADE_data *movingStructure, proshade_double *eulA, proshade_double *eulB, proshade_double *eulG)
 This function finds the optimal rotation between two structures as described by the settings object. More...
 
void getOptimalTranslation (ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *staticStructure, ProSHADE_internal_data::ProSHADE_data *movingStructure, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ, proshade_double eulA, proshade_double eulB, proshade_double eulG)
 This function finds the optimal translation between two structures as described by the settings object given a rotation between the two objects. More...
 
void computeBeforeAfterZeroCounts (proshade_unsign *addXPre, proshade_unsign *addYPre, proshade_unsign *addZPre, proshade_unsign *addXPost, proshade_unsign *addYPost, proshade_unsign *addZPost, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_unsign xDimIndices, proshade_unsign yDimIndices, proshade_unsign zDimIndices)
 This function finds the number of zeroes to be added after and before the structure along each dimension. More...
 
void paddMapWithZeroes (proshade_double *oldMap, proshade_double *&newMap, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_unsign xDimIndices, proshade_unsign yDimIndices, proshade_unsign zDimIndices, proshade_unsign addXPre, proshade_unsign addYPre, proshade_unsign addZPre)
 This function adds zeroes before and after the central map and copies the central map values into a new map. More...
 
void allocateTranslationFunctionMemory (fftw_complex *&tmpIn1, fftw_complex *&tmpOut1, fftw_complex *&tmpIn2, fftw_complex *&tmpOut2, fftw_complex *&resIn, fftw_complex *&resOut, fftw_plan &forwardFourierObj1, fftw_plan &forwardFourierObj2, fftw_plan &inverseFourierCombo, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD)
 This function allocates the memory for the Fourier transforms required for translation function computation. More...
 
void combineFourierForTranslation (fftw_complex *tmpOut1, fftw_complex *tmpOut2, fftw_complex *&resOut, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD)
 This function combines Fourier coefficients of two structures in a way, so that inverse Fourier of the combination will be the translation function. More...
 
void findHighestValueInMap (fftw_complex *resIn, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ, proshade_double *mapPeak)
 This function simply finds the highest value in fftw_complex map and returns its position and value. More...
 
void freeTranslationFunctionMemory (fftw_complex *&tmpIn1, fftw_complex *&tmpOut1, fftw_complex *&tmpIn2, fftw_complex *&tmpOut2, fftw_complex *&resOut, fftw_plan &forwardFourierObj1, fftw_plan &forwardFourierObj2, fftw_plan &inverseFourierCombo)
 This function releases the memory for the Fourier transforms required for translation function computation. More...
 
void computeAngularThreshold (std::vector< proshade_double > *lonCO, std::vector< proshade_double > *latCO, proshade_unsign angRes)
 This function computes the angular thresholds for longitude and lattitude angles. More...
 
void initialiseInverseSHComputation (proshade_unsign shBand, double *&sigR, double *&sigI, double *&rcoeffs, double *&icoeffs, double *&weights, double *&workspace, fftw_plan &idctPlan, fftw_plan &ifftPlan)
 This function initialises internal variables for inverse Spherical Harmonics computation. More...
 

Detailed Description

This namespace contains all the functions required for map overlays.

The map overlay task does have some specific functions that are betters groupped together for higher clarity of the code and also simple modifications later. This is where these live.

Function Documentation

◆ allocateTranslationFunctionMemory()

void ProSHADE_internal_overlay::allocateTranslationFunctionMemory ( fftw_complex *&  tmpIn1,
fftw_complex *&  tmpOut1,
fftw_complex *&  tmpIn2,
fftw_complex *&  tmpOut2,
fftw_complex *&  resIn,
fftw_complex *&  resOut,
fftw_plan &  forwardFourierObj1,
fftw_plan &  forwardFourierObj2,
fftw_plan &  inverseFourierCombo,
proshade_unsign  xD,
proshade_unsign  yD,
proshade_unsign  zD 
)

This function allocates the memory for the Fourier transforms required for translation function computation.

Parameters
[in]tmpIn1Array to hold the static structure Fourier inputs.
[in]tmpOut1Array to hold the static structure Fourier outputs.
[in]tmpIn2Array to hold the moving structure Fourier inputs.
[in]tmpOut2Array to hold the moving structure Fourier outputs.
[in]resInArray to hold the final translation function values.
[in]resOutArray to hold the combined Fourier coefficients of both structures.
[in]forwardFourierObj1FFTW plan for the forward Fourier of the static structure.
[in]forwardFourierObj2FFTW plan for the forward Fourier of the moving structure.
[in]inverseFourierComboFFTW plan for the backward Fourier of the combined Fourier factors of both structures.
[in]xDThe dimension of the X axis of the structures (assumes both structures have the same sizes and sampling).
[in]yDThe dimension of the Y axis of the structures (assumes both structures have the same sizes and sampling).
[in]zDThe dimension of the Z axis of the structures (assumes both structures have the same sizes and sampling).

Definition at line 432 of file ProSHADE_overlay.cpp.

433 {
434  //================================================ Allocate memory
435  tmpIn1 = new fftw_complex[xD * yD * zD];
436  tmpOut1 = new fftw_complex[xD * yD * zD];
437  tmpIn2 = new fftw_complex[xD * yD * zD];
438  tmpOut2 = new fftw_complex[xD * yD * zD];
439  resIn = new fftw_complex[xD * yD * zD];
440  resOut = new fftw_complex[xD * yD * zD];
441 
442  //================================================ Check memory allocation
443  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn1 , __FILE__, __LINE__, __func__ );
444  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut1, __FILE__, __LINE__, __func__ );
445  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn2, __FILE__, __LINE__, __func__ );
446  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut2, __FILE__, __LINE__, __func__ );
447  ProSHADE_internal_misc::checkMemoryAllocation ( resIn, __FILE__, __LINE__, __func__ );
448  ProSHADE_internal_misc::checkMemoryAllocation ( resOut, __FILE__, __LINE__, __func__ );
449 
450  //================================================ Get Fourier transforms of the maps
451  forwardFourierObj1 = fftw_plan_dft_3d ( xD, yD, zD, tmpIn1, tmpOut1, FFTW_FORWARD , FFTW_ESTIMATE );
452  forwardFourierObj2 = fftw_plan_dft_3d ( xD, yD, zD, tmpIn2, tmpOut2, FFTW_FORWARD , FFTW_ESTIMATE );
453  inverseFourierCombo = fftw_plan_dft_3d ( xD, yD, zD, resOut, resIn , FFTW_BACKWARD, FFTW_ESTIMATE );
454 
455  //================================================ Done
456  return ;
457 
458 }

◆ combineFourierForTranslation()

void ProSHADE_internal_overlay::combineFourierForTranslation ( fftw_complex *  tmpOut1,
fftw_complex *  tmpOut2,
fftw_complex *&  resOut,
proshade_unsign  xD,
proshade_unsign  yD,
proshade_unsign  zD 
)

This function combines Fourier coefficients of two structures in a way, so that inverse Fourier of the combination will be the translation function.

Parameters
[in]tmpOut1Array holding the static structure Fourier outputs.
[in]tmpOut2Array holding the moving structure Fourier outputs.
[in]resOutArray to hold the combined Fourier coefficients of both structures.
[in]xDThe dimension of the X axis of the structures (assumes both structures have the same sizes and sampling).
[in]yDThe dimension of the Y axis of the structures (assumes both structures have the same sizes and sampling).
[in]zDThe dimension of the Z axis of the structures (assumes both structures have the same sizes and sampling).

Definition at line 497 of file ProSHADE_overlay.cpp.

498 {
499  //================================================ Initialise local variables
500  double normFactor = static_cast<double> ( xD * yD * zD );
501  proshade_signed h1, k1, l1, hlpPos, arrPos;
502  proshade_signed uo2 = xD / 2;
503  proshade_signed vo2 = yD / 2;
504  proshade_signed wo2 = zD / 2;
505 
506  //================================================ Combine the coefficients
507  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> ( xD ); uIt++ )
508  {
509  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> ( yD ); vIt++ )
510  {
511  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> ( zD ); wIt++ )
512  {
513  //==================================== Convert to HKL
514  if ( uIt > uo2 ) { h1 = uIt - xD; } else { h1 = uIt; }
515  if ( vIt > vo2 ) { k1 = vIt - yD; } else { k1 = vIt; }
516  if ( wIt > wo2 ) { l1 = wIt - zD; } else { l1 = wIt; }
517 
518  //==================================== Make HKL into indexable numbers
519  if ( h1 < 0 ) { h1 += xD; }
520  if ( k1 < 0 ) { k1 += yD; }
521  if ( l1 < 0 ) { l1 += zD; }
522 
523  //==================================== Find indices
524  hlpPos = l1 + zD * ( k1 + yD * h1 );
525  arrPos = wIt + zD * ( vIt + yD * uIt );
526 
527  //==================================== Combine
529  &tmpOut1[hlpPos][1],
530  &tmpOut2[arrPos][0],
531  &tmpOut2[arrPos][1],
532  &resOut[hlpPos][0],
533  &resOut[hlpPos][1] );
534 
535  //==================================== Save
536  resOut[hlpPos][0] /= normFactor;
537  resOut[hlpPos][1] /= normFactor;
538  }
539  }
540  }
541 
542  //======================================== Done
543  return ;
544 
545 }

◆ computeAngularThreshold()

void ProSHADE_internal_overlay::computeAngularThreshold ( std::vector< proshade_double > *  lonCO,
std::vector< proshade_double > *  latCO,
proshade_unsign  angRes 
)

This function computes the angular thresholds for longitude and lattitude angles.

Parameters
[in]lonCOPointer to vector where longitude thresholds are to be stored at.
[in]latCOPointer to vector where lattitude thresholds are to be stored at.
[in]angResThe angular resolution of the computation.

Definition at line 1018 of file ProSHADE_overlay.cpp.

1019 {
1020  //================================================ Longitude angular thresholds
1021  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1022  {
1023  ProSHADE_internal_misc::addToDoubleVector ( lonCO, static_cast<proshade_double> ( iter ) * ( ( static_cast<proshade_double> ( M_PI ) * 2.0 ) / static_cast<proshade_double> ( angRes ) ) - ( static_cast<double> ( M_PI ) ) );
1024  }
1025 
1026  //================================================ Lattitude angular thresholds
1027  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1028  {
1029  ProSHADE_internal_misc::addToDoubleVector ( latCO, static_cast<proshade_double> ( iter ) * ( static_cast<proshade_double> ( M_PI ) / static_cast<proshade_double> ( angRes ) ) - ( static_cast<proshade_double> ( M_PI ) / 2.0 ) );
1030  }
1031 
1032  //================================================ Done
1033  return ;
1034 
1035 }

◆ computeBeforeAfterZeroCounts()

void ProSHADE_internal_overlay::computeBeforeAfterZeroCounts ( proshade_unsign *  addXPre,
proshade_unsign *  addYPre,
proshade_unsign *  addZPre,
proshade_unsign *  addXPost,
proshade_unsign *  addYPost,
proshade_unsign *  addZPost,
proshade_unsign  xDim,
proshade_unsign  yDim,
proshade_unsign  zDim,
proshade_unsign  xDimIndices,
proshade_unsign  yDimIndices,
proshade_unsign  zDimIndices 
)

This function finds the number of zeroes to be added after and before the structure along each dimension.

Parameters
[in]addXPreVariable pointer where the number of zeroes to be added before the map along X axis is saved.
[in]addYPreVariable pointer where the number of zeroes to be added before the map along Y axis is saved.
[in]addZPreVariable pointer where the number of zeroes to be added before the map along Z axis is saved.
[in]addXPostVariable pointer where the number of zeroes to be added after the map along X axis is saved.
[in]addYPostVariable pointer where the number of zeroes to be added after the map along Y axis is saved.
[in]addZPostVariable pointer where the number of zeroes to be added after the map along Z axis is saved.
[in]xDimThe X dimension size in indices of the new map.
[in]yDimThe Y dimension size in indices of the new map.
[in]zDimThe Z dimension size in indices of the new map.
[in]xDimIndicesThe X dimension size in indices of the old map.
[in]yDimIndicesThe Y dimension size in indices of the old map.
[in]zDimIndicesThe Z dimension size in indices of the old map.

Definition at line 657 of file ProSHADE_overlay.cpp.

658 {
659  //================================================ Compute
660  *addXPre = ( xDim - xDimIndices ) / 2;
661  *addYPre = ( yDim - yDimIndices ) / 2;
662  *addZPre = ( zDim - zDimIndices ) / 2;
663  *addXPost = *addXPre; if ( ( xDim - xDimIndices ) % 2 == 1 ) { *addXPost += 1; }
664  *addYPost = *addYPre; if ( ( yDim - yDimIndices ) % 2 == 1 ) { *addYPost += 1; }
665  *addZPost = *addZPre; if ( ( zDim - zDimIndices ) % 2 == 1 ) { *addZPost += 1; }
666 
667  //================================================ Done
668  return ;
669 
670 }

◆ findHighestValueInMap()

void ProSHADE_internal_overlay::findHighestValueInMap ( fftw_complex *  resIn,
proshade_unsign  xD,
proshade_unsign  yD,
proshade_unsign  zD,
proshade_double *  trsX,
proshade_double *  trsY,
proshade_double *  trsZ,
proshade_double *  mapPeak 
)

This function simply finds the highest value in fftw_complex map and returns its position and value.

Parameters
[in]resInArray holding the translation function values.
[in]xDThe dimension of the X axis of the structures (assumes both structures have the same sizes and sampling).
[in]yDThe dimension of the Y axis of the structures (assumes both structures have the same sizes and sampling).
[in]zDThe dimension of the Z axis of the structures (assumes both structures have the same sizes and sampling).
[in]trsXVariable to which the X axis translation function peak position will be saved to.
[in]trsYVariable to which the Y axis translation function peak position will be saved to.
[in]trsZVariable to which the Z axis translation function peak position will be saved to.
[in]mapPeakVariable to which the height of the translation function peak will be saved to.

Definition at line 558 of file ProSHADE_overlay.cpp.

559 {
560  //================================================ Initialise variables
561  proshade_signed arrPos;
562  *mapPeak = 0.0;
563 
564  //================================================ Search the map
565  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> ( xD ); uIt++ )
566  {
567  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> ( yD ); vIt++ )
568  {
569  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> ( zD ); wIt++ )
570  {
571  arrPos = wIt + zD * ( vIt + yD * uIt );
572  if ( resIn[arrPos][0] > *mapPeak )
573  {
574  *mapPeak = resIn[arrPos][0];
575  *trsX = uIt;
576  *trsY = vIt;
577  *trsZ = wIt;
578  }
579  }
580  }
581  }
582 
583  //================================================ Done
584  return ;
585 
586 }

◆ freeTranslationFunctionMemory()

void ProSHADE_internal_overlay::freeTranslationFunctionMemory ( fftw_complex *&  tmpIn1,
fftw_complex *&  tmpOut1,
fftw_complex *&  tmpIn2,
fftw_complex *&  tmpOut2,
fftw_complex *&  resOut,
fftw_plan &  forwardFourierObj1,
fftw_plan &  forwardFourierObj2,
fftw_plan &  inverseFourierCombo 
)

This function releases the memory for the Fourier transforms required for translation function computation.

Parameters
[in]tmpIn1Array to hold the static structure Fourier inputs.
[in]tmpOut1Array to hold the static structure Fourier outputs.
[in]tmpIn2Array to hold the moving structure Fourier inputs.
[in]tmpOut2Array to hold the moving structure Fourier outputs.
[in]resOutArray to hold the combined Fourier coefficients of both structures.
[in]forwardFourierObj1FFTW plan for the forward Fourier of the static structure.
[in]forwardFourierObj2FFTW plan for the forward Fourier of the moving structure.
[in]inverseFourierComboFFTW plan for the backward Fourier of the combined Fourier factors of both structures.

Definition at line 471 of file ProSHADE_overlay.cpp.

472 {
473  //================================================ Release memory
474  fftw_destroy_plan ( forwardFourierObj1 );
475  fftw_destroy_plan ( forwardFourierObj2 );
476  fftw_destroy_plan ( inverseFourierCombo );
477  delete[] tmpIn1;
478  delete[] tmpIn2;
479  delete[] tmpOut1;
480  delete[] tmpOut2;
481  delete[] resOut;
482 
483  //======================================== Done
484  return ;
485 
486 }

◆ getOptimalRotation()

void ProSHADE_internal_overlay::getOptimalRotation ( ProSHADE_settings settings,
ProSHADE_internal_data::ProSHADE_data staticStructure,
ProSHADE_internal_data::ProSHADE_data movingStructure,
proshade_double *  eulA,
proshade_double *  eulB,
proshade_double *  eulG 
)

This function finds the optimal rotation between two structures as described by the settings object.

This function takes the settings and two structure classes. It then reads in and processes both structures so that the globally optimal rotation overlay Euler angles are detected. However, this is only the case for rotation along the centre of the map of the second structure; therefore, either use Patterson data (usePhase = false), or be aware that better rotation may exist for different centre of rotation.

Parameters
[in]settingsA pointer to settings class containing all the information required for map symmetry detection.
[in]obj1A pointer to the data class object of the other ( static ) structure.
[in]obj2A pointer to the data class object of the first ( moving ) structure.
[in]eulAThe variable to which the best Euler alpha angle value will be saved to.
[in]eulBThe variable to which the best Euler beta angle value will be saved to.
[in]eulGThe variable to which the best Euler gamma angle value will be saved to.

Definition at line 154 of file ProSHADE_overlay.cpp.

155 {
156  //================================================ Read in the structures
157  staticStructure->readInStructure ( settings->inputFiles.at(0), 0, settings );
158  movingStructure->readInStructure ( settings->inputFiles.at(1), 1, settings );
159 
160  //================================================ Internal data processing (COM, norm, mask, extra space)
161  staticStructure->processInternalMap ( settings );
162  movingStructure->processInternalMap ( settings );
163 
164  //================================================ Map to sphere
165  staticStructure->mapToSpheres ( settings );
166  movingStructure->mapToSpheres ( settings );
167 
168  //================================================ Get spherical harmonics
169  staticStructure->computeSphericalHarmonics ( settings );
170  movingStructure->computeSphericalHarmonics ( settings );
171 
172  //================================================ Get the rotation function of the pair
173  movingStructure->getOverlayRotationFunction ( settings, staticStructure );
174 
175  //================================================ Get inverse SO(3) map top peak Euler angle values
177  std::min ( staticStructure->getMaxBand(), movingStructure->getMaxBand() ) * 2,
178  eulA, eulB, eulG, settings );
179 
180  //================================================ Done
181  return ;
182 
183 }

◆ getOptimalTranslation()

void ProSHADE_internal_overlay::getOptimalTranslation ( ProSHADE_settings settings,
ProSHADE_internal_data::ProSHADE_data staticStructure,
ProSHADE_internal_data::ProSHADE_data movingStructure,
proshade_double *  trsX,
proshade_double *  trsY,
proshade_double *  trsZ,
proshade_double  eulA,
proshade_double  eulB,
proshade_double  eulG 
)

This function finds the optimal translation between two structures as described by the settings object given a rotation between the two objects.

This function starts by loading and processing the structures according to the settings object (keeping phase is assumed, but callers responsibility). It then applies the required rotation to the second (moing) strucutre and then it follows with zero padding to make sure the structures have the same dimensions (again, it assumes map re-sampling was done, but setting to is callers responsibility). It then computes the translation function, finds the highest peak and returns the positions as well as height of this peak.

Parameters
[in]settingsA pointer to settings class containing all the information required for map overlay computation.
[in]staticStructureA pointer to the data class object of the other ( static ) structure.
[in]movingStructureA pointer to the data class object of the first ( moving ) structure.
[in]trsXThe variable to which the best X-axis position value will be saved to.
[in]trsYThe variable to which the best Y-axis position value will be saved to.
[in]trsZThe variable to which the best Z-axis position value will be saved to.
[in]eulAThe Euler alpha angle value, by which the moving structure is to be rotated by.
[in]eulBThe Euler beta angle value, by which the moving structure is to be rotated by.
[in]eulGThe Euler gamma angle value, by which the moving structure is to be rotated by.

Definition at line 203 of file ProSHADE_overlay.cpp.

204 {
205  //================================================ Read in the structures
206  staticStructure->readInStructure ( settings->inputFiles.at(0), 0, settings );
207  movingStructure->readInStructure ( settings->inputFiles.at(1), 1, settings );
208 
209  //================================================ Internal data processing (COM, norm, mask, extra space)
210  staticStructure->processInternalMap ( settings );
211  movingStructure->processInternalMap ( settings );
212 
213  //================================================ Compute spherical harmonics to allow Fourier space rotation
214  movingStructure->mapToSpheres ( settings );
215  movingStructure->computeSphericalHarmonics ( settings );
216 
217  //================================================ Rotate map
218  movingStructure->rotateMap ( settings, eulA, eulB, eulG );
219 
220  //================================================ Zero padding for smaller structure
221  staticStructure->zeroPaddToDims ( std::max ( staticStructure->getXDim(), movingStructure->getXDim() ),
222  std::max ( staticStructure->getYDim(), movingStructure->getYDim() ),
223  std::max ( staticStructure->getZDim(), movingStructure->getZDim() ) );
224  movingStructure->zeroPaddToDims ( std::max ( staticStructure->getXDim(), movingStructure->getXDim() ),
225  std::max ( staticStructure->getYDim(), movingStructure->getYDim() ),
226  std::max ( staticStructure->getZDim(), movingStructure->getZDim() ) );
227 
228  //================================================ Report progress
229  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting translation function computation." );
230 
231  //================================================ Compute the translation function map
232  movingStructure->computeTranslationMap ( staticStructure );
233 
234  //================================================ Find highest peak in translation function
235  proshade_double mapPeak = 0.0;
236  proshade_unsign xDimS = staticStructure->getXDim();
237  proshade_unsign yDimS = staticStructure->getYDim();
238  proshade_unsign zDimS = staticStructure->getZDim();
239  ProSHADE_internal_overlay::findHighestValueInMap ( movingStructure->getTranslationFnPointer(), xDimS, yDimS, zDimS, trsX, trsY, trsZ, &mapPeak );
240 
241  //================================================ Dont translate over half
242  if ( *trsX > ( xDimS / 2 ) ) { *trsX = *trsX - xDimS; }
243  if ( *trsY > ( yDimS / 2 ) ) { *trsY = *trsY - yDimS; }
244  if ( *trsZ > ( zDimS / 2 ) ) { *trsZ = *trsZ - zDimS; }
245 
246  //================================================ Move map
247  proshade_single xCor = ( staticStructure->xFrom - movingStructure->xFrom ) *
248  ( static_cast<proshade_double> ( staticStructure->getXDimSize() ) / staticStructure->getXDim() );
249  proshade_single yCor = ( staticStructure->yFrom - movingStructure->yFrom ) *
250  ( static_cast<proshade_double> ( staticStructure->getYDimSize() ) / staticStructure->getYDim() );
251  proshade_single zCor = ( staticStructure->zFrom - movingStructure->zFrom ) *
252  ( static_cast<proshade_double> ( staticStructure->getZDimSize() ) / staticStructure->getZDim() );
253  proshade_single xMov = staticStructure->mapCOMProcessChangeX - xCor -
254  ( *trsX * static_cast<proshade_double> ( staticStructure->getXDimSize() ) / staticStructure->getXDim() );
255  proshade_single yMov = staticStructure->mapCOMProcessChangeY - yCor -
256  ( *trsY * static_cast<proshade_double> ( staticStructure->getYDimSize() ) / staticStructure->getYDim() );
257  proshade_single zMov = staticStructure->mapCOMProcessChangeZ - zCor -
258  ( *trsZ * static_cast<proshade_double> ( staticStructure->getZDimSize() ) / staticStructure->getZDim() );
259 
260  //================================================ Save translation vector back
261  *trsX = -xMov;
262  *trsY = -yMov;
263  *trsZ = -zMov;
264 
265  //================================================ Report progress
266  std::stringstream hlpSS;
267  hlpSS << "Optimal map translation distances are " << *trsX << " ; " << *trsY << " ; " << *trsZ << " Angstroms with peak height " << mapPeak / ( xDimS * yDimS * zDimS );
268 
269  //================================================ Save original from variables for PDB output
270  movingStructure->mapMovFromsChangeX = movingStructure->xFrom;
271  movingStructure->mapMovFromsChangeY = movingStructure->yFrom;
272  movingStructure->mapMovFromsChangeZ = movingStructure->zFrom;
273 
274  //================================================ Move the map
275  ProSHADE_internal_mapManip::moveMapByIndices ( &xMov, &yMov, &zMov, movingStructure->getXDimSize(), movingStructure->getYDimSize(), movingStructure->getZDimSize(),
276  movingStructure->getXFromPtr(), movingStructure->getXToPtr(),
277  movingStructure->getYFromPtr(), movingStructure->getYToPtr(),
278  movingStructure->getZFromPtr(), movingStructure->getZToPtr(),
279  movingStructure->getXAxisOrigin(), movingStructure->getYAxisOrigin(), movingStructure->getZAxisOrigin() );
280 
281  ProSHADE_internal_mapManip::moveMapByFourier ( movingStructure->getInternalMap(), xMov, yMov, zMov,
282  movingStructure->getXDimSize(), movingStructure->getYDimSize(), movingStructure->getZDimSize(),
283  movingStructure->getXDim(), movingStructure->getYDim(), movingStructure->getZDim() );
284 
285  //================================================ Report progress
286  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, hlpSS.str() );
287  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Translation function computation complete." );
288 
289  //================================================ Keep only the change in from and to variables
290  movingStructure->mapMovFromsChangeX = movingStructure->xFrom - movingStructure->mapMovFromsChangeX;
291  movingStructure->mapMovFromsChangeY = movingStructure->yFrom - movingStructure->mapMovFromsChangeY;
292  movingStructure->mapMovFromsChangeZ = movingStructure->zFrom - movingStructure->mapMovFromsChangeZ;
293 
294  //================================================ Compute the optimal rotation centre for co-ordinates
295  movingStructure->computePdbRotationCentre ( );
296  movingStructure->computeOptimalTranslation ( eulA, eulB, eulG, *trsX, *trsY, *trsZ );
297 
298  //================================================ Done
299  return ;
300 
301 }

◆ initialiseInverseSHComputation()

void ProSHADE_internal_overlay::initialiseInverseSHComputation ( proshade_unsign  shBand,
double *&  sigR,
double *&  sigI,
double *&  rcoeffs,
double *&  icoeffs,
double *&  weights,
double *&  workspace,
fftw_plan &  idctPlan,
fftw_plan &  ifftPlan 
)

This function initialises internal variables for inverse Spherical Harmonics computation.

Parameters
[in]shBandThe bandwidth for this particular shell.
[in]sigRPointer to be initialised for the real signal values.
[in]sigIPointer to be initialised for the imaginary signal values.
[in]rcoeffsPointer to be initialised for the real coefficient values.
[in]icoeffsPointer to be initialised for the imaginary coefficient values.
[in]weightsPointer to be initialised for the transform weight values.
[in]workspacePointer to be initialised for the computation screatch space.
[in]idctPlanPointer reference to the cosine/sine transform plan to be created.
[in]ifftPlanPointer reference to the discrete 3D Fourier transform plan to be created.

Definition at line 902 of file ProSHADE_overlay.cpp.

903 {
904  //================================================ Initialise internal variables
905  proshade_unsign oneDim = shBand * 2;
906 
907  //================================================ Allocate memory
908  sigR = new double [(oneDim * oneDim * 4)];
909  sigI = sigR + (oneDim * oneDim * 2);
910  rcoeffs = new double [(oneDim * oneDim * 2)];
911  icoeffs = rcoeffs + (oneDim * oneDim);
912  weights = new double [4 * shBand];
913  workspace = new double [( 20 * shBand * shBand ) + ( 24 * shBand )];
914 
915  //================================================ Check memory allocation
916  ProSHADE_internal_misc::checkMemoryAllocation ( sigR, __FILE__, __LINE__, __func__ );
917  ProSHADE_internal_misc::checkMemoryAllocation ( rcoeffs, __FILE__, __LINE__, __func__ );
918  ProSHADE_internal_misc::checkMemoryAllocation ( weights, __FILE__, __LINE__, __func__ );
919  ProSHADE_internal_misc::checkMemoryAllocation ( workspace, __FILE__, __LINE__, __func__ );
920 
921  //================================================ Create the cosine/sine transform plan
922  idctPlan = fftw_plan_r2r_1d ( oneDim, weights, workspace, FFTW_REDFT01, FFTW_ESTIMATE );
923 
924  //================================================ Set up the discrete Fourier transform
925  int rank, howmany_rank;
926  fftw_iodim dims[1], howmany_dims[1];
927 
928  rank = 1;
929  dims[0].n = 2 * shBand;
930  dims[0].is = 2 * shBand;
931  dims[0].os = 1;
932  howmany_rank = 1;
933  howmany_dims[0].n = 2 * shBand;
934  howmany_dims[0].is = 1;
935  howmany_dims[0].os = 2 * shBand;
936 
937  //================================================ Create the discrete Fourier transform
938  ifftPlan = fftw_plan_guru_split_dft ( rank, dims, howmany_rank, howmany_dims, sigR, sigI, rcoeffs, icoeffs, FFTW_ESTIMATE );
939 
940  //================================================ Done
941  return ;
942 
943 }

◆ paddMapWithZeroes()

void ProSHADE_internal_overlay::paddMapWithZeroes ( proshade_double *  oldMap,
proshade_double *&  newMap,
proshade_unsign  xDim,
proshade_unsign  yDim,
proshade_unsign  zDim,
proshade_unsign  xDimIndices,
proshade_unsign  yDimIndices,
proshade_unsign  zDimIndices,
proshade_unsign  addXPre,
proshade_unsign  addYPre,
proshade_unsign  addZPre 
)

This function adds zeroes before and after the central map and copies the central map values into a new map.

Parameters
[in]oldMapThe original, unpadded map.
[in]newMapThe array to which the new map will be saved into.
[in]xDimThe X dimension size of the new map.
[in]yDimThe Y dimension size of the new map.
[in]zDimThe Z dimension size of the new map.
[in]xDimIndicesThe X dimension size in indices of the old map.
[in]yDimIndicesThe Y dimension size in indices of the old map.
[in]zDimIndicesThe Z dimension size in indices of the old map.
[in]addXPreHow many zeroes are to be added before the central map along the X axis?
[in]addYPreHow many zeroes are to be added before the central map along the Y axis?
[in]addZPreHow many zeroes are to be added before the central map along the Z axis?

Definition at line 686 of file ProSHADE_overlay.cpp.

687 {
688  //================================================ Initialise local variables
689  proshade_unsign newMapIndex = 0;
690  proshade_unsign oldMapIndex = 0;
691 
692  //================================================ Copy and padd as appropriate
693  for ( proshade_unsign xIt = 0; xIt < xDim; xIt++ )
694  {
695  for ( proshade_unsign yIt = 0; yIt < yDim; yIt++ )
696  {
697  for ( proshade_unsign zIt = 0; zIt < zDim; zIt++ )
698  {
699  //==================================== Find map location
700  newMapIndex = zIt + zDim * ( yIt + yDim * xIt );
701 
702  //==================================== If less than needed, add zero
703  if ( xIt < addXPre ) { newMap[newMapIndex] = 0.0; continue; }
704  if ( yIt < addYPre ) { newMap[newMapIndex] = 0.0; continue; }
705  if ( zIt < addZPre ) { newMap[newMapIndex] = 0.0; continue; }
706 
707  //==================================== If more than needed, add zero
708  if ( xIt >= ( addXPre + xDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
709  if ( yIt >= ( addYPre + yDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
710  if ( zIt >= ( addZPre + zDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
711 
712  //==================================== If not padding, copy original values
713  oldMapIndex = (zIt-addZPre) + zDimIndices * ( (yIt-addYPre) + yDimIndices * (xIt-addXPre) );
714  newMap[newMapIndex] = oldMap[oldMapIndex];
715  }
716  }
717  }
718 
719  //======================================== Done
720  return ;
721 
722 }
ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics
void computeSphericalHarmonics(ProSHADE_settings *settings)
This function computes the spherical harmonics decomposition for the whole structure.
Definition: ProSHADE_data.cpp:1608
ProSHADE_internal_data::ProSHADE_data::mapMovFromsChangeX
proshade_double mapMovFromsChangeX
When the map is translated, the xFrom and xTo values are changed. This variable holds how much they h...
Definition: ProSHADE_data.hpp:94
ProSHADE_internal_data::ProSHADE_data::getOverlayRotationFunction
void getOverlayRotationFunction(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *obj2)
This function computes the overlay rotation function (i.e. the correlation function in SO(3) space).
Definition: ProSHADE_overlay.cpp:35
ProSHADE_internal_data::ProSHADE_data::zFrom
proshade_signed zFrom
This is the starting index along the z axis.
Definition: ProSHADE_data.hpp:112
ProSHADE_internal_data::ProSHADE_data::getZDimSize
proshade_single getZDimSize(void)
This function allows access to the map size in angstroms along the Z axis.
Definition: ProSHADE_data.cpp:3303
ProSHADE_internal_data::ProSHADE_data::computeOptimalTranslation
void computeOptimalTranslation(proshade_double eulA, proshade_double eulB, proshade_double eulG, proshade_double trsX, proshade_double trsY, proshade_double trsZ)
This function computes and saves the optimal translation vector from the already determined translati...
Definition: ProSHADE_overlay.cpp:109
ProSHADE_internal_data::ProSHADE_data::getXAxisOrigin
proshade_signed * getXAxisOrigin(void)
This function allows access to the map X axis origin value.
Definition: ProSHADE_data.cpp:3403
ProSHADE_internal_data::ProSHADE_data::getXDimSize
proshade_single getXDimSize(void)
This function allows access to the map size in angstroms along the X axis.
Definition: ProSHADE_data.cpp:3283
ProSHADE_internal_data::ProSHADE_data::getInternalMap
proshade_double *& getInternalMap(void)
This function allows access to the first map array value address.
Definition: ProSHADE_data.cpp:3433
ProSHADE_internal_data::ProSHADE_data::getYDimSize
proshade_single getYDimSize(void)
This function allows access to the map size in angstroms along the Y axis.
Definition: ProSHADE_data.cpp:3293
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeX
proshade_double mapCOMProcessChangeX
The change in X axis between the creation of the structure (originalMapXCom) and just before rotation...
Definition: ProSHADE_data.hpp:97
ProSHADE_internal_mapManip::moveMapByFourier
void moveMapByFourier(proshade_double *&map, proshade_single xMov, proshade_single yMov, proshade_single zMov, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
Function for moving map back to original PDB location by using Fourier transformation.
Definition: ProSHADE_mapManip.cpp:792
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeZ
proshade_double mapCOMProcessChangeZ
The change in Z axis between the creation of the structure (originalMapZCom) and just before rotation...
Definition: ProSHADE_data.hpp:99
ProSHADE_internal_data::ProSHADE_data::computePdbRotationCentre
void computePdbRotationCentre(void)
This function computes the optimal rotation centre for co-ordinates.
Definition: ProSHADE_overlay.cpp:68
ProSHADE_internal_data::ProSHADE_data::getXFromPtr
proshade_signed * getXFromPtr(void)
This function allows access to the map start along the X axis.
Definition: ProSHADE_data.cpp:3343
ProSHADE_settings::verbose
proshade_signed verbose
Should the software report on the progress, or just be quiet? Value between -1 (nothing) and 4 (loud)
Definition: ProSHADE_settings.hpp:188
ProSHADE_internal_misc::addToDoubleVector
void addToDoubleVector(std::vector< proshade_double > *vecToAddTo, proshade_double elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:77
ProSHADE_internal_data::ProSHADE_data::getXDim
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
Definition: ProSHADE_data.cpp:3313
ProSHADE_internal_data::ProSHADE_data::rotateMap
void rotateMap(ProSHADE_settings *settings, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function rotates a map based on the given Euler angles.
Definition: ProSHADE_overlay.cpp:736
ProSHADE_internal_data::ProSHADE_data::mapToSpheres
void mapToSpheres(ProSHADE_settings *settings)
This function converts the internal map onto a set of concentric spheres.
Definition: ProSHADE_data.cpp:1564
ProSHADE_internal_data::ProSHADE_data::mapMovFromsChangeY
proshade_double mapMovFromsChangeY
When the map is translated, the yFrom and yTo values are changed. This variable holds how much they h...
Definition: ProSHADE_data.hpp:95
ProSHADE_internal_peakSearch::getBestPeakEulerAngsNaive
void getBestPeakEulerAngsNaive(proshade_complex *map, proshade_unsign dim, proshade_double *eulA, proshade_double *eulB, proshade_double *eulG, ProSHADE_settings *settings)
This function finds the highest peaks optimised Euler angles using the "naive" approach.
Definition: ProSHADE_peakSearch.cpp:351
ProSHADE_internal_data::ProSHADE_data::computeTranslationMap
void computeTranslationMap(ProSHADE_internal_data::ProSHADE_data *obj1)
This function does the computation of the translation map and saves results internally.
Definition: ProSHADE_overlay.cpp:389
ProSHADE_internal_data::ProSHADE_data::xFrom
proshade_signed xFrom
This is the starting index along the x axis.
Definition: ProSHADE_data.hpp:110
ProSHADE_internal_data::ProSHADE_data::getMaxBand
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
Definition: ProSHADE_data.cpp:3004
ProSHADE_internal_data::ProSHADE_data::getYDim
proshade_unsign getYDim(void)
This function allows access to the map size in indices along the Y axis.
Definition: ProSHADE_data.cpp:3323
ProSHADE_internal_data::ProSHADE_data::processInternalMap
void processInternalMap(ProSHADE_settings *settings)
This function simply clusters several other functions which should be called together.
Definition: ProSHADE_data.cpp:1464
ProSHADE_internal_data::ProSHADE_data::getTranslationFnPointer
proshade_complex * getTranslationFnPointer(void)
This function allows access to the translation function through a pointer.
Definition: ProSHADE_data.cpp:3443
ProSHADE_internal_data::ProSHADE_data::getZDim
proshade_unsign getZDim(void)
This function allows access to the map size in indices along the Z axis.
Definition: ProSHADE_data.cpp:3333
ProSHADE_internal_overlay::findHighestValueInMap
void findHighestValueInMap(fftw_complex *resIn, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ, proshade_double *mapPeak)
This function simply finds the highest value in fftw_complex map and returns its position and value.
Definition: ProSHADE_overlay.cpp:558
ProSHADE_internal_data::ProSHADE_data::mapMovFromsChangeZ
proshade_double mapMovFromsChangeZ
When the map is translated, the zFrom and zTo values are changed. This variable holds how much they h...
Definition: ProSHADE_data.hpp:96
ProSHADE_internal_data::ProSHADE_data::getYToPtr
proshade_signed * getYToPtr(void)
This function allows access to the map last position along the Y axis.
Definition: ProSHADE_data.cpp:3383
ProSHADE_internal_data::ProSHADE_data::getZAxisOrigin
proshade_signed * getZAxisOrigin(void)
This function allows access to the map Z axis origin value.
Definition: ProSHADE_data.cpp:3423
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:65
ProSHADE_internal_data::ProSHADE_data::getXToPtr
proshade_signed * getXToPtr(void)
This function allows access to the map last position along the X axis.
Definition: ProSHADE_data.cpp:3373
ProSHADE_internal_data::ProSHADE_data::getYFromPtr
proshade_signed * getYFromPtr(void)
This function allows access to the map start along the Y axis.
Definition: ProSHADE_data.cpp:3353
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeY
proshade_double mapCOMProcessChangeY
The change in Y axis between the creation of the structure (originalMapYCom) and just before rotation...
Definition: ProSHADE_data.hpp:98
ProSHADE_internal_data::ProSHADE_data::getInvSO3Coeffs
proshade_complex * getInvSO3Coeffs(void)
This function allows access to the inverse SO(3) coefficients array.
Definition: ProSHADE_data.cpp:3231
ProSHADE_internal_data::ProSHADE_data::yFrom
proshade_signed yFrom
This is the starting index along the y axis.
Definition: ProSHADE_data.hpp:111
ProSHADE_settings::inputFiles
std::vector< std::string > inputFiles
This vector contains the filenames of all input structure files.
Definition: ProSHADE_settings.hpp:95
ProSHADE_internal_maths::complexMultiplicationConjug
void complexMultiplicationConjug(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2, proshade_double *retReal, proshade_double *retImag)
Function to multiply two complex numbers by using the second number's conjugate.
Definition: ProSHADE_maths.cpp:62
ProSHADE_internal_data::ProSHADE_data::getYAxisOrigin
proshade_signed * getYAxisOrigin(void)
This function allows access to the map Y axis origin value.
Definition: ProSHADE_data.cpp:3413
ProSHADE_internal_messages::printProgressMessage
void printProgressMessage(proshade_signed verbose, proshade_signed messageLevel, std::string message)
General stdout message printing.
Definition: ProSHADE_messages.cpp:70
ProSHADE_internal_mapManip::moveMapByIndices
void moveMapByIndices(proshade_single *xMov, proshade_single *yMov, proshade_single *zMov, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed *xFrom, proshade_signed *xTo, proshade_signed *yFrom, proshade_signed *yTo, proshade_signed *zFrom, proshade_signed *zTo, proshade_signed *xOrigin, proshade_signed *yOrigin, proshade_signed *zOrigin)
Function for moving map back to original PDB location by changing the indices.
Definition: ProSHADE_mapManip.cpp:744
ProSHADE_internal_data::ProSHADE_data::zeroPaddToDims
void zeroPaddToDims(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim)
This function changes the size of a structure to fit the supplied new limits.
Definition: ProSHADE_overlay.cpp:598
ProSHADE_internal_data::ProSHADE_data::readInStructure
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
Definition: ProSHADE_data.cpp:447
ProSHADE_internal_data::ProSHADE_data::getZFromPtr
proshade_signed * getZFromPtr(void)
This function allows access to the map start along the Z axis.
Definition: ProSHADE_data.cpp:3363
ProSHADE_internal_data::ProSHADE_data::getZToPtr
proshade_signed * getZToPtr(void)
This function allows access to the map last position along the Z axis.
Definition: ProSHADE_data.cpp:3393