ProSHADE  0.7.6.0 (JUL 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 441 of file ProSHADE_overlay.cpp.

442 {
443  //================================================ Allocate memory
444  tmpIn1 = new fftw_complex[xD * yD * zD];
445  tmpOut1 = new fftw_complex[xD * yD * zD];
446  tmpIn2 = new fftw_complex[xD * yD * zD];
447  tmpOut2 = new fftw_complex[xD * yD * zD];
448  resIn = new fftw_complex[xD * yD * zD];
449  resOut = new fftw_complex[xD * yD * zD];
450 
451  //================================================ Check memory allocation
452  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn1 , __FILE__, __LINE__, __func__ );
453  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut1, __FILE__, __LINE__, __func__ );
454  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn2, __FILE__, __LINE__, __func__ );
455  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut2, __FILE__, __LINE__, __func__ );
456  ProSHADE_internal_misc::checkMemoryAllocation ( resIn, __FILE__, __LINE__, __func__ );
457  ProSHADE_internal_misc::checkMemoryAllocation ( resOut, __FILE__, __LINE__, __func__ );
458 
459  //================================================ Get Fourier transforms of the maps
460  forwardFourierObj1 = fftw_plan_dft_3d ( static_cast< int > ( xD ), static_cast< int > ( yD ), static_cast< int > ( zD ), tmpIn1, tmpOut1, FFTW_FORWARD , FFTW_ESTIMATE );
461  forwardFourierObj2 = fftw_plan_dft_3d ( static_cast< int > ( xD ), static_cast< int > ( yD ), static_cast< int > ( zD ), tmpIn2, tmpOut2, FFTW_FORWARD , FFTW_ESTIMATE );
462  inverseFourierCombo = fftw_plan_dft_3d ( static_cast< int > ( xD ), static_cast< int > ( yD ), static_cast< int > ( zD ), resOut, resIn , FFTW_BACKWARD, FFTW_ESTIMATE );
463 
464  //================================================ Done
465  return ;
466 
467 }

◆ 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 506 of file ProSHADE_overlay.cpp.

507 {
508  //================================================ Initialise local variables
509  double normFactor = static_cast<double> ( xD * yD * zD );
510  proshade_signed arrPos;
511 
512  //================================================ Allocate local memory
513  proshade_single *mins = new proshade_single[3];
514  proshade_single *maxs = new proshade_single[3];
515 
516  //================================================ Check local memory
517  ProSHADE_internal_misc::checkMemoryAllocation ( mins, __FILE__, __LINE__, __func__ );
518  ProSHADE_internal_misc::checkMemoryAllocation ( maxs, __FILE__, __LINE__, __func__ );
519 
520  //================================================ Determine reciprocal space indexing
521  mins[0] = std::floor ( static_cast< proshade_single > ( xD ) / -2.0f );
522  mins[1] = std::floor ( static_cast< proshade_single > ( yD ) / -2.0f );
523  mins[2] = std::floor ( static_cast< proshade_single > ( zD ) / -2.0f );
524 
525  maxs[0] = -mins[0];
526  maxs[1] = -mins[1];
527  maxs[2] = -mins[2];
528 
529  if ( xD % 2 == 0 ) { mins[0] += 1.0f; }
530  if ( yD % 2 == 0 ) { mins[1] += 1.0f; }
531  if ( zD % 2 == 0 ) { mins[2] += 1.0f; }
532 
533  //================================================ Combine the coefficients
534  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xD ); xIt++ )
535  {
536  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yD ); yIt++ )
537  {
538  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zD / 2 ) + 1; zIt++ )
539  {
540  //==================================== Find indices
541  arrPos = zIt + static_cast< proshade_signed > ( zD ) * ( yIt + static_cast< proshade_signed > ( yD ) * xIt );
542 
543  //==================================== Combine
545  &tmpOut1[arrPos][1],
546  &tmpOut2[arrPos][0],
547  &tmpOut2[arrPos][1],
548  &resOut[arrPos][0],
549  &resOut[arrPos][1] );
550 
551  //==================================== Save
552  resOut[arrPos][0] /= normFactor;
553  resOut[arrPos][1] /= normFactor;
554  }
555  }
556  }
557 
558  //================================================ Delete local memory
559  delete[] mins;
560  delete[] maxs;
561 
562  //================================================ Done
563  return ;
564 
565 }

◆ 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 1263 of file ProSHADE_overlay.cpp.

1264 {
1265  //================================================ Longitude angular thresholds
1266  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1267  {
1268  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 ) ) );
1269  }
1270 
1271  //================================================ Lattitude angular thresholds
1272  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1273  {
1274  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 ) );
1275  }
1276 
1277  //================================================ Done
1278  return ;
1279 
1280 }

◆ 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 677 of file ProSHADE_overlay.cpp.

678 {
679  //================================================ Compute
680  *addXPre = ( xDim - xDimIndices ) / 2;
681  *addYPre = ( yDim - yDimIndices ) / 2;
682  *addZPre = ( zDim - zDimIndices ) / 2;
683  *addXPost = *addXPre; if ( ( xDim - xDimIndices ) % 2 == 1 ) { *addXPost += 1; }
684  *addYPost = *addYPre; if ( ( yDim - yDimIndices ) % 2 == 1 ) { *addYPost += 1; }
685  *addZPost = *addZPre; if ( ( zDim - zDimIndices ) % 2 == 1 ) { *addZPost += 1; }
686 
687  //================================================ Done
688  return ;
689 
690 }

◆ 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 578 of file ProSHADE_overlay.cpp.

579 {
580  //================================================ Initialise variables
581  proshade_signed arrPos;
582  *mapPeak = 0.0;
583 
584  //================================================ Search the map
585  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> ( xD ); uIt++ )
586  {
587  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> ( yD ); vIt++ )
588  {
589  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> ( zD ); wIt++ )
590  {
591  arrPos = wIt + static_cast< proshade_signed > ( zD ) * ( vIt + static_cast< proshade_signed > ( yD ) * uIt );
592  if ( resIn[arrPos][0] > *mapPeak )
593  {
594  *mapPeak = resIn[arrPos][0];
595  *trsX = static_cast< proshade_double > ( uIt );
596  *trsY = static_cast< proshade_double > ( vIt );
597  *trsZ = static_cast< proshade_double > ( wIt );
598  }
599  }
600  }
601  }
602 
603  //================================================ Done
604  return ;
605 
606 }

◆ 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 480 of file ProSHADE_overlay.cpp.

481 {
482  //================================================ Release memory
483  fftw_destroy_plan ( forwardFourierObj1 );
484  fftw_destroy_plan ( forwardFourierObj2 );
485  fftw_destroy_plan ( inverseFourierCombo );
486  delete[] tmpIn1;
487  delete[] tmpIn2;
488  delete[] tmpOut1;
489  delete[] tmpOut2;
490  delete[] resOut;
491 
492  //======================================== Done
493  return ;
494 
495 }

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

◆ 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 1147 of file ProSHADE_overlay.cpp.

1148 {
1149  //================================================ Initialise internal variables
1150  proshade_unsign oneDim = shBand * 2;
1151 
1152  //================================================ Allocate memory
1153  sigR = new double [(oneDim * oneDim * 4)];
1154  sigI = sigR + (oneDim * oneDim * 2);
1155  rcoeffs = new double [(oneDim * oneDim * 2)];
1156  icoeffs = rcoeffs + (oneDim * oneDim);
1157  weights = new double [4 * shBand];
1158  workspace = new double [( 20 * shBand * shBand ) + ( 24 * shBand )];
1159 
1160  //================================================ Check memory allocation
1161  ProSHADE_internal_misc::checkMemoryAllocation ( sigR, __FILE__, __LINE__, __func__ );
1162  ProSHADE_internal_misc::checkMemoryAllocation ( rcoeffs, __FILE__, __LINE__, __func__ );
1163  ProSHADE_internal_misc::checkMemoryAllocation ( weights, __FILE__, __LINE__, __func__ );
1164  ProSHADE_internal_misc::checkMemoryAllocation ( workspace, __FILE__, __LINE__, __func__ );
1165 
1166  //================================================ Create the cosine/sine transform plan
1167  idctPlan = fftw_plan_r2r_1d ( static_cast< int > ( oneDim ), weights, workspace, FFTW_REDFT01, FFTW_ESTIMATE );
1168 
1169  //================================================ Set up the discrete Fourier transform
1170  int rank, howmany_rank;
1171  fftw_iodim dims[1], howmany_dims[1];
1172 
1173  rank = 1;
1174  dims[0].n = 2 * static_cast< int > ( shBand );
1175  dims[0].is = 2 * static_cast< int > ( shBand );
1176  dims[0].os = 1;
1177  howmany_rank = 1;
1178  howmany_dims[0].n = 2 * static_cast< int > ( shBand );
1179  howmany_dims[0].is = 1;
1180  howmany_dims[0].os = 2 * static_cast< int > ( shBand );
1181 
1182  //================================================ Create the discrete Fourier transform
1183  ifftPlan = fftw_plan_guru_split_dft ( rank, dims, howmany_rank, howmany_dims, sigR, sigI, rcoeffs, icoeffs, FFTW_ESTIMATE );
1184 
1185  //================================================ Done
1186  return ;
1187 
1188 }

◆ 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 706 of file ProSHADE_overlay.cpp.

707 {
708  //================================================ Initialise local variables
709  proshade_unsign newMapIndex = 0;
710  proshade_unsign oldMapIndex = 0;
711 
712  //================================================ Copy and padd as appropriate
713  for ( proshade_unsign xIt = 0; xIt < xDim; xIt++ )
714  {
715  for ( proshade_unsign yIt = 0; yIt < yDim; yIt++ )
716  {
717  for ( proshade_unsign zIt = 0; zIt < zDim; zIt++ )
718  {
719  //==================================== Find map location
720  newMapIndex = zIt + zDim * ( yIt + yDim * xIt );
721 
722  //==================================== If less than needed, add zero
723  if ( xIt < addXPre ) { newMap[newMapIndex] = 0.0; continue; }
724  if ( yIt < addYPre ) { newMap[newMapIndex] = 0.0; continue; }
725  if ( zIt < addZPre ) { newMap[newMapIndex] = 0.0; continue; }
726 
727  //==================================== If more than needed, add zero
728  if ( xIt >= ( addXPre + xDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
729  if ( yIt >= ( addYPre + yDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
730  if ( zIt >= ( addZPre + zDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
731 
732  //==================================== If not padding, copy original values
733  oldMapIndex = (zIt-addZPre) + zDimIndices * ( (yIt-addYPre) + yDimIndices * (xIt-addXPre) );
734  newMap[newMapIndex] = oldMap[oldMapIndex];
735  }
736  }
737  }
738 
739  //======================================== Done
740  return ;
741 
742 }
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:1756
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:3929
ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpaceInPlace
void rotateMapRealSpaceInPlace(proshade_double eulA, proshade_double eulB, proshade_double eulG)
This function rotates a map based on the given Euler angles in place.
Definition: ProSHADE_overlay.cpp:982
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:4029
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:3909
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:4059
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:3919
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:811
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:3969
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:140
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:3939
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:1712
ProSHADE_internal_data::ProSHADE_data::computeOptimalTranslation
void computeOptimalTranslation(proshade_double eulA, proshade_double eulB, proshade_double eulG, proshade_single trsX, proshade_single trsY, proshade_single 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::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:352
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:398
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:3630
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:3949
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:1612
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:4069
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:3959
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:578
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:4009
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:4049
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:67
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:3999
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:3979
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:3857
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:43
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:4039
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:763
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:618
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:501
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:3989
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:4019