ProSHADE  0.7.6.1 (AUG 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 computeTranslationsFromPeak (ProSHADE_internal_data::ProSHADE_data *staticStructure, ProSHADE_internal_data::ProSHADE_data *movingStructure, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ)
 This function computes the translation in Angstroms that corresponds to the translation function peak. 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 367 of file ProSHADE_overlay.cpp.

368 {
369  //================================================ Allocate memory
370  tmpIn1 = new fftw_complex[xD * yD * zD];
371  tmpOut1 = new fftw_complex[xD * yD * zD];
372  tmpIn2 = new fftw_complex[xD * yD * zD];
373  tmpOut2 = new fftw_complex[xD * yD * zD];
374  resIn = new fftw_complex[xD * yD * zD];
375  resOut = new fftw_complex[xD * yD * zD];
376 
377  //================================================ Check memory allocation
378  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn1 , __FILE__, __LINE__, __func__ );
379  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut1, __FILE__, __LINE__, __func__ );
380  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn2, __FILE__, __LINE__, __func__ );
381  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut2, __FILE__, __LINE__, __func__ );
382  ProSHADE_internal_misc::checkMemoryAllocation ( resIn, __FILE__, __LINE__, __func__ );
383  ProSHADE_internal_misc::checkMemoryAllocation ( resOut, __FILE__, __LINE__, __func__ );
384 
385  //================================================ Get Fourier transforms of the maps
386  forwardFourierObj1 = fftw_plan_dft_3d ( static_cast< int > ( xD ), static_cast< int > ( yD ), static_cast< int > ( zD ), tmpIn1, tmpOut1, FFTW_FORWARD , FFTW_ESTIMATE );
387  forwardFourierObj2 = fftw_plan_dft_3d ( static_cast< int > ( xD ), static_cast< int > ( yD ), static_cast< int > ( zD ), tmpIn2, tmpOut2, FFTW_FORWARD , FFTW_ESTIMATE );
388  inverseFourierCombo = fftw_plan_dft_3d ( static_cast< int > ( xD ), static_cast< int > ( yD ), static_cast< int > ( zD ), resOut, resIn , FFTW_BACKWARD, FFTW_ESTIMATE );
389 
390  //================================================ Done
391  return ;
392 
393 }

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

433 {
434  //================================================ Initialise local variables
435  double normFactor = static_cast<double> ( xD * yD * zD );
436  proshade_signed arrPos;
437 
438  //================================================ Combine the coefficients
439  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xD ); xIt++ )
440  {
441  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yD ); yIt++ )
442  {
443  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zD ); zIt++ )
444  {
445  //==================================== Find indices
446  arrPos = zIt + static_cast< proshade_signed > ( zD ) * ( yIt + static_cast< proshade_signed > ( yD ) * xIt );
447 
448  //==================================== Combine
450  &tmpOut1[arrPos][1],
451  &tmpOut2[arrPos][0],
452  &tmpOut2[arrPos][1],
453  &resOut[arrPos][0],
454  &resOut[arrPos][1] );
455 
456  //==================================== Save
457  resOut[arrPos][0] /= normFactor;
458  resOut[arrPos][1] /= normFactor;
459  }
460  }
461  }
462 
463  //================================================ Done
464  return ;
465 
466 }

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

1178 {
1179  //================================================ Longitude angular thresholds
1180  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1181  {
1182  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 ) ) );
1183  }
1184 
1185  //================================================ Lattitude angular thresholds
1186  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1187  {
1188  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 ) );
1189  }
1190 
1191  //================================================ Done
1192  return ;
1193 
1194 }

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

579 {
580  //================================================ Compute
581  *addXPre = ( xDim - xDimIndices ) / 2;
582  *addYPre = ( yDim - yDimIndices ) / 2;
583  *addZPre = ( zDim - zDimIndices ) / 2;
584  *addXPost = *addXPre; if ( ( xDim - xDimIndices ) % 2 == 1 ) { *addXPost += 1; }
585  *addYPost = *addYPre; if ( ( yDim - yDimIndices ) % 2 == 1 ) { *addYPost += 1; }
586  *addZPost = *addZPre; if ( ( zDim - zDimIndices ) % 2 == 1 ) { *addZPost += 1; }
587 
588  //================================================ Done
589  return ;
590 
591 }

◆ computeTranslationsFromPeak()

void ProSHADE_internal_overlay::computeTranslationsFromPeak ( ProSHADE_internal_data::ProSHADE_data staticStructure,
ProSHADE_internal_data::ProSHADE_data movingStructure,
proshade_double *  trsX,
proshade_double *  trsY,
proshade_double *  trsZ 
)

This function computes the translation in Angstroms that corresponds to the translation function peak.

Parameters
[in]staticStructureA pointer to the data class object of the static structure.
[in]movingStructureA pointer to the data class object of the 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.

Definition at line 220 of file ProSHADE_overlay.cpp.

221 {
222  //================================================ Compute map movement
223  proshade_single xSamplRate = ( static_cast< proshade_single > ( staticStructure->getXDimSize() ) / static_cast< proshade_single > ( staticStructure->getXDim() ) );
224  proshade_single ySamplRate = ( static_cast< proshade_single > ( staticStructure->getYDimSize() ) / static_cast< proshade_single > ( staticStructure->getYDim() ) );
225  proshade_single zSamplRate = ( static_cast< proshade_single > ( staticStructure->getZDimSize() ) / static_cast< proshade_single > ( staticStructure->getZDim() ) );
226 
227  proshade_single xCor = static_cast< proshade_single > ( static_cast< proshade_single > ( staticStructure->xFrom - movingStructure->xFrom ) * xSamplRate );
228  proshade_single yCor = static_cast< proshade_single > ( static_cast< proshade_single > ( staticStructure->yFrom - movingStructure->yFrom ) * ySamplRate );
229  proshade_single zCor = static_cast< proshade_single > ( static_cast< proshade_single > ( staticStructure->zFrom - movingStructure->zFrom ) * zSamplRate );
230 
231  proshade_single xMov = static_cast< proshade_single > ( 0.0 - static_cast< proshade_double > ( xCor ) - ( *trsX * static_cast< proshade_double > ( xSamplRate ) ) );
232  proshade_single yMov = static_cast< proshade_single > ( 0.0 - static_cast< proshade_double > ( yCor ) - ( *trsY * static_cast< proshade_double > ( ySamplRate ) ) );
233  proshade_single zMov = static_cast< proshade_single > ( 0.0 - static_cast< proshade_double > ( zCor ) - ( *trsZ * static_cast< proshade_double > ( zSamplRate ) ) );
234 
235  //================================================ Save translation vector back
236  *trsX = static_cast< proshade_double > ( -xMov );
237  *trsY = static_cast< proshade_double > ( -yMov );
238  *trsZ = static_cast< proshade_double > ( -zMov );
239 
240  //================================================ Done
241  return ;
242 
243 }

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

480 {
481  //================================================ Initialise variables
482  proshade_signed arrPos;
483  *mapPeak = 0.0;
484 
485  //================================================ Search the map
486  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> ( xD ); uIt++ )
487  {
488  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> ( yD ); vIt++ )
489  {
490  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> ( zD ); wIt++ )
491  {
492  arrPos = wIt + static_cast< proshade_signed > ( zD ) * ( vIt + static_cast< proshade_signed > ( yD ) * uIt );
493  if ( resIn[arrPos][0] > *mapPeak )
494  {
495  *mapPeak = resIn[arrPos][0];
496  *trsX = static_cast< proshade_double > ( uIt );
497  *trsY = static_cast< proshade_double > ( vIt );
498  *trsZ = static_cast< proshade_double > ( wIt );
499  }
500  }
501  }
502  }
503 
504  //================================================ Done
505  return ;
506 
507 }

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

407 {
408  //================================================ Release memory
409  fftw_destroy_plan ( forwardFourierObj1 );
410  fftw_destroy_plan ( forwardFourierObj2 );
411  fftw_destroy_plan ( inverseFourierCombo );
412  delete[] tmpIn1;
413  delete[] tmpIn2;
414  delete[] tmpOut1;
415  delete[] tmpOut2;
416  delete[] resOut;
417 
418  //======================================== Done
419  return ;
420 
421 }

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

75 {
76  //================================================ Read in the structures
77  staticStructure->readInStructure ( settings->inputFiles.at(0), 0, settings );
78  movingStructure->readInStructure ( settings->inputFiles.at(1), 1, settings );
79 
80  //================================================ Internal data processing (COM, norm, mask, extra space)
81  staticStructure->processInternalMap ( settings );
82  movingStructure->processInternalMap ( settings );
83 
84  //================================================ Map to sphere
85  staticStructure->mapToSpheres ( settings );
86  movingStructure->mapToSpheres ( settings );
87 
88  //================================================ Get spherical harmonics
89  staticStructure->computeSphericalHarmonics ( settings );
90  movingStructure->computeSphericalHarmonics ( settings );
91 
92  //================================================ Get the rotation function of the pair
93  movingStructure->getOverlayRotationFunction ( settings, staticStructure );
94 
95  //================================================ Get inverse SO(3) map top peak Euler angle values
97  std::min ( staticStructure->getMaxBand(), movingStructure->getMaxBand() ) * 2,
98  eulA, eulB, eulG, settings );
99 
100  //================================================ Done
101  return ;
102 
103 }

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

124 {
125  //================================================ Read in the structures
126  staticStructure->readInStructure ( settings->inputFiles.at(0), 0, settings );
127  movingStructure->readInStructure ( settings->inputFiles.at(1), 1, settings );
128 
129  //================================================ Determine spherical harmonics variables from the first structure (otherwise, they would be determined from the second structure)
130  settings->determineAllSHValues ( staticStructure->xDimIndices, staticStructure->yDimIndices,
131  staticStructure->xDimSize, staticStructure->yDimSize, staticStructure->zDimSize );
132 
133  //================================================ Internal data processing (COM, norm, mask, extra space)
134  staticStructure->processInternalMap ( settings );
135  movingStructure->processInternalMap ( settings );
136 
137  //================================================ Rotate map
138  std::vector< proshade_double > rotCen = movingStructure->rotateMapRealSpaceInPlace ( eulA, eulB, eulG );
139 
140  //================================================ Zero padding for smaller structure
141  staticStructure->zeroPaddToDims ( std::max ( staticStructure->getXDim(), movingStructure->getXDim() ),
142  std::max ( staticStructure->getYDim(), movingStructure->getYDim() ),
143  std::max ( staticStructure->getZDim(), movingStructure->getZDim() ) );
144  movingStructure->zeroPaddToDims ( std::max ( staticStructure->getXDim(), movingStructure->getXDim() ),
145  std::max ( staticStructure->getYDim(), movingStructure->getYDim() ),
146  std::max ( staticStructure->getZDim(), movingStructure->getZDim() ) );
147 
148  //================================================ Report progress
149  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting translation function computation." );
150 
151  //================================================ Compute the translation function map
152  movingStructure->computeTranslationMap ( staticStructure );
153 
154  //================================================ Find highest peak in translation function
155  proshade_double mapPeak = 0.0;
156  proshade_unsign xDimS = staticStructure->getXDim();
157  proshade_unsign yDimS = staticStructure->getYDim();
158  proshade_unsign zDimS = staticStructure->getZDim();
159  ProSHADE_internal_overlay::findHighestValueInMap ( movingStructure->getTranslationFnPointer(), xDimS, yDimS, zDimS, trsX, trsY, trsZ, &mapPeak );
160 
161  //================================================ Dont translate over half
162  if ( *trsX > ( static_cast< proshade_double > ( xDimS ) / 2.0 ) ) { *trsX = *trsX - static_cast< proshade_double > ( xDimS ); }
163  if ( *trsY > ( static_cast< proshade_double > ( yDimS ) / 2.0 ) ) { *trsY = *trsY - static_cast< proshade_double > ( yDimS ); }
164  if ( *trsZ > ( static_cast< proshade_double > ( zDimS ) / 2.0 ) ) { *trsZ = *trsZ - static_cast< proshade_double > ( zDimS ); }
165 
166  //================================================ Convert the map positions onto translation in Angstroms and save the results
167  computeTranslationsFromPeak ( staticStructure, movingStructure, trsX, trsY, trsZ );
168 
169  movingStructure->originalPdbTransX = *trsX;
170  movingStructure->originalPdbTransY = *trsY;
171  movingStructure->originalPdbTransZ = *trsZ;
172 
173  proshade_single xMov = static_cast< proshade_single > ( -(*trsX) );
174  proshade_single yMov = static_cast< proshade_single > ( -(*trsY) );
175  proshade_single zMov = static_cast< proshade_single > ( -(*trsZ) );
176 
177  //================================================ Report progress
178  std::stringstream hlpSS;
179  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 ) );
180 
181  //================================================ Save original from variables for PDB output
182  movingStructure->mapMovFromsChangeX = static_cast< proshade_double > ( movingStructure->xFrom );
183  movingStructure->mapMovFromsChangeY = static_cast< proshade_double > ( movingStructure->yFrom );
184  movingStructure->mapMovFromsChangeZ = static_cast< proshade_double > ( movingStructure->zFrom );
185 
186  //================================================ Move the map
187  ProSHADE_internal_mapManip::moveMapByIndices ( &xMov, &yMov, &zMov, movingStructure->getXDimSize(), movingStructure->getYDimSize(), movingStructure->getZDimSize(),
188  movingStructure->getXFromPtr(), movingStructure->getXToPtr(),
189  movingStructure->getYFromPtr(), movingStructure->getYToPtr(),
190  movingStructure->getZFromPtr(), movingStructure->getZToPtr(),
191  movingStructure->getXAxisOrigin(), movingStructure->getYAxisOrigin(), movingStructure->getZAxisOrigin() );
192 
193  ProSHADE_internal_mapManip::moveMapByFourier ( movingStructure->getInternalMap(), xMov, yMov, zMov,
194  movingStructure->getXDimSize(), movingStructure->getYDimSize(), movingStructure->getZDimSize(),
195  static_cast< proshade_signed > ( movingStructure->getXDim() ), static_cast< proshade_signed > ( movingStructure->getYDim() ),
196  static_cast< proshade_signed > ( movingStructure->getZDim() ) );
197 
198  //================================================ Report progress
199  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, hlpSS.str() );
200  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Translation function computation complete." );
201 
202  //================================================ Keep only the change in from and to variables
203  movingStructure->mapMovFromsChangeX = static_cast< proshade_double > ( movingStructure->xFrom ) - movingStructure->mapMovFromsChangeX;
204  movingStructure->mapMovFromsChangeY = static_cast< proshade_double > ( movingStructure->yFrom ) - movingStructure->mapMovFromsChangeY;
205  movingStructure->mapMovFromsChangeZ = static_cast< proshade_double > ( movingStructure->zFrom ) - movingStructure->mapMovFromsChangeZ;
206 
207  //================================================ Done
208  return ;
209 
210 }

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

1062 {
1063  //================================================ Initialise internal variables
1064  proshade_unsign oneDim = shBand * 2;
1065 
1066  //================================================ Allocate memory
1067  sigR = new double [(oneDim * oneDim * 4)];
1068  sigI = sigR + (oneDim * oneDim * 2);
1069  rcoeffs = new double [(oneDim * oneDim * 2)];
1070  icoeffs = rcoeffs + (oneDim * oneDim);
1071  weights = new double [4 * shBand];
1072  workspace = new double [( 20 * shBand * shBand ) + ( 24 * shBand )];
1073 
1074  //================================================ Check memory allocation
1075  ProSHADE_internal_misc::checkMemoryAllocation ( sigR, __FILE__, __LINE__, __func__ );
1076  ProSHADE_internal_misc::checkMemoryAllocation ( rcoeffs, __FILE__, __LINE__, __func__ );
1077  ProSHADE_internal_misc::checkMemoryAllocation ( weights, __FILE__, __LINE__, __func__ );
1078  ProSHADE_internal_misc::checkMemoryAllocation ( workspace, __FILE__, __LINE__, __func__ );
1079 
1080  //================================================ Create the cosine/sine transform plan
1081  idctPlan = fftw_plan_r2r_1d ( static_cast< int > ( oneDim ), weights, workspace, FFTW_REDFT01, FFTW_ESTIMATE );
1082 
1083  //================================================ Set up the discrete Fourier transform
1084  int rank, howmany_rank;
1085  fftw_iodim dims[1], howmany_dims[1];
1086 
1087  rank = 1;
1088  dims[0].n = 2 * static_cast< int > ( shBand );
1089  dims[0].is = 2 * static_cast< int > ( shBand );
1090  dims[0].os = 1;
1091  howmany_rank = 1;
1092  howmany_dims[0].n = 2 * static_cast< int > ( shBand );
1093  howmany_dims[0].is = 1;
1094  howmany_dims[0].os = 2 * static_cast< int > ( shBand );
1095 
1096  //================================================ Create the discrete Fourier transform
1097  ifftPlan = fftw_plan_guru_split_dft ( rank, dims, howmany_rank, howmany_dims, sigR, sigI, rcoeffs, icoeffs, FFTW_ESTIMATE );
1098 
1099  //================================================ Done
1100  return ;
1101 
1102 }

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

608 {
609  //================================================ Initialise local variables
610  proshade_unsign newMapIndex = 0;
611  proshade_unsign oldMapIndex = 0;
612 
613  //================================================ Copy and padd as appropriate
614  for ( proshade_unsign xIt = 0; xIt < xDim; xIt++ )
615  {
616  for ( proshade_unsign yIt = 0; yIt < yDim; yIt++ )
617  {
618  for ( proshade_unsign zIt = 0; zIt < zDim; zIt++ )
619  {
620  //==================================== Find map location
621  newMapIndex = zIt + zDim * ( yIt + yDim * xIt );
622 
623  //==================================== If less than needed, add zero
624  if ( xIt < addXPre ) { newMap[newMapIndex] = 0.0; continue; }
625  if ( yIt < addYPre ) { newMap[newMapIndex] = 0.0; continue; }
626  if ( zIt < addZPre ) { newMap[newMapIndex] = 0.0; continue; }
627 
628  //==================================== If more than needed, add zero
629  if ( xIt >= ( addXPre + xDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
630  if ( yIt >= ( addYPre + yDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
631  if ( zIt >= ( addZPre + zDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
632 
633  //==================================== If not padding, copy original values
634  oldMapIndex = (zIt-addZPre) + zDimIndices * ( (yIt-addYPre) + yDimIndices * (xIt-addXPre) );
635  newMap[newMapIndex] = oldMap[oldMapIndex];
636  }
637  }
638  }
639 
640  //======================================== Done
641  return ;
642 
643 }
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:1843
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::xDimIndices
proshade_unsign xDimIndices
This is the size of the map cell x dimension in indices.
Definition: ProSHADE_data.hpp:65
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:3898
ProSHADE_settings::determineAllSHValues
void determineAllSHValues(proshade_unsign xDim, proshade_unsign yDim, proshade_single xDimAngs, proshade_single yDimAngs, proshade_single zDimAngs)
This function determines all the required values for spherical harmonics computation.
Definition: ProSHADE.cpp:1615
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:3998
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:3878
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:4028
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:3888
ProSHADE_internal_data::ProSHADE_data::zDimSize
proshade_single zDimSize
This is the size of the map cell z dimension in Angstroms.
Definition: ProSHADE_data.hpp:61
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::getXFromPtr
proshade_signed * getXFromPtr(void)
This function allows access to the map start along the X axis.
Definition: ProSHADE_data.cpp:3938
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:142
ProSHADE_internal_data::ProSHADE_data::yDimIndices
proshade_unsign yDimIndices
This is the size of the map cell y dimension in indices.
Definition: ProSHADE_data.hpp:66
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:3908
ProSHADE_internal_data::ProSHADE_data::originalPdbTransX
proshade_double originalPdbTransX
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
Definition: ProSHADE_data.hpp:105
ProSHADE_internal_data::ProSHADE_data::readInStructure
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings, proshade_double *maskArr=nullptr, proshade_unsign maskXDim=0, proshade_unsign maskYDim=0, proshade_unsign maskZDim=0, proshade_double *weightsArr=nullptr, proshade_unsign weigXDim=0, proshade_unsign weigYDim=0, proshade_unsign weigZDim=0)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
Definition: ProSHADE_data.cpp:509
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:1799
ProSHADE_internal_data::ProSHADE_data::originalPdbTransY
proshade_double originalPdbTransY
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
Definition: ProSHADE_data.hpp:106
ProSHADE_internal_data::ProSHADE_data::xDimSize
proshade_single xDimSize
This is the size of the map cell x dimension in Angstroms.
Definition: ProSHADE_data.hpp:59
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:324
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:3599
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:3918
ProSHADE_internal_data::ProSHADE_data::originalPdbTransZ
proshade_double originalPdbTransZ
The optimal translation vector as it relates to the original PDB positions (and not the ProSHADE inte...
Definition: ProSHADE_data.hpp:107
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:1699
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:4038
ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpaceInPlace
std::vector< proshade_double > 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:891
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:3928
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:479
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:3978
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:4018
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:3968
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:3948
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:3826
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:4008
ProSHADE_internal_overlay::computeTranslationsFromPeak
void computeTranslationsFromPeak(ProSHADE_internal_data::ProSHADE_data *staticStructure, ProSHADE_internal_data::ProSHADE_data *movingStructure, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ)
This function computes the translation in Angstroms that corresponds to the translation function peak...
Definition: ProSHADE_overlay.cpp:220
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_data::ProSHADE_data::yDimSize
proshade_single yDimSize
This is the size of the map cell y dimension in Angstroms.
Definition: ProSHADE_data.hpp:60
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:519
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:3958
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:3988