ProSHADE  0.7.6.0 (JUL 2021)
Protein Shape Detection
ProSHADE_overlay.cpp
Go to the documentation of this file.
1 
22 //==================================================== ProSHADE
23 #include "ProSHADE_overlay.hpp"
24 
36 {
37  //================================================ Report progress
38  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting rotation function computation." );
39 
40  //================================================ Compute un-weighted E matrices and their weights
41  ProSHADE_internal_distances::computeEMatrices ( obj2, this, settings );
42 
43  //================================================ Normalise E matrices by the magnitudes
44  ProSHADE_internal_distances::normaliseEMatrices ( obj2, this, settings );
45 
46  //================================================ Generate SO(3) coefficients
48 
49  //================================================ Compute the inverse SO(3) Fourier Transform (SOFT) on the newly computed coefficients
51 
52  //================================================ Report completion
53  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Rotation function obtained." );
54 
55  //================================================ Done
56  return ;
57 
58 }
59 
69 {
70  //================================================ First. determine the sampling rates (value to multiply indices with to get Angstroms)
71  proshade_double xSamplRate = static_cast<proshade_double> ( this->xDimSizeOriginal / static_cast< proshade_single > ( this->xDimIndicesOriginal ) );
72  proshade_double ySamplRate = static_cast<proshade_double> ( this->yDimSizeOriginal / static_cast< proshade_single > ( this->yDimIndicesOriginal ) );
73  proshade_double zSamplRate = static_cast<proshade_double> ( this->zDimSizeOriginal / static_cast< proshade_single > ( this->zDimIndicesOriginal ) );
74 
75  //================================================ Compute the rotation centre for the co-ordinates
76  proshade_double xRotPos = ( static_cast<proshade_double> ( this->xFrom ) - this->mapMovFromsChangeX ) * xSamplRate + // Corner X position in Angstroms
77  ( ( ( static_cast<proshade_double> ( this->xTo ) - static_cast<proshade_double> ( this->xFrom ) ) / 2.0 ) * xSamplRate ); // Half of box X size
78 
79  proshade_double yRotPos = ( static_cast<proshade_double> ( this->yFrom ) - this->mapMovFromsChangeY ) * ySamplRate + // Corner Y position in Angstroms
80  ( ( ( static_cast<proshade_double> ( this->yTo ) - static_cast<proshade_double> ( this->yFrom ) ) / 2.0 ) * ySamplRate ); // Half of box Y size
81 
82  proshade_double zRotPos = ( static_cast<proshade_double> ( this->zFrom ) - this->mapMovFromsChangeZ ) * zSamplRate + // Corner Z position in Angstroms
83  ( ( ( static_cast<proshade_double> ( this->zTo ) - static_cast<proshade_double> ( this->zFrom ) ) / 2.0 ) * zSamplRate ); // Half of box Z size
84 
85  //============================================ Modify by change during ProSHADE map processing
86  this->originalPdbRotCenX = xRotPos - ( this->mapCOMProcessChangeX / 2.0 );
87  this->originalPdbRotCenY = yRotPos - ( this->mapCOMProcessChangeY / 2.0 );
88  this->originalPdbRotCenZ = zRotPos - ( this->mapCOMProcessChangeZ / 2.0 );
89 
90  //================================================ Done
91  return ;
92 
93 }
94 
109 void ProSHADE_internal_data::ProSHADE_data::computeOptimalTranslation ( proshade_double eulA, proshade_double eulB, proshade_double eulG, proshade_single trsX, proshade_single trsY, proshade_single trsZ )
110 {
111  //================================================ Reset class variables
112  this->originalPdbTransX = 0.0;
113  this->originalPdbTransY = 0.0;
114  this->originalPdbTransZ = 0.0;
115 
116  //================================================ Correctly apply any map modifications that ProSHADE may have done to the map to make sure map matches co-ordinates.
117  if ( ( eulA != 0.0 ) || ( eulB != 0.0 ) || ( eulG != 0.0 ) )
118  {
119  //============================================ If rotation is to be done, then ProSHADE processing map changes are already dealt with
120  ;
121  }
122  else
123  {
124  //============================================ In not, then they need to be added
125  this->originalPdbTransX = this->mapCOMProcessChangeX;
126  this->originalPdbTransY = this->mapCOMProcessChangeY;
127  this->originalPdbTransZ = this->mapCOMProcessChangeZ;
128  }
129 
130  //================================================ Save the values
131  this->originalPdbTransX += static_cast< proshade_double > ( trsX );
132  this->originalPdbTransY += static_cast< proshade_double > ( trsY );
133  this->originalPdbTransZ += static_cast< proshade_double > ( trsZ );
134 
135  //================================================ Done
136  return ;
137 
138 }
139 
154 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 )
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 }
184 
203 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 )
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 }
305 
311 std::vector< proshade_double > ProSHADE_internal_data::ProSHADE_data::getBestTranslationMapPeaksAngstrom ( ProSHADE_internal_data::ProSHADE_data* staticStructure, proshade_double eulA, proshade_double eulB, proshade_double eulG )
312 {
313  //================================================ Initialise local variables
314  std::vector< proshade_double > ret;
315  proshade_double mapPeak = 0.0;
316  proshade_double trsX = 0.0, trsY = 0.0, trsZ = 0.0;
317  ProSHADE_internal_overlay::findHighestValueInMap ( this->getTranslationFnPointer(),
318  staticStructure->getXDim(),
319  staticStructure->getYDim(),
320  staticStructure->getZDim(),
321  &trsX,
322  &trsY,
323  &trsZ,
324  &mapPeak );
325 
326  //================================================ Dont translate over half
327  if ( trsX > ( static_cast< proshade_double > ( staticStructure->getXDim() ) / 2.0 ) ) { trsX = trsX - static_cast< proshade_double > ( this->getXDim() ); }
328  if ( trsY > ( static_cast< proshade_double > ( staticStructure->getYDim() ) / 2.0 ) ) { trsY = trsY - static_cast< proshade_double > ( this->getYDim() ); }
329  if ( trsZ > ( static_cast< proshade_double > ( staticStructure->getZDim() ) / 2.0 ) ) { trsZ = trsZ - static_cast< proshade_double > ( this->getZDim() ); }
330 
331  //================================================ Move map
332  proshade_single xCor = static_cast< proshade_single > ( static_cast< proshade_double > ( staticStructure->xFrom - this->xFrom ) *
333  ( static_cast< proshade_double > ( staticStructure->getXDimSize() ) /
334  static_cast< proshade_double > ( staticStructure->getXDim() ) ) );
335  proshade_single yCor = static_cast< proshade_single > ( static_cast< proshade_double > ( staticStructure->yFrom - this->yFrom ) *
336  ( static_cast< proshade_double > ( staticStructure->getYDimSize() ) /
337  static_cast< proshade_double > ( staticStructure->getYDim() ) ) );
338  proshade_single zCor = static_cast< proshade_single > ( static_cast< proshade_double > ( staticStructure->zFrom - this->zFrom ) *
339  ( static_cast< proshade_double > ( staticStructure->getZDimSize() ) /
340  static_cast< proshade_double > ( staticStructure->getZDim() ) ) );
341  proshade_single xMov = static_cast< proshade_single > ( staticStructure->mapCOMProcessChangeX - static_cast< proshade_double > ( xCor ) -
342  ( trsX * static_cast<proshade_double> ( staticStructure->getXDimSize() ) /
343  static_cast<proshade_double> ( staticStructure->getXDim() ) ) );
344  proshade_single yMov = static_cast< proshade_single > ( staticStructure->mapCOMProcessChangeY - static_cast< proshade_double > ( yCor ) -
345  ( trsY * static_cast<proshade_double> ( staticStructure->getYDimSize() ) /
346  static_cast<proshade_double> ( staticStructure->getYDim() ) ) );
347  proshade_single zMov = static_cast< proshade_single > ( staticStructure->mapCOMProcessChangeZ - static_cast< proshade_double > ( zCor ) -
348  ( trsZ * static_cast<proshade_double> ( staticStructure->getZDimSize() ) /
349  static_cast<proshade_double> ( staticStructure->getZDim() ) ) );
350  proshade_single modXMov = xMov;
351  proshade_single modYMov = yMov;
352  proshade_single modZMov = zMov;
353 
354  //================================================ Save results as vector
355  ProSHADE_internal_misc::addToDoubleVector ( &ret, static_cast<proshade_double> ( -xMov ) );
356  ProSHADE_internal_misc::addToDoubleVector ( &ret, static_cast<proshade_double> ( -yMov ) );
357  ProSHADE_internal_misc::addToDoubleVector ( &ret, static_cast<proshade_double> ( -zMov ) );
358 
359  //================================================ Save original from variables for PDB output
360  this->mapMovFromsChangeX = static_cast<proshade_double> ( this->xFrom );
361  this->mapMovFromsChangeY = static_cast<proshade_double> ( this->yFrom );
362  this->mapMovFromsChangeZ = static_cast<proshade_double> ( this->zFrom );
363 
364  //================================================ Move the map
365  ProSHADE_internal_mapManip::moveMapByIndices ( &modXMov, &modYMov, &modZMov, this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
366  this->getXFromPtr(), this->getXToPtr(),
367  this->getYFromPtr(), this->getYToPtr(),
368  this->getZFromPtr(), this->getZToPtr(),
369  this->getXAxisOrigin(), this->getYAxisOrigin(), this->getZAxisOrigin() );
370 
371  ProSHADE_internal_mapManip::moveMapByFourier ( this->getInternalMap(), modXMov, modYMov, modZMov,
372  this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
373  static_cast< proshade_signed > ( this->getXDim() ), static_cast< proshade_signed > ( this->getYDim() ), static_cast< proshade_signed > ( this->getZDim() ) );
374 
375  //================================================ Keep only the change in from and to variables
376  this->mapMovFromsChangeX = static_cast<proshade_double> ( this->xFrom ) - this->mapMovFromsChangeX;
377  this->mapMovFromsChangeY = static_cast<proshade_double> ( this->yFrom ) - this->mapMovFromsChangeY;
378  this->mapMovFromsChangeZ = static_cast<proshade_double> ( this->zFrom ) - this->mapMovFromsChangeZ;
379 
380  //================================================ Compute the optimal rotation centre for co-ordinates
381  this->computePdbRotationCentre ( );
382  this->computeOptimalTranslation ( eulA, eulB, eulG, -xMov, -yMov, -zMov );
383 
384  //================================================ Done
385  return ( ret );
386 
387 }
388 
399 {
400  //================================================ Do this using Fourier!
401  fftw_complex *tmpIn1 = nullptr, *tmpOut1 = nullptr, *tmpIn2 = nullptr, *tmpOut2 = nullptr, *resOut = nullptr;
402  fftw_plan forwardFourierObj1, forwardFourierObj2, inverseFourierCombo;
403  proshade_unsign dimMult = staticStructure->getXDim() * staticStructure->getYDim() * staticStructure->getZDim();
404  ProSHADE_internal_overlay::allocateTranslationFunctionMemory ( tmpIn1, tmpOut1, tmpIn2, tmpOut2, this->translationMap, resOut, forwardFourierObj1, forwardFourierObj2, inverseFourierCombo, staticStructure->getXDim(), staticStructure->getYDim(), staticStructure->getZDim() );
405 
406  //================================================ Fill in input data
407  for ( proshade_unsign iter = 0; iter < dimMult; iter++ ) { tmpIn1[iter][0] = staticStructure->getMapValue ( iter ); tmpIn1[iter][1] = 0.0; }
408  for ( proshade_unsign iter = 0; iter < dimMult; iter++ ) { tmpIn2[iter][0] = this->getMapValue ( iter ); tmpIn2[iter][1] = 0.0; }
409 
410  //================================================ Calculate Fourier
411  fftw_execute ( forwardFourierObj1 );
412  fftw_execute ( forwardFourierObj2 );
413 
414  //================================================ Combine Fourier coeffs and invert
415  ProSHADE_internal_overlay::combineFourierForTranslation ( tmpOut1, tmpOut2, resOut, staticStructure->getXDim(), staticStructure->getYDim(), staticStructure->getZDim() );
416  fftw_execute ( inverseFourierCombo );
417 
418  //================================================ Free memory
419  ProSHADE_internal_overlay::freeTranslationFunctionMemory ( tmpIn1, tmpOut1, tmpIn2, tmpOut2, resOut, forwardFourierObj1, forwardFourierObj2, inverseFourierCombo );
420 
421  //================================================ Done
422  return ;
423 
424 }
425 
441 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 )
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 }
468 
480 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 )
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 }
496 
506 void ProSHADE_internal_overlay::combineFourierForTranslation ( fftw_complex* tmpOut1, fftw_complex* tmpOut2, fftw_complex*& resOut, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD )
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 }
566 
578 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 )
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 }
607 
618 void ProSHADE_internal_data::ProSHADE_data::zeroPaddToDims ( proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim )
619 {
620  //================================================ Sanity check
621  if ( ( this->xDimIndices > xDim ) || ( this->yDimIndices > yDim ) || ( this->zDimIndices > zDim ) )
622  {
623  throw ProSHADE_exception ( "Cannot zero-pad in negative direction.", "EO00034", __FILE__, __LINE__, __func__, "The requested padded size of a structure is smaller than\n : the current size. If the user sees this error, there is\n : likely a considerable bug. Please report this error." );
624  }
625 
626  //================================================ If done, do nothing
627  if ( ( this->xDimIndices == xDim ) && ( this->yDimIndices == yDim ) && ( this->zDimIndices == zDim ) ) { return ; }
628 
629  //================================================ Find out how many zeroes need to be added before and after
630  proshade_unsign addXPre, addYPre, addZPre, addXPost, addYPost, addZPost;
631  ProSHADE_internal_overlay::computeBeforeAfterZeroCounts ( &addXPre, &addYPre, &addZPre, &addXPost, &addYPost, &addZPost, xDim, yDim, zDim, this->xDimIndices, this->yDimIndices, this->zDimIndices );
632 
633  //================================================ Create a new map
634  proshade_double* newMap = new proshade_double [xDim * yDim * zDim];
635 
636  //================================================ Do the hard work
637  ProSHADE_internal_overlay::paddMapWithZeroes ( this->internalMap, newMap, xDim, yDim, zDim, this->xDimIndices, this->yDimIndices, this->zDimIndices, addXPre, addYPre, addZPre );
638 
639  //================================================ Create a new internal map and copy
640  delete[] this->internalMap;
641  this->internalMap = new proshade_double [xDim * yDim * zDim];
642  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( xDim * yDim * zDim ); iter++ ) { this->internalMap[iter] = newMap[iter]; }
643 
644  //================================================ Release memory
645  delete[] newMap;
646 
647  //================================================ Update map related variables
648  this->xDimSize = static_cast< proshade_single > ( xDim ) * ( this->xDimSize / static_cast< proshade_single > ( this->xDimIndices ) );
649  this->yDimSize = static_cast< proshade_single > ( yDim ) * ( this->yDimSize / static_cast< proshade_single > ( this->yDimIndices ) );
650  this->zDimSize = static_cast< proshade_single > ( zDim ) * ( this->zDimSize / static_cast< proshade_single > ( this->zDimIndices ) );
651  this->xDimIndices = xDim ; this->yDimIndices = yDim ; this->zDimIndices = zDim;
652  this->xGridIndices = xDim ; this->yGridIndices = yDim ; this->zGridIndices = zDim;
653  this->xFrom -= addXPre ; this->yFrom -= addYPre ; this->zFrom -= addZPre;
654  this->xTo += addXPost; this->yTo += addYPost; this->zTo += addZPost;
655  this->xAxisOrigin -= addXPre ; this->yAxisOrigin -= addYPre ; this->zAxisOrigin -= addZPre ;
656 
657  //================================================ Done
658  return ;
659 
660 }
661 
677 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 )
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 }
691 
706 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 )
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 }
743 
756 void ProSHADE_internal_data::ProSHADE_data::rotateMapReciprocalSpace ( ProSHADE_settings* settings, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma )
757 {
758  //================================================ Set maximum comparison bandwidth to maximum object bandwidth
759  this->maxCompBand = this->spheres[this->noSpheres-1]->getLocalBandwidth();
760 
761  //================================================ Save map COM after processing but before rotation
762  this->findMapCOM ( );
763  this->mapCOMProcessChangeX += ( this->xCom - this->originalMapXCom );
764  this->mapCOMProcessChangeY += ( this->yCom - this->originalMapYCom );
765  this->mapCOMProcessChangeZ += ( this->zCom - this->originalMapZCom );
766 
767  //================================================ Compute the Wigner D matrices for the Euler angles
768  ProSHADE_internal_wigner::computeWignerMatricesForRotation ( settings, this, -eulerAlpha, eulerBeta, -eulerGamma );
769 
770  //================================================ Initialise rotated Spherical Harmonics memory
771  this->allocateRotatedSHMemory ( );
772 
773  //================================================ Multiply SH coeffs by Wigner
774  this->computeRotatedSH ( );
775 
776  //================================================ Inverse the SH coeffs to shells
777  this->invertSHCoefficients ( );
778 
779  //================================================ Find spherical cut-offs
780  std::vector<proshade_double> lonCO, latCO;
781  ProSHADE_internal_overlay::computeAngularThreshold ( &lonCO, &latCO, settings->maxBandwidth * 2 );
782 
783  //================================================ Allocate memory for the rotated map
784  proshade_double *densityMapRotated = new proshade_double [this->xDimIndices * this->yDimIndices * this->zDimIndices];
785  ProSHADE_internal_misc::checkMemoryAllocation ( densityMapRotated, __FILE__, __LINE__, __func__ );
786  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( this->xDimIndices * this->yDimIndices * this->zDimIndices ); iter++ ) { densityMapRotated[iter] = 0.0; }
787 
788  //================================================ Interpolate onto cartesian grid
789  this->interpolateMapFromSpheres ( densityMapRotated );
790 
791  //================================================ Copy map
792  for ( proshade_unsign iter = 0; iter < ( this->xDimIndices * this->yDimIndices * this->zDimIndices ); iter++ )
793  {
794  this->internalMap[iter] = densityMapRotated[iter];
795  }
796 
797  //================================================ Release rotated map (original is now rotated)
798  delete[] densityMapRotated;
799 
800  //================================================ Done
801  return ;
802 
803 }
804 
816 void ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpace ( proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, proshade_double*& map )
817 {
818  //================================================ Initialise local variables
819  bool withinBounds = true;
820  proshade_double c000, c001, c010, c011, c100, c101, c110, c111, c00, c01, c10, c11, c0, c1;
821  size_t arrPos = 0;
822  proshade_double xCOM, yCOM, zCOM;
823 
824  //================================================ Store sampling rates
825  proshade_single xSampRate = this->xDimSize / static_cast< proshade_single > ( this->xTo - this->xFrom );
826  proshade_single ySampRate = this->yDimSize / static_cast< proshade_single > ( this->yTo - this->yFrom );
827  proshade_single zSampRate = this->zDimSize / static_cast< proshade_single > ( this->zTo - this->zFrom );
828 
829  //================================================ Compute map COM
830  ProSHADE_internal_mapManip::findMAPCOMValues ( this->internalMap, &xCOM, &yCOM, &zCOM, this->xDimSize, this->yDimSize, this->zDimSize, this->xFrom, this->xTo, this->yFrom, this->yTo, this->zFrom, this->zTo );
831 
832  //================================================ Allocate local variables
833  proshade_single *mins = new proshade_single[3];
834  proshade_single *maxs = new proshade_single[3];
835  proshade_single *rotMat = new proshade_single[9];
836  proshade_single *rotVec;
837  proshade_single *interpMins = new proshade_single[3];
838  proshade_single *interpMaxs = new proshade_single[3];
839  proshade_single *interpDiff = new proshade_single[3];
840  proshade_single *movs = new proshade_single[3];
841  map = new proshade_double[ this->xDimIndices * this->yDimIndices * this->zDimIndices ];
842 
843  //================================================ Check memory allocation
844  ProSHADE_internal_misc::checkMemoryAllocation ( mins, __FILE__, __LINE__, __func__ );
845  ProSHADE_internal_misc::checkMemoryAllocation ( maxs, __FILE__, __LINE__, __func__ );
846  ProSHADE_internal_misc::checkMemoryAllocation ( rotMat, __FILE__, __LINE__, __func__ );
847  ProSHADE_internal_misc::checkMemoryAllocation ( interpMins, __FILE__, __LINE__, __func__ );
848  ProSHADE_internal_misc::checkMemoryAllocation ( interpMaxs, __FILE__, __LINE__, __func__ );
849  ProSHADE_internal_misc::checkMemoryAllocation ( interpDiff, __FILE__, __LINE__, __func__ );
850  ProSHADE_internal_misc::checkMemoryAllocation ( movs, __FILE__, __LINE__, __func__ );
851  ProSHADE_internal_misc::checkMemoryAllocation ( map, __FILE__, __LINE__, __func__ );
852 
853  //================================================ Fill map with zeroes
854  for ( size_t iter = 0; iter < static_cast< size_t > ( this->xDimIndices * this->yDimIndices * this->zDimIndices ); iter++ ) { map[iter] = 0.0; }
855 
856  //================================================ Determine map max's and min's in terms of the hkl system
857  mins[0] = std::floor ( static_cast< proshade_single > ( this->xDimIndices ) / -2.0f );
858  mins[1] = std::floor ( static_cast< proshade_single > ( this->yDimIndices ) / -2.0f );
859  mins[2] = std::floor ( static_cast< proshade_single > ( this->zDimIndices ) / -2.0f );
860 
861  maxs[0] = -mins[0];
862  maxs[1] = -mins[1];
863  maxs[2] = -mins[2];
864 
865  if ( this->xDimIndices % 2 == 0 ) { maxs[0] -= 1.0f; }
866  if ( this->yDimIndices % 2 == 0 ) { maxs[1] -= 1.0f; }
867  if ( this->zDimIndices % 2 == 0 ) { maxs[2] -= 1.0f; }
868 
869  //================================================ Rotate about COM instead of map midpoint
870  movs[0] = ( static_cast< proshade_single > ( xCOM ) / xSampRate ) + ( mins[0] - static_cast< proshade_single > ( this->xFrom ) );
871  movs[1] = ( static_cast< proshade_single > ( yCOM ) / ySampRate ) + ( mins[1] - static_cast< proshade_single > ( this->yFrom ) );
872  movs[2] = ( static_cast< proshade_single > ( zCOM ) / zSampRate ) + ( mins[2] - static_cast< proshade_single > ( this->zFrom ) );
873 
874  //================================================ Get rotation matrix from Euler angles
875  ProSHADE_internal_maths::getRotationMatrixFromAngleAxis ( rotMat, axX, axY, axZ, axAng );
876 
877  //================================================ For each point
878  for ( proshade_single xIt = mins[0]; xIt <= maxs[0]; xIt += 1.0f )
879  {
880  for ( proshade_single yIt = mins[1]; yIt <= maxs[1]; yIt += 1.0f )
881  {
882  for ( proshade_single zIt = mins[2]; zIt <= maxs[2]; zIt += 1.0f )
883  {
884  //==================================== Compute new point position
885  rotVec = ProSHADE_internal_maths::compute3x3MatrixVectorMultiplication ( rotMat, xIt - movs[0], yIt - movs[1], zIt - movs[2] );
886 
887  //==================================== Find surrounding grid points indices and check for boundaries
888  withinBounds = true;
889  for ( size_t posIt = 0; posIt < 3; posIt++ )
890  {
891  //================================ Determine surrounding points indices in this dimension
892  interpMins[posIt] = std::floor ( rotVec[posIt] );
893  interpMaxs[posIt] = interpMins[posIt] + 1.0f;
894 
895  //================================ Check for boundaries
896  if ( ( maxs[posIt] < interpMins[posIt] ) || ( interpMins[posIt] < mins[posIt] ) || ( maxs[posIt] < interpMaxs[posIt] ) || ( interpMaxs[posIt] < mins[posIt] ) )
897  {
898  withinBounds = false;
899  break;
900  }
901 
902  //================================ Compute the difference from position to min index along this axis
903  interpDiff[posIt] = rotVec[posIt] - interpMins[posIt];
904  }
905  if ( !withinBounds ) { continue; }
906 
907  //==================================== Make sure interpolation max's are within bounds
908  for ( size_t posIt = 0; posIt < 3; posIt++ )
909  {
910  interpMaxs[posIt] = std::min ( maxs[posIt], std::max ( mins[posIt], interpMaxs[posIt] ) );
911  }
912 
913  //==================================== Release memory
914  delete[] rotVec;
915 
916  //==================================== Find the surrounding points values from their indices
917  arrPos = static_cast< size_t > ( ( interpMins[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMins[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMins[0] - mins[0] ) ) );
918  c000 = this->internalMap[arrPos];
919 
920  arrPos = static_cast< size_t > ( ( interpMaxs[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMins[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMins[0] - mins[0] ) ) );
921  c001 = this->internalMap[arrPos];
922 
923  arrPos = static_cast< size_t > ( ( interpMins[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMaxs[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMins[0] - mins[0] ) ) );
924  c010 = this->internalMap[arrPos];
925 
926  arrPos = static_cast< size_t > ( ( interpMaxs[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMaxs[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMins[0] - mins[0] ) ) );
927  c011 = this->internalMap[arrPos];
928 
929  arrPos = static_cast< size_t > ( ( interpMins[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMins[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMaxs[0] - mins[0] ) ) );
930  c100 = this->internalMap[arrPos];
931 
932  arrPos = static_cast< size_t > ( ( interpMaxs[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMins[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMaxs[0] - mins[0] ) ) );
933  c101 = this->internalMap[arrPos];
934 
935  arrPos = static_cast< size_t > ( ( interpMins[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMaxs[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMaxs[0] - mins[0] ) ) );
936  c110 = this->internalMap[arrPos];
937 
938  arrPos = static_cast< size_t > ( ( interpMaxs[2] - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( interpMaxs[1] - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( interpMaxs[0] - mins[0] ) ) );
939  c111 = this->internalMap[arrPos];
940 
941  //==================================== Interpolate along x-axis
942  c00 = ( c000 * ( 1.0 - static_cast< proshade_double > ( interpDiff[0] ) ) ) + ( c100 * static_cast< proshade_double > ( interpDiff[0] ) );
943  c01 = ( c001 * ( 1.0 - static_cast< proshade_double > ( interpDiff[0] ) ) ) + ( c101 * static_cast< proshade_double > ( interpDiff[0] ) );
944  c10 = ( c010 * ( 1.0 - static_cast< proshade_double > ( interpDiff[0] ) ) ) + ( c110 * static_cast< proshade_double > ( interpDiff[0] ) );
945  c11 = ( c011 * ( 1.0 - static_cast< proshade_double > ( interpDiff[0] ) ) ) + ( c111 * static_cast< proshade_double > ( interpDiff[0] ) );
946 
947  //==================================== Interpolate along y-axis
948  c0 = ( c00 * ( 1.0 - static_cast< proshade_double > ( interpDiff[1] ) ) ) + ( c10 * static_cast< proshade_double > ( interpDiff[1] ) );
949  c1 = ( c01 * ( 1.0 - static_cast< proshade_double > ( interpDiff[1] ) ) ) + ( c11 * static_cast< proshade_double > ( interpDiff[1] ) );
950 
951  //==================================== Interpolate along z-axis
952  arrPos = static_cast< size_t > ( ( zIt - mins[2] ) + static_cast< proshade_single > ( this->zDimIndices ) * ( ( yIt - mins[1] ) + static_cast< proshade_single > ( this->yDimIndices ) * ( xIt - mins[0] ) ) );
953  map[arrPos] = ( c0 * ( 1.0 - static_cast< proshade_double > ( interpDiff[2] ) ) ) + ( c1 * static_cast< proshade_double > ( interpDiff[2] ) );
954  }
955  }
956  }
957 
958  //================================================ Release memory
959  delete[] mins;
960  delete[] maxs;
961  delete[] rotMat;
962  delete[] interpMins;
963  delete[] interpMaxs;
964  delete[] interpDiff;
965  delete[] movs;
966 
967  //================================================ Done
968  return ;
969 
970 }
971 
982 void ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpaceInPlace ( proshade_double eulA, proshade_double eulB, proshade_double eulG )
983 {
984  //================================================ Initialise local variables
985  proshade_double axX, axY, axZ, axAng, tmp, *rMat, *map;
986 
987  //================================================ Allocate local memory
988  rMat = new proshade_double[9];
989 
990  //================================================ Check local memory allocation
991  ProSHADE_internal_misc::checkMemoryAllocation ( rMat, __FILE__, __LINE__, __func__ );
992 
993  //================================================ Convert Euler angles to rotation matrix
995 
996  //================================================ Transpose the rotation matrix
997  tmp = rMat[1];
998  rMat[1] = rMat[3];
999  rMat[3] = tmp;
1000 
1001  tmp = rMat[2];
1002  rMat[2] = rMat[6];
1003  rMat[6] = tmp;
1004 
1005  tmp = rMat[5];
1006  rMat[5] = rMat[7];
1007  rMat[7] = tmp;
1008 
1009  //================================================ Convert rotation matrix to angle-axis
1010  ProSHADE_internal_maths::getAxisAngleFromRotationMatrix ( rMat, &axX, &axY, &axZ, &axAng );
1011 
1012  //================================================ Rotate the internal map
1013  this->rotateMapRealSpace ( axX, axY, axZ, axAng, map );
1014 
1015  //================================================ Copy the rotated map in place of the internal map
1016  for ( size_t iter = 0; iter < static_cast< size_t > ( this->xDimIndices * this->yDimIndices * this->zDimIndices ); iter++ )
1017  {
1018  this->internalMap[iter] = map[iter];
1019  }
1020 
1021  //================================================ Release local memory
1022  delete[] rMat;
1023  delete[] map;
1024 
1025  //================================================ Done
1026  return ;
1027 
1028 }
1029 
1040 void ProSHADE_internal_data::ProSHADE_data::translateMap ( proshade_double trsX, proshade_double trsY, proshade_double trsZ )
1041 {
1042  //================================================ Initialise local variables
1043  proshade_single xMov = static_cast< proshade_single > ( -trsX );
1044  proshade_single yMov = static_cast< proshade_single > ( -trsY );
1045  proshade_single zMov = static_cast< proshade_single > ( -trsZ );
1046 
1047  //================================================ Move the whole map frame to minimise the Fourier movement
1048  ProSHADE_internal_mapManip::moveMapByIndices ( &xMov, &yMov, &zMov, this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
1049  this->getXFromPtr(), this->getXToPtr(), this->getYFromPtr(), this->getYToPtr(),
1050  this->getZFromPtr(), this->getZToPtr(), this->getXAxisOrigin(), this->getYAxisOrigin(), this->getZAxisOrigin() );
1051 
1052  //================================================ Finalise the movement by in-frame Fourier movement
1053  ProSHADE_internal_mapManip::moveMapByFourier ( this->getInternalMap(), xMov, yMov, zMov, this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
1054  static_cast< proshade_signed > ( this->getXDim() ), static_cast< proshade_signed > ( this->getYDim() ),
1055  static_cast< proshade_signed > ( this->getZDim() ) );
1056 
1057  //================================================ Done
1058  return ;
1059 
1060 }
1061 
1065 {
1066  //================================================ Allocate the main pointer and check
1067  this->rotSphericalHarmonics = new proshade_complex* [this->noSpheres];
1068  ProSHADE_internal_misc::checkMemoryAllocation ( this->rotSphericalHarmonics, __FILE__, __LINE__, __func__ );
1069 
1070  //================================================ For each sphere
1071  for ( proshade_unsign iter = 0; iter < this->noSpheres; iter++ )
1072  {
1073  //============================================ Allocate the sphere storage space
1074  this->rotSphericalHarmonics[iter] = new proshade_complex [static_cast<proshade_unsign> ( pow ( (this->spheres[iter]->getLocalBandwidth() * 2), 2) )];
1075  ProSHADE_internal_misc::checkMemoryAllocation ( this->rotSphericalHarmonics[iter], __FILE__, __LINE__, __func__ );
1076 
1077  //============================================ Set values to zeroes
1078  for ( proshade_unsign it = 0; it < static_cast<proshade_unsign> ( pow ( (this->spheres[iter]->getLocalBandwidth() * 2), 2) ); it++ )
1079  {
1080  this->rotSphericalHarmonics[iter][it][0] = 0.0;
1081  this->rotSphericalHarmonics[iter][it][1] = 0.0;
1082  }
1083  }
1084 
1085  //================================================ Done
1086  return ;
1087 
1088 }
1089 
1093 {
1094  //================================================ Initialise variables
1095  proshade_double WigDR, WigDI, *ShR, *ShI, retR, retI;
1096  proshade_unsign arrPos;
1097 
1098  //================================================ Compute
1099  for ( proshade_signed shell = 0; shell < static_cast<proshade_signed> ( this->noSpheres ); shell++ )
1100  {
1101  //============================================ For each band
1102  for ( proshade_signed bandIter = 0; bandIter < static_cast<proshade_signed> ( this->spheres[shell]->getLocalBandwidth() ); bandIter++ )
1103  {
1104  //======================================== For each order1
1105  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
1106  {
1107  //==================================== Get Spherical Harmonics value
1108  ShR = getRealSphHarmValue ( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( order1 ), static_cast< proshade_unsign > ( shell ) );
1109  ShI = getImagSphHarmValue ( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( order1 ), static_cast< proshade_unsign > ( shell ) );
1110 
1111  //==================================== For each order2
1112  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
1113  {
1114  //================================ Get Wigner D values
1115  this->getWignerMatrixValue ( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( order1 ), static_cast< proshade_unsign > ( order2 ), &WigDR, &WigDI );
1116 
1117  //================================ Multiply SH and Wigner
1118  ProSHADE_internal_maths::complexMultiplication ( ShR, ShI, &WigDR, &WigDI, &retR, &retI );
1119 
1120  //================================ Save
1121  arrPos = static_cast<proshade_unsign> ( seanindex ( static_cast< int > ( order2-bandIter ), static_cast< int > ( bandIter ),
1122  static_cast< int > ( this->spheres[shell]->getLocalBandwidth() ) ) );
1123  this->rotSphericalHarmonics[shell][arrPos][0] += retR;
1124  this->rotSphericalHarmonics[shell][arrPos][1] += retI;
1125  }
1126  }
1127  }
1128  }
1129 
1130  //================================================ Done
1131  return ;
1132 
1133 }
1134 
1147 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 )
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 }
1189 
1194 {
1195  //================================================ Initialise local variables
1196  double *sigR = nullptr, *sigI = nullptr, *rcoeffs = nullptr, *icoeffs = nullptr, *weights = nullptr, *workspace = nullptr;
1197  fftw_plan idctPlan, ifftPlan;
1198 
1199  //================================================ For each shell
1200  for ( int shell = 0; shell < static_cast<int> ( this->noSpheres ); shell++ )
1201  {
1202  //=========================================== Initialise internal variables
1203  proshade_unsign oneDim = this->spheres[shell]->getLocalBandwidth() * 2;
1204 
1205  //=========================================== Allocate memory
1206  ProSHADE_internal_overlay::initialiseInverseSHComputation ( this->spheres[shell]->getLocalBandwidth(), sigR, sigI, rcoeffs, icoeffs, weights, workspace, idctPlan, ifftPlan );
1207 
1208  //=========================================== Compute weights for the transform using the appropriate shell related band
1209  makeweights ( static_cast< int > ( this->spheres[shell]->getLocalBandwidth() ), weights );
1210 
1211  //============================================ Allocate rotated shell mapped data memory
1212  this->spheres[shell]->allocateRotatedMap ( );
1213 
1214  //============================================ Load SH coeffs to arrays
1215  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( oneDim * oneDim ); iter++ )
1216  {
1217  rcoeffs[iter] = this->rotSphericalHarmonics[shell][iter][0];
1218  icoeffs[iter] = this->rotSphericalHarmonics[shell][iter][1];
1219  sigR[iter] = 0.0;
1220  sigI[iter] = 0.0;
1221  }
1222 
1223  //============================================ Get inverse spherical harmonics transform for the shell
1224  InvFST_semi_fly ( rcoeffs,
1225  icoeffs,
1226  sigR,
1227  sigI,
1228  static_cast< int > ( this->spheres[shell]->getLocalBandwidth() ),
1229  workspace,
1230  0,
1231  static_cast< int > ( this->spheres[shell]->getLocalBandwidth() ),
1232  &idctPlan,
1233  &ifftPlan );
1234 
1235  //=========================================== Copy the results to the rotated shells array
1236  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( oneDim * oneDim ); iter++ )
1237  {
1238  this->spheres[shell]->setRotatedMappedData ( iter, sigR[iter] );
1239  }
1240 
1241  //=========================================== Release the plans
1242  fftw_destroy_plan ( idctPlan );
1243  fftw_destroy_plan ( ifftPlan );
1244 
1245  //=========================================== Release the memory
1246  delete[] sigR;
1247  delete[] rcoeffs;
1248  delete[] weights;
1249  delete[] workspace;
1250  }
1251 
1252  //================================================ Done
1253  return ;
1254 
1255 }
1256 
1263 void ProSHADE_internal_overlay::computeAngularThreshold ( std::vector<proshade_double>* lonCO, std::vector<proshade_double>* latCO, proshade_unsign angRes )
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 }
1281 
1286 void ProSHADE_internal_data::ProSHADE_data::interpolateMapFromSpheres ( proshade_double*& densityMapRotated )
1287 {
1288  //================================================ Initialise variables
1289  proshade_double rad = 0.0, lon = 0.0, lat = 0.0, newU = 0.0, newV = 0.0, newW = 0.0;
1290  proshade_unsign lowerLonL = 0, upperLonL = 0, lowerLonU = 0, upperLonU = 0, lowerLatL = 0, upperLatL = 0, lowerLatU = 0, upperLatU = 0, lowerShell = 0, upperShell = 0;
1291  proshade_double x00 = 0.0, x01 = 0.0, x10 = 0.0, x11 = 0.0, distLLon = 0.0, distLLat = 0.0, distLRad = 0.0, valLLon = 0.0, valULon = 0.0;
1292  proshade_double lowerShellValue = 0.0, upperShellValue = 0.0;
1293  proshade_double xSamplingRate = static_cast<proshade_double> ( this->xDimSize ) / static_cast<proshade_double> ( this->xDimIndices );
1294  proshade_double ySamplingRate = static_cast<proshade_double> ( this->yDimSize ) / static_cast<proshade_double> ( this->yDimIndices );
1295  proshade_double zSamplingRate = static_cast<proshade_double> ( this->zDimSize ) / static_cast<proshade_double> ( this->zDimIndices );
1296  proshade_signed arrPos;
1297  std::vector<proshade_double> lonCOU, latCOU, lonCOL, latCOL;
1298 
1299  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> (this->xDimIndices); uIt++ )
1300  {
1301  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> (this->yDimIndices); vIt++ )
1302  {
1303  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> (this->zDimIndices); wIt++ )
1304  {
1305  //==================================== Convert to centered coords
1306  newU = static_cast<proshade_double> ( uIt - ( static_cast<proshade_signed> (this->xDimIndices) / 2 ) );
1307  newV = static_cast<proshade_double> ( vIt - ( static_cast<proshade_signed> (this->yDimIndices) / 2 ) );
1308  newW = static_cast<proshade_double> ( wIt - ( static_cast<proshade_signed> (this->zDimIndices) / 2 ) );
1309 
1310  //==================================== Deal with 0 ; 0 ; 0
1311  if ( ( newU == 0.0 ) && ( newV == 0.0 ) && ( newW == 0.0 ) )
1312  {
1313  arrPos = wIt + static_cast< proshade_signed > ( this->zDimIndices ) * ( vIt + static_cast< proshade_signed > ( this->yDimIndices ) * uIt );
1314  densityMapRotated[arrPos] = this->internalMap[arrPos];
1315  continue;
1316  }
1317 
1318  //==================================== Convert to spherical coords
1319  rad = sqrt ( pow( ( newU * xSamplingRate ), 2.0 ) +
1320  pow( ( newV * ySamplingRate ), 2.0 ) +
1321  pow( ( newW * zSamplingRate ), 2.0 ) );
1322  lon = atan2 ( ( newV * ySamplingRate ), ( newU * xSamplingRate ) );
1323  lat = asin ( ( newW * zSamplingRate ) / rad );
1324 
1325  //==================================== Deal with nan's
1326  if ( rad != rad ) { rad = 0.0; }
1327  if ( lon != lon ) { lon = 0.0; }
1328  if ( lat != lat ) { lat = 0.0; }
1329 
1330  //==================================== Find shells above and below
1331  lowerShell = 0;
1332  upperShell = 0;
1333  for ( proshade_unsign iter = 0; iter < (this->noSpheres-1); iter++ )
1334  {
1335  if ( ( static_cast< proshade_double > ( this->spherePos.at(iter) ) <= rad ) && ( static_cast< proshade_double > ( this->spherePos.at(iter+1) ) > rad ) )
1336  {
1337  lowerShell = iter;
1338  upperShell = iter+1;
1339  break;
1340  }
1341  }
1342 
1343  if ( upperShell == 0 )
1344  {
1345  arrPos = wIt + static_cast< proshade_signed > ( this->zDimIndices ) * ( vIt + static_cast< proshade_signed > ( this->yDimIndices ) * uIt );
1346  densityMapRotated[arrPos] = 0.0;
1347  continue;
1348  }
1349 
1350  //==================================== Get the longitude and lattitude cut-offs for this shell resolution
1351  lonCOL.clear(); latCOL.clear(); lonCOU.clear(); latCOU.clear();
1352  ProSHADE_internal_overlay::computeAngularThreshold ( &lonCOL, &latCOL, this->spheres[lowerShell]->getLocalAngRes() );
1353  ProSHADE_internal_overlay::computeAngularThreshold ( &lonCOU, &latCOU, this->spheres[upperShell]->getLocalAngRes() );
1354 
1355  //==================================== Find the angle cutoffs around the point
1356  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( lonCOL.size() ); iter++ )
1357  {
1358  if ( iter == ( static_cast<proshade_unsign> ( lonCOL.size() ) - 1 ) )
1359  {
1360  lowerLonL = 0;
1361  upperLonL = 1;
1362  break;
1363  }
1364  if ( ( std::floor(10000. * lonCOL.at(iter)) <= std::floor(10000. * lon) ) && ( std::floor(10000. * lonCOL.at(iter+1)) > std::floor(10000. * lon) ) )
1365  {
1366  lowerLonL = iter;
1367  upperLonL = iter+1;
1368  break;
1369  }
1370  }
1371  if ( upperLonL == this->spheres[lowerShell]->getLocalAngRes() ) { upperLonL = 0; }
1372 
1373  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( lonCOU.size() ); iter++ )
1374  {
1375  if ( iter == ( static_cast<proshade_unsign> ( lonCOU.size() ) - 1 ) )
1376  {
1377  lowerLonU = 0;
1378  upperLonU = 1;
1379  break;
1380  }
1381  if ( ( std::floor(10000. * lonCOU.at(iter)) <= std::floor(10000. * lon) ) && ( std::floor(10000. * lonCOU.at(iter+1)) > std::floor(10000. * lon) ) )
1382  {
1383  lowerLonU = iter;
1384  upperLonU = iter+1;
1385  break;
1386  }
1387  }
1388  if ( upperLonU == this->spheres[upperShell]->getLocalAngRes() ) { upperLonU = 0; }
1389 
1390  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( latCOL.size() ); iter++ )
1391  {
1392  if ( iter == ( static_cast<proshade_unsign> ( latCOL.size() ) - 1 ) )
1393  {
1394  lowerLatL = 0;
1395  upperLatL = 1;
1396  break;
1397  }
1398  if ( ( std::floor(10000. * latCOL.at(iter)) <= std::floor(10000. * lat) ) && ( std::floor(10000. * latCOL.at(iter+1)) > std::floor(10000. * lat) ) )
1399  {
1400  lowerLatL = iter;
1401  upperLatL = iter+1;
1402  break;
1403  }
1404  }
1405  if ( upperLatL == this->spheres[lowerShell]->getLocalAngRes() ) { upperLatL = 0; }
1406 
1407  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( latCOU.size() ); iter++ )
1408  {
1409  if ( iter == ( static_cast<proshade_unsign> ( latCOU.size() ) - 1 ) )
1410  {
1411  lowerLatU = 0;
1412  upperLatU = 1;
1413  break;
1414  }
1415  if ( ( std::floor(10000. * latCOU.at(iter)) <= std::floor(10000. * lat) ) && ( std::floor(10000. * latCOU.at(iter+1)) > std::floor(10000. * lat) ) )
1416  {
1417  lowerLatU = iter;
1418  upperLatU = iter+1;
1419  break;
1420  }
1421  }
1422  if ( upperLatU == this->spheres[upperShell]->getLocalAngRes() ) { upperLatU = 0; }
1423 
1424  //==================================== Interpolate lower shell
1425  x00 = this->spheres[lowerShell]->getRotatedMappedData ( lowerLatL * this->spheres[lowerShell]->getLocalAngRes() + lowerLonL );
1426  x01 = this->spheres[lowerShell]->getRotatedMappedData ( lowerLatL * this->spheres[lowerShell]->getLocalAngRes() + upperLonL );
1427  x10 = this->spheres[lowerShell]->getRotatedMappedData ( upperLatL * this->spheres[lowerShell]->getLocalAngRes() + lowerLonL );
1428  x11 = this->spheres[lowerShell]->getRotatedMappedData ( upperLatL * this->spheres[lowerShell]->getLocalAngRes() + upperLonL );
1429 
1430  distLLon = std::abs ( lon - lonCOL.at(lowerLonL) ) / ( std::abs( lon - lonCOL.at(lowerLonL) ) + std::abs( lon - lonCOL.at(upperLonL) ) );
1431  valLLon = ( ( 1.0 - distLLon ) * x00 ) + ( distLLon * x01 );
1432  valULon = ( ( 1.0 - distLLon ) * x10 ) + ( distLLon * x11 );
1433 
1434  distLLat = std::abs ( lat - latCOL.at(lowerLatL) ) / ( std::abs( lat - latCOL.at(lowerLatL) ) + std::abs( lat - latCOL.at(upperLatL) ) );
1435  lowerShellValue = ( ( 1.0 - distLLat ) * valLLon ) + ( distLLat * valULon );
1436 
1437  //==================================== Interpolate upper shell
1438  x00 = this->spheres[upperShell]->getRotatedMappedData ( lowerLatU * this->spheres[upperShell]->getLocalAngRes() + lowerLonU );
1439  x01 = this->spheres[upperShell]->getRotatedMappedData ( lowerLatU * this->spheres[upperShell]->getLocalAngRes() + upperLonU );
1440  x10 = this->spheres[upperShell]->getRotatedMappedData ( upperLatU * this->spheres[upperShell]->getLocalAngRes() + lowerLonU );
1441  x11 = this->spheres[upperShell]->getRotatedMappedData ( upperLatU * this->spheres[upperShell]->getLocalAngRes() + upperLonU );
1442 
1443  distLLon = std::abs ( lon - lonCOU.at(lowerLonU) ) / ( std::abs( lon - lonCOU.at(lowerLonU) ) + std::abs( lon - lonCOU.at(upperLonU) ) );
1444  valLLon = ( ( 1.0 - distLLon ) * x00 ) + ( distLLon * x01 );
1445  valULon = ( ( 1.0 - distLLon ) * x10 ) + ( distLLon * x11 );
1446 
1447  distLLat = std::abs ( lat - latCOU.at(lowerLatU) ) / ( std::abs( lat - latCOU.at(lowerLatU) ) + std::abs( lat - latCOU.at(upperLatU) ) );
1448  upperShellValue = ( ( 1.0 - distLLat ) * valLLon ) + ( distLLat * valULon );
1449 
1450  //==================================== Interpolate between shells
1451  distLRad = std::abs ( rad - static_cast< proshade_double > ( this->spherePos.at(lowerShell) ) ) / ( std::abs( rad - static_cast< proshade_double > ( this->spherePos.at(lowerShell) ) ) +
1452  std::abs( rad - static_cast< proshade_double > ( this->spherePos.at(upperShell) ) ) );
1453 
1454  arrPos = wIt + static_cast< proshade_signed > ( this->zDimIndices ) * ( vIt + static_cast< proshade_signed > ( this->yDimIndices ) * uIt );
1455  densityMapRotated[arrPos] = ( ( 1.0 - distLRad ) * lowerShellValue ) + ( distLRad * upperShellValue );
1456  }
1457 
1458  }
1459 
1460  }
1461 
1462  //================================================ Done
1463  return ;
1464 
1465 }
1466 
1473 {
1474  //================================================ Initialise local variables
1475  std::vector < proshade_double > ret;
1476  proshade_double eulA, eulB, eulG;
1477 
1478  //================================================ Get inverse SO(3) map top peak Euler angle values
1480  this->getMaxBand() * 2,
1481  &eulA, &eulB, &eulG, settings );
1482 
1483  //================================================ Re-format to vector
1487 
1488  //================================================ Done
1489  return ( ret );
1490 
1491 }
ProSHADE_internal_distances::normaliseEMatrices
void normaliseEMatrices(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function normalises the E matrices.
Definition: ProSHADE_distances.cpp:572
ProSHADE_settings::maxBandwidth
proshade_unsign maxBandwidth
The bandwidth of spherical harmonics decomposition for the largest sphere.
Definition: ProSHADE_settings.hpp:58
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::allocateRotatedSHMemory
void allocateRotatedSHMemory(void)
This function allocates the memory required for storing the rotated Spherical Harmonics coefficients.
Definition: ProSHADE_overlay.cpp:1064
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::getBestTranslationMapPeaksAngstrom
std::vector< proshade_double > getBestTranslationMapPeaksAngstrom(ProSHADE_internal_data::ProSHADE_data *staticStructure, proshade_double eulA, proshade_double eulB, proshade_double eulG)
This function gets the optimal translation vector and returns it as a standard library vector....
Definition: ProSHADE_overlay.cpp:311
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::interpolateMapFromSpheres
void interpolateMapFromSpheres(proshade_double *&densityMapRotated)
This function interpolates the density map from the sphere mapped data.
Definition: ProSHADE_overlay.cpp:1286
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_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
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
This class contains all inputed and derived data for a single structure.
Definition: ProSHADE_data.hpp:49
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::findMAPCOMValues
void findMAPCOMValues(proshade_double *map, proshade_double *xCom, proshade_double *yCom, proshade_double *zCom, 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)
This function finds the Centre of Mass for a map.
Definition: ProSHADE_mapManip.cpp:240
ProSHADE_internal_data::ProSHADE_data::getBestRotationMapPeaksEulerAngles
std::vector< proshade_double > getBestRotationMapPeaksEulerAngles(ProSHADE_settings *settings)
This function returns a vector of three floats, the three Euler angles of the best peak in the rotati...
Definition: ProSHADE_overlay.cpp:1472
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_maths::getRotationMatrixFromEulerZXZAngles
void getRotationMatrixFromEulerZXZAngles(proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma, proshade_double *matrix)
Function to find the rotation matrix from Euler angles (ZXZ convention).
Definition: ProSHADE_maths.cpp:1011
ProSHADE_internal_overlay::computeBeforeAfterZeroCounts
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 dimens...
Definition: ProSHADE_overlay.cpp:677
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_maths::complexMultiplication
void complexMultiplication(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2, proshade_double *retReal, proshade_double *retImag)
Function to multiply two complex numbers.
Definition: ProSHADE_maths.cpp:38
ProSHADE_internal_overlay::getOptimalRotation
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.
Definition: ProSHADE_overlay.cpp:154
ProSHADE_internal_overlay::freeTranslationFunctionMemory
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 comput...
Definition: ProSHADE_overlay.cpp:480
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::rotateMapReciprocalSpace
void rotateMapReciprocalSpace(ProSHADE_settings *settings, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function rotates a map based on the given Euler angles.
Definition: ProSHADE_overlay.cpp:756
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_wigner::computeWignerMatricesForRotation
void computeWignerMatricesForRotation(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *obj, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function computes the Wigner D matrices for a particular set of Euler angles.
Definition: ProSHADE_wignerMatrices.cpp:258
ProSHADE_overlay.hpp
This header file declares the functions required for the structure overlay computation.
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::computeRotatedSH
void computeRotatedSH(void)
This function multiplies the objects spherical harmonics with the Wigner D matrices,...
Definition: ProSHADE_overlay.cpp:1092
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_maths::getRotationMatrixFromAngleAxis
void getRotationMatrixFromAngleAxis(proshade_double *rotMat, proshade_double x, proshade_double y, proshade_double z, proshade_double ang)
This function converts the axis-angle representation to the rotation matrix representation.
Definition: ProSHADE_maths.cpp:1448
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::translateMap
void translateMap(proshade_double trsX, proshade_double trsY, proshade_double trsZ)
This function simply translates the map by a given number of Angstroms along the three axes.
Definition: ProSHADE_overlay.cpp:1040
ProSHADE_internal_data::ProSHADE_data::invertSHCoefficients
void invertSHCoefficients(void)
This function computes the shell mapped data from inverting the Spherical Harmonics coefficients.
Definition: ProSHADE_overlay.cpp:1193
ProSHADE_internal_overlay::getOptimalTranslation
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 objec...
Definition: ProSHADE_overlay.cpp:203
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_settings
This class stores all the settings and is passed to the executive classes instead of a multitude of p...
Definition: ProSHADE_settings.hpp:37
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_distances::computeEMatrices
void computeEMatrices(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the complete E matrices and their weights between any two objects.
Definition: ProSHADE_distances.cpp:510
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::computeAngularThreshold
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.
Definition: ProSHADE_overlay.cpp:1263
ProSHADE_internal_overlay::initialiseInverseSHComputation
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.
Definition: ProSHADE_overlay.cpp:1147
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_distances::computeInverseSOFTTransform
void computeInverseSOFTTransform(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the inverse SO(3) transform.
Definition: ProSHADE_distances.cpp:854
ProSHADE_internal_data::ProSHADE_data::getMapValue
proshade_double getMapValue(proshade_unsign pos)
This function returns the internal map representation value of a particular array position.
Definition: ProSHADE_data.cpp:3620
ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpace
void rotateMapRealSpace(proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, proshade_double *&map)
This function rotates a map based on the given angle-axis rotation.
Definition: ProSHADE_overlay.cpp:816
ProSHADE_internal_maths::getAxisAngleFromRotationMatrix
void getAxisAngleFromRotationMatrix(proshade_double *rotMat, proshade_double *x, proshade_double *y, proshade_double *z, proshade_double *ang)
This function converts rotation matrix to the axis-angle representation.
Definition: ProSHADE_maths.cpp:1128
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_maths::compute3x3MatrixVectorMultiplication
proshade_double * compute3x3MatrixVectorMultiplication(proshade_double *mat, proshade_double x, proshade_double y, proshade_double z)
Function for computing a 3x3 matrix to 3x1 vector multiplication.
Definition: ProSHADE_maths.cpp:1861
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_overlay::allocateTranslationFunctionMemory
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 compu...
Definition: ProSHADE_overlay.cpp:441
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_overlay::paddMapWithZeroes
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 n...
Definition: ProSHADE_overlay.cpp:706
ProSHADE_internal_overlay::combineFourierForTranslation
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 th...
Definition: ProSHADE_overlay.cpp:506
ProSHADE_internal_distances::generateSO3CoeffsFromEMatrices
void generateSO3CoeffsFromEMatrices(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function converts the E matrices to SO(3) coefficients.
Definition: ProSHADE_distances.cpp:701
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