ProSHADE  0.7.5.4 (MAR 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 = ( this->xDimSizeOriginal / static_cast<proshade_double> ( this->xDimIndicesOriginal ) );
72  proshade_double ySamplRate = ( this->yDimSizeOriginal / static_cast<proshade_double> ( this->yDimIndicesOriginal ) );
73  proshade_double zSamplRate = ( this->zDimSizeOriginal / static_cast<proshade_double> ( 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_double trsX, proshade_double trsY, proshade_double 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 += trsX;
132  this->originalPdbTransY += trsY;
133  this->originalPdbTransZ += 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  //================================================ Compute spherical harmonics to allow Fourier space rotation
214  movingStructure->mapToSpheres ( settings );
215  movingStructure->computeSphericalHarmonics ( settings );
216 
217  //================================================ Rotate map
218  movingStructure->rotateMap ( settings, eulA, eulB, eulG );
219 
220  //================================================ Zero padding for smaller structure
221  staticStructure->zeroPaddToDims ( std::max ( staticStructure->getXDim(), movingStructure->getXDim() ),
222  std::max ( staticStructure->getYDim(), movingStructure->getYDim() ),
223  std::max ( staticStructure->getZDim(), movingStructure->getZDim() ) );
224  movingStructure->zeroPaddToDims ( std::max ( staticStructure->getXDim(), movingStructure->getXDim() ),
225  std::max ( staticStructure->getYDim(), movingStructure->getYDim() ),
226  std::max ( staticStructure->getZDim(), movingStructure->getZDim() ) );
227 
228  //================================================ Report progress
229  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting translation function computation." );
230 
231  //================================================ Compute the translation function map
232  movingStructure->computeTranslationMap ( staticStructure );
233 
234  //================================================ Find highest peak in translation function
235  proshade_double mapPeak = 0.0;
236  proshade_unsign xDimS = staticStructure->getXDim();
237  proshade_unsign yDimS = staticStructure->getYDim();
238  proshade_unsign zDimS = staticStructure->getZDim();
239  ProSHADE_internal_overlay::findHighestValueInMap ( movingStructure->getTranslationFnPointer(), xDimS, yDimS, zDimS, trsX, trsY, trsZ, &mapPeak );
240 
241  //================================================ Dont translate over half
242  if ( *trsX > ( xDimS / 2 ) ) { *trsX = *trsX - xDimS; }
243  if ( *trsY > ( yDimS / 2 ) ) { *trsY = *trsY - yDimS; }
244  if ( *trsZ > ( zDimS / 2 ) ) { *trsZ = *trsZ - zDimS; }
245 
246  //================================================ Move map
247  proshade_single xCor = ( staticStructure->xFrom - movingStructure->xFrom ) *
248  ( static_cast<proshade_double> ( staticStructure->getXDimSize() ) / staticStructure->getXDim() );
249  proshade_single yCor = ( staticStructure->yFrom - movingStructure->yFrom ) *
250  ( static_cast<proshade_double> ( staticStructure->getYDimSize() ) / staticStructure->getYDim() );
251  proshade_single zCor = ( staticStructure->zFrom - movingStructure->zFrom ) *
252  ( static_cast<proshade_double> ( staticStructure->getZDimSize() ) / staticStructure->getZDim() );
253  proshade_single xMov = staticStructure->mapCOMProcessChangeX - xCor -
254  ( *trsX * static_cast<proshade_double> ( staticStructure->getXDimSize() ) / staticStructure->getXDim() );
255  proshade_single yMov = staticStructure->mapCOMProcessChangeY - yCor -
256  ( *trsY * static_cast<proshade_double> ( staticStructure->getYDimSize() ) / staticStructure->getYDim() );
257  proshade_single zMov = staticStructure->mapCOMProcessChangeZ - zCor -
258  ( *trsZ * static_cast<proshade_double> ( staticStructure->getZDimSize() ) / staticStructure->getZDim() );
259 
260  //================================================ Save translation vector back
261  *trsX = -xMov;
262  *trsY = -yMov;
263  *trsZ = -zMov;
264 
265  //================================================ Report progress
266  std::stringstream hlpSS;
267  hlpSS << "Optimal map translation distances are " << *trsX << " ; " << *trsY << " ; " << *trsZ << " Angstroms with peak height " << mapPeak / ( xDimS * yDimS * zDimS );
268 
269  //================================================ Save original from variables for PDB output
270  movingStructure->mapMovFromsChangeX = movingStructure->xFrom;
271  movingStructure->mapMovFromsChangeY = movingStructure->yFrom;
272  movingStructure->mapMovFromsChangeZ = movingStructure->zFrom;
273 
274  //================================================ Move the map
275  ProSHADE_internal_mapManip::moveMapByIndices ( &xMov, &yMov, &zMov, movingStructure->getXDimSize(), movingStructure->getYDimSize(), movingStructure->getZDimSize(),
276  movingStructure->getXFromPtr(), movingStructure->getXToPtr(),
277  movingStructure->getYFromPtr(), movingStructure->getYToPtr(),
278  movingStructure->getZFromPtr(), movingStructure->getZToPtr(),
279  movingStructure->getXAxisOrigin(), movingStructure->getYAxisOrigin(), movingStructure->getZAxisOrigin() );
280 
281  ProSHADE_internal_mapManip::moveMapByFourier ( movingStructure->getInternalMap(), xMov, yMov, zMov,
282  movingStructure->getXDimSize(), movingStructure->getYDimSize(), movingStructure->getZDimSize(),
283  movingStructure->getXDim(), movingStructure->getYDim(), movingStructure->getZDim() );
284 
285  //================================================ Report progress
286  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, hlpSS.str() );
287  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Translation function computation complete." );
288 
289  //================================================ Keep only the change in from and to variables
290  movingStructure->mapMovFromsChangeX = movingStructure->xFrom - movingStructure->mapMovFromsChangeX;
291  movingStructure->mapMovFromsChangeY = movingStructure->yFrom - movingStructure->mapMovFromsChangeY;
292  movingStructure->mapMovFromsChangeZ = movingStructure->zFrom - movingStructure->mapMovFromsChangeZ;
293 
294  //================================================ Compute the optimal rotation centre for co-ordinates
295  movingStructure->computePdbRotationCentre ( );
296  movingStructure->computeOptimalTranslation ( eulA, eulB, eulG, *trsX, *trsY, *trsZ );
297 
298  //================================================ Done
299  return ;
300 
301 }
302 
308 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 )
309 {
310  //================================================ Initialise local variables
311  std::vector< proshade_double > ret;
312  proshade_double mapPeak = 0.0;
313  proshade_double trsX = 0.0, trsY = 0.0, trsZ = 0.0;
314  ProSHADE_internal_overlay::findHighestValueInMap ( this->getTranslationFnPointer(),
315  staticStructure->getXDim(),
316  staticStructure->getYDim(),
317  staticStructure->getZDim(),
318  &trsX,
319  &trsY,
320  &trsZ,
321  &mapPeak );
322 
323  //================================================ Dont translate over half
324  if ( trsX > ( staticStructure->getXDim() / 2 ) ) { trsX = trsX - this->getXDim(); }
325  if ( trsY > ( staticStructure->getYDim() / 2 ) ) { trsY = trsY - this->getYDim(); }
326  if ( trsZ > ( staticStructure->getZDim() / 2 ) ) { trsZ = trsZ - this->getZDim(); }
327 
328  //================================================ Move map
329  proshade_single xCor = ( staticStructure->xFrom - this->xFrom ) *
330  ( static_cast<proshade_double> ( staticStructure->getXDimSize() ) / staticStructure->getXDim() );
331  proshade_single yCor = ( staticStructure->yFrom - this->yFrom ) *
332  ( static_cast<proshade_double> ( staticStructure->getYDimSize() ) / staticStructure->getYDim() );
333  proshade_single zCor = ( staticStructure->zFrom - this->zFrom ) *
334  ( static_cast<proshade_double> ( staticStructure->getZDimSize() ) / staticStructure->getZDim() );
335  proshade_single xMov = staticStructure->mapCOMProcessChangeX - xCor -
336  ( trsX * static_cast<proshade_double> ( staticStructure->getXDimSize() ) / staticStructure->getXDim() );
337  proshade_single yMov = staticStructure->mapCOMProcessChangeY - yCor -
338  ( trsY * static_cast<proshade_double> ( staticStructure->getYDimSize() ) / staticStructure->getYDim() );
339  proshade_single zMov = staticStructure->mapCOMProcessChangeZ - zCor -
340  ( trsZ * static_cast<proshade_double> ( staticStructure->getZDimSize() ) / staticStructure->getZDim() );
341  proshade_single modXMov = xMov;
342  proshade_single modYMov = yMov;
343  proshade_single modZMov = zMov;
344 
345  //================================================ Save results as vector
346  ProSHADE_internal_misc::addToDoubleVector ( &ret, static_cast<proshade_double> ( -xMov ) );
347  ProSHADE_internal_misc::addToDoubleVector ( &ret, static_cast<proshade_double> ( -yMov ) );
348  ProSHADE_internal_misc::addToDoubleVector ( &ret, static_cast<proshade_double> ( -zMov ) );
349 
350  //================================================ Save original from variables for PDB output
351  this->mapMovFromsChangeX = this->xFrom;
352  this->mapMovFromsChangeY = this->yFrom;
353  this->mapMovFromsChangeZ = this->zFrom;
354 
355  //================================================ Move the map
356  ProSHADE_internal_mapManip::moveMapByIndices ( &modXMov, &modYMov, &modZMov, this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
357  this->getXFromPtr(), this->getXToPtr(),
358  this->getYFromPtr(), this->getYToPtr(),
359  this->getZFromPtr(), this->getZToPtr(),
360  this->getXAxisOrigin(), this->getYAxisOrigin(), this->getZAxisOrigin() );
361 
362  ProSHADE_internal_mapManip::moveMapByFourier ( this->getInternalMap(), modXMov, modYMov, modZMov,
363  this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
364  this->getXDim(), this->getYDim(), this->getZDim() );
365 
366  //================================================ Keep only the change in from and to variables
367  this->mapMovFromsChangeX = this->xFrom - this->mapMovFromsChangeX;
368  this->mapMovFromsChangeY = this->yFrom - this->mapMovFromsChangeY;
369  this->mapMovFromsChangeZ = this->zFrom - this->mapMovFromsChangeZ;
370 
371  //================================================ Compute the optimal rotation centre for co-ordinates
372  this->computePdbRotationCentre ( );
373  this->computeOptimalTranslation ( eulA, eulB, eulG, -xMov, -yMov, -zMov );
374 
375  //================================================ Done
376  return ( ret );
377 
378 }
379 
390 {
391  //================================================ Do this using Fourier!
392  fftw_complex *tmpIn1 = NULL, *tmpOut1 = NULL, *tmpIn2 = NULL, *tmpOut2 = NULL, *resOut = NULL;
393  fftw_plan forwardFourierObj1, forwardFourierObj2, inverseFourierCombo;
394  proshade_unsign dimMult = staticStructure->getXDim() * staticStructure->getYDim() * staticStructure->getZDim();
395  ProSHADE_internal_overlay::allocateTranslationFunctionMemory ( tmpIn1, tmpOut1, tmpIn2, tmpOut2, this->translationMap, resOut, forwardFourierObj1, forwardFourierObj2, inverseFourierCombo, staticStructure->getXDim(), staticStructure->getYDim(), staticStructure->getZDim() );
396 
397  //================================================ Fill in input data
398  for ( proshade_unsign iter = 0; iter < dimMult; iter++ ) { tmpIn1[iter][0] = staticStructure->getMapValue ( iter ); tmpIn1[iter][1] = 0.0; }
399  for ( proshade_unsign iter = 0; iter < dimMult; iter++ ) { tmpIn2[iter][0] = this->getMapValue ( iter ); tmpIn2[iter][1] = 0.0; }
400 
401  //================================================ Calculate Fourier
402  fftw_execute ( forwardFourierObj1 );
403  fftw_execute ( forwardFourierObj2 );
404 
405  //================================================ Combine Fourier coeffs and invert
406  ProSHADE_internal_overlay::combineFourierForTranslation ( tmpOut1, tmpOut2, resOut, staticStructure->getXDim(), staticStructure->getYDim(), staticStructure->getZDim() );
407  fftw_execute ( inverseFourierCombo );
408 
409  //================================================ Free memory
410  ProSHADE_internal_overlay::freeTranslationFunctionMemory ( tmpIn1, tmpOut1, tmpIn2, tmpOut2, resOut, forwardFourierObj1, forwardFourierObj2, inverseFourierCombo );
411 
412  //================================================ Done
413  return ;
414 
415 }
416 
432 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 )
433 {
434  //================================================ Allocate memory
435  tmpIn1 = new fftw_complex[xD * yD * zD];
436  tmpOut1 = new fftw_complex[xD * yD * zD];
437  tmpIn2 = new fftw_complex[xD * yD * zD];
438  tmpOut2 = new fftw_complex[xD * yD * zD];
439  resIn = new fftw_complex[xD * yD * zD];
440  resOut = new fftw_complex[xD * yD * zD];
441 
442  //================================================ Check memory allocation
443  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn1 , __FILE__, __LINE__, __func__ );
444  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut1, __FILE__, __LINE__, __func__ );
445  ProSHADE_internal_misc::checkMemoryAllocation ( tmpIn2, __FILE__, __LINE__, __func__ );
446  ProSHADE_internal_misc::checkMemoryAllocation ( tmpOut2, __FILE__, __LINE__, __func__ );
447  ProSHADE_internal_misc::checkMemoryAllocation ( resIn, __FILE__, __LINE__, __func__ );
448  ProSHADE_internal_misc::checkMemoryAllocation ( resOut, __FILE__, __LINE__, __func__ );
449 
450  //================================================ Get Fourier transforms of the maps
451  forwardFourierObj1 = fftw_plan_dft_3d ( xD, yD, zD, tmpIn1, tmpOut1, FFTW_FORWARD , FFTW_ESTIMATE );
452  forwardFourierObj2 = fftw_plan_dft_3d ( xD, yD, zD, tmpIn2, tmpOut2, FFTW_FORWARD , FFTW_ESTIMATE );
453  inverseFourierCombo = fftw_plan_dft_3d ( xD, yD, zD, resOut, resIn , FFTW_BACKWARD, FFTW_ESTIMATE );
454 
455  //================================================ Done
456  return ;
457 
458 }
459 
471 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 )
472 {
473  //================================================ Release memory
474  fftw_destroy_plan ( forwardFourierObj1 );
475  fftw_destroy_plan ( forwardFourierObj2 );
476  fftw_destroy_plan ( inverseFourierCombo );
477  delete[] tmpIn1;
478  delete[] tmpIn2;
479  delete[] tmpOut1;
480  delete[] tmpOut2;
481  delete[] resOut;
482 
483  //======================================== Done
484  return ;
485 
486 }
487 
497 void ProSHADE_internal_overlay::combineFourierForTranslation ( fftw_complex* tmpOut1, fftw_complex* tmpOut2, fftw_complex*& resOut, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD )
498 {
499  //================================================ Initialise local variables
500  double normFactor = static_cast<double> ( xD * yD * zD );
501  proshade_signed h1, k1, l1, hlpPos, arrPos;
502  proshade_signed uo2 = xD / 2;
503  proshade_signed vo2 = yD / 2;
504  proshade_signed wo2 = zD / 2;
505 
506  //================================================ Combine the coefficients
507  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> ( xD ); uIt++ )
508  {
509  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> ( yD ); vIt++ )
510  {
511  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> ( zD ); wIt++ )
512  {
513  //==================================== Convert to HKL
514  if ( uIt > uo2 ) { h1 = uIt - xD; } else { h1 = uIt; }
515  if ( vIt > vo2 ) { k1 = vIt - yD; } else { k1 = vIt; }
516  if ( wIt > wo2 ) { l1 = wIt - zD; } else { l1 = wIt; }
517 
518  //==================================== Make HKL into indexable numbers
519  if ( h1 < 0 ) { h1 += xD; }
520  if ( k1 < 0 ) { k1 += yD; }
521  if ( l1 < 0 ) { l1 += zD; }
522 
523  //==================================== Find indices
524  hlpPos = l1 + zD * ( k1 + yD * h1 );
525  arrPos = wIt + zD * ( vIt + yD * uIt );
526 
527  //==================================== Combine
529  &tmpOut1[hlpPos][1],
530  &tmpOut2[arrPos][0],
531  &tmpOut2[arrPos][1],
532  &resOut[hlpPos][0],
533  &resOut[hlpPos][1] );
534 
535  //==================================== Save
536  resOut[hlpPos][0] /= normFactor;
537  resOut[hlpPos][1] /= normFactor;
538  }
539  }
540  }
541 
542  //======================================== Done
543  return ;
544 
545 }
546 
558 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 )
559 {
560  //================================================ Initialise variables
561  proshade_signed arrPos;
562  *mapPeak = 0.0;
563 
564  //================================================ Search the map
565  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> ( xD ); uIt++ )
566  {
567  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> ( yD ); vIt++ )
568  {
569  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> ( zD ); wIt++ )
570  {
571  arrPos = wIt + zD * ( vIt + yD * uIt );
572  if ( resIn[arrPos][0] > *mapPeak )
573  {
574  *mapPeak = resIn[arrPos][0];
575  *trsX = uIt;
576  *trsY = vIt;
577  *trsZ = wIt;
578  }
579  }
580  }
581  }
582 
583  //================================================ Done
584  return ;
585 
586 }
587 
598 void ProSHADE_internal_data::ProSHADE_data::zeroPaddToDims ( proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim )
599 {
600  //================================================ Sanity check
601  if ( ( this->xDimIndices > xDim ) || ( this->yDimIndices > yDim ) || ( this->zDimIndices > zDim ) )
602  {
603  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." );
604  }
605 
606  //================================================ If done, do nothing
607  if ( ( this->xDimIndices == xDim ) && ( this->yDimIndices == yDim ) && ( this->zDimIndices == zDim ) ) { return ; }
608 
609  //================================================ Find out how many zeroes need to be added before and after
610  proshade_unsign addXPre, addYPre, addZPre, addXPost, addYPost, addZPost;
611  ProSHADE_internal_overlay::computeBeforeAfterZeroCounts ( &addXPre, &addYPre, &addZPre, &addXPost, &addYPost, &addZPost, xDim, yDim, zDim, this->xDimIndices, this->yDimIndices, this->zDimIndices );
612 
613  //================================================ Create a new map
614  proshade_double* newMap = new proshade_double [xDim * yDim * zDim];
615 
616  //================================================ Do the hard work
617  ProSHADE_internal_overlay::paddMapWithZeroes ( this->internalMap, newMap, xDim, yDim, zDim, this->xDimIndices, this->yDimIndices, this->zDimIndices, addXPre, addYPre, addZPre );
618 
619  //================================================ Create a new internal map and copy
620  delete[] this->internalMap;
621  this->internalMap = new proshade_double [xDim * yDim * zDim];
622  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( xDim * yDim * zDim ); iter++ ) { this->internalMap[iter] = newMap[iter]; }
623 
624  //================================================ Release memory
625  delete[] newMap;
626 
627  //================================================ Update map related variables
628  this->xDimSize = xDim * ( this->xDimSize / this->xDimIndices );
629  this->yDimSize = yDim * ( this->yDimSize / this->yDimIndices );
630  this->zDimSize = zDim * ( this->zDimSize / this->zDimIndices );
631  this->xDimIndices = xDim ; this->yDimIndices = yDim ; this->zDimIndices = zDim;
632  this->xGridIndices = xDim ; this->yGridIndices = yDim ; this->zGridIndices = zDim;
633  this->xFrom -= addXPre ; this->yFrom -= addYPre ; this->zFrom -= addZPre;
634  this->xTo += addXPost; this->yTo += addYPost; this->zTo += addZPost;
635  this->xAxisOrigin -= addXPre ; this->yAxisOrigin -= addYPre ; this->zAxisOrigin -= addZPre ;
636 
637  //================================================ Done
638  return ;
639 
640 }
641 
657 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 )
658 {
659  //================================================ Compute
660  *addXPre = ( xDim - xDimIndices ) / 2;
661  *addYPre = ( yDim - yDimIndices ) / 2;
662  *addZPre = ( zDim - zDimIndices ) / 2;
663  *addXPost = *addXPre; if ( ( xDim - xDimIndices ) % 2 == 1 ) { *addXPost += 1; }
664  *addYPost = *addYPre; if ( ( yDim - yDimIndices ) % 2 == 1 ) { *addYPost += 1; }
665  *addZPost = *addZPre; if ( ( zDim - zDimIndices ) % 2 == 1 ) { *addZPost += 1; }
666 
667  //================================================ Done
668  return ;
669 
670 }
671 
686 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 )
687 {
688  //================================================ Initialise local variables
689  proshade_unsign newMapIndex = 0;
690  proshade_unsign oldMapIndex = 0;
691 
692  //================================================ Copy and padd as appropriate
693  for ( proshade_unsign xIt = 0; xIt < xDim; xIt++ )
694  {
695  for ( proshade_unsign yIt = 0; yIt < yDim; yIt++ )
696  {
697  for ( proshade_unsign zIt = 0; zIt < zDim; zIt++ )
698  {
699  //==================================== Find map location
700  newMapIndex = zIt + zDim * ( yIt + yDim * xIt );
701 
702  //==================================== If less than needed, add zero
703  if ( xIt < addXPre ) { newMap[newMapIndex] = 0.0; continue; }
704  if ( yIt < addYPre ) { newMap[newMapIndex] = 0.0; continue; }
705  if ( zIt < addZPre ) { newMap[newMapIndex] = 0.0; continue; }
706 
707  //==================================== If more than needed, add zero
708  if ( xIt >= ( addXPre + xDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
709  if ( yIt >= ( addYPre + yDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
710  if ( zIt >= ( addZPre + zDimIndices ) ) { newMap[newMapIndex] = 0.0; continue; }
711 
712  //==================================== If not padding, copy original values
713  oldMapIndex = (zIt-addZPre) + zDimIndices * ( (yIt-addYPre) + yDimIndices * (xIt-addXPre) );
714  newMap[newMapIndex] = oldMap[oldMapIndex];
715  }
716  }
717  }
718 
719  //======================================== Done
720  return ;
721 
722 }
723 
736 void ProSHADE_internal_data::ProSHADE_data::rotateMap ( ProSHADE_settings* settings, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma )
737 {
738  //================================================ Set maximum comparison bandwidth to maximum object bandwidth
739  this->maxCompBand = this->spheres[this->noSpheres-1]->getLocalBandwidth();
740 
741  //================================================ Save map COM after processing but before rotation
742  this->findMapCOM ( );
743  this->mapCOMProcessChangeX = this->xCom - this->originalMapXCom;
744  this->mapCOMProcessChangeY = this->yCom - this->originalMapYCom;
745  this->mapCOMProcessChangeZ = this->zCom - this->originalMapZCom;
746 
747  //================================================ Compute the Wigner D matrices for the Euler angles
748  ProSHADE_internal_wigner::computeWignerMatricesForRotation ( settings, this, -eulerAlpha, eulerBeta, -eulerGamma );
749 
750  //================================================ Initialise rotated Spherical Harmonics memory
751  this->allocateRotatedSHMemory ( settings );
752 
753  //================================================ Multiply SH coeffs by Wigner
754  this->computeRotatedSH ( settings );
755 
756  //================================================ Inverse the SH coeffs to shells
757  this->invertSHCoefficients ( );
758 
759  //================================================ Find spherical cut-offs
760  std::vector<proshade_double> lonCO, latCO;
761  ProSHADE_internal_overlay::computeAngularThreshold ( &lonCO, &latCO, settings->maxBandwidth * 2 );
762 
763  //================================================ Allocate memory for the rotated map
764  proshade_double *densityMapRotated = new proshade_double [this->xDimIndices * this->yDimIndices * this->zDimIndices];
765  ProSHADE_internal_misc::checkMemoryAllocation ( densityMapRotated, __FILE__, __LINE__, __func__ );
766  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( this->xDimIndices * this->yDimIndices * this->zDimIndices ); iter++ ) { densityMapRotated[iter] = 0.0; }
767 
768  //================================================ Interpolate onto cartesian grid
769  this->interpolateMapFromSpheres ( settings, densityMapRotated );
770 
771  //================================================ Copy map
772  for ( proshade_unsign iter = 0; iter < ( this->xDimIndices * this->yDimIndices * this->zDimIndices ); iter++ )
773  {
774  this->internalMap[iter] = densityMapRotated[iter];
775  }
776 
777  //================================================ Release rotated map (original is now rotated)
778  delete[] densityMapRotated;
779 
780  //================================================ Done
781  return ;
782 
783 }
784 
795 void ProSHADE_internal_data::ProSHADE_data::translateMap ( ProSHADE_settings* settings, proshade_double trsX, proshade_double trsY, proshade_double trsZ )
796 {
797  //================================================ Initialise local variables
798  proshade_single xMov = -trsX, yMov = -trsY, zMov = -trsZ;
799 
800  //================================================ Move the whole map frame to minimise the Fourier movement
801  ProSHADE_internal_mapManip::moveMapByIndices ( &xMov, &yMov, &zMov, this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
802  this->getXFromPtr(), this->getXToPtr(), this->getYFromPtr(), this->getYToPtr(),
803  this->getZFromPtr(), this->getZToPtr(), this->getXAxisOrigin(), this->getYAxisOrigin(), this->getZAxisOrigin() );
804 
805  //================================================ Finalise the movement by in-frame Fourier movement
806  ProSHADE_internal_mapManip::moveMapByFourier ( this->getInternalMap(), xMov, yMov, zMov, this->getXDimSize(), this->getYDimSize(), this->getZDimSize(),
807  this->getXDim(), this->getYDim(), this->getZDim() );
808 
809  //================================================ Done
810  return ;
811 
812 }
813 
819 {
820  //================================================ Allocate the main pointer and check
821  this->rotSphericalHarmonics = new proshade_complex* [this->noSpheres];
822  ProSHADE_internal_misc::checkMemoryAllocation ( this->rotSphericalHarmonics, __FILE__, __LINE__, __func__ );
823 
824  //================================================ For each sphere
825  for ( proshade_unsign iter = 0; iter < this->noSpheres; iter++ )
826  {
827  //============================================ Allocate the sphere storage space
828  this->rotSphericalHarmonics[iter] = new proshade_complex [static_cast<proshade_unsign> ( pow ( (this->spheres[iter]->getLocalBandwidth() * 2), 2) )];
829  ProSHADE_internal_misc::checkMemoryAllocation ( this->rotSphericalHarmonics[iter], __FILE__, __LINE__, __func__ );
830 
831  //============================================ Set values to zeroes
832  for ( proshade_unsign it = 0; it < static_cast<proshade_unsign> ( pow ( (this->spheres[iter]->getLocalBandwidth() * 2), 2) ); it++ )
833  {
834  this->rotSphericalHarmonics[iter][it][0] = 0.0;
835  this->rotSphericalHarmonics[iter][it][1] = 0.0;
836  }
837  }
838 
839  //================================================ Done
840  return ;
841 
842 }
843 
849 {
850  //================================================ Initialise variables
851  proshade_double WigDR, WigDI, *ShR, *ShI, retR, retI;
852  proshade_unsign arrPos;
853 
854  //================================================ Compute
855  for ( proshade_signed shell = 0; shell < static_cast<proshade_signed> ( this->noSpheres ); shell++ )
856  {
857  //============================================ For each band
858  for ( proshade_signed bandIter = 0; bandIter < static_cast<proshade_signed> ( this->spheres[shell]->getLocalBandwidth() ); bandIter++ )
859  {
860  //======================================== For each order1
861  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
862  {
863  //==================================== Get Spherical Harmonics value
864  ShR = getRealSphHarmValue ( bandIter, order1, shell );
865  ShI = getImagSphHarmValue ( bandIter, order1, shell );
866 
867  //==================================== For each order2
868  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
869  {
870  //================================ Get Wigner D values
871  this->getWignerMatrixValue ( bandIter, order1, order2, &WigDR, &WigDI );
872 
873  //================================ Multiply SH and Wigner
874  ProSHADE_internal_maths::complexMultiplication ( ShR, ShI, &WigDR, &WigDI, &retR, &retI );
875 
876  //================================ Save
877  arrPos = static_cast<proshade_unsign> ( seanindex ( order2-bandIter, bandIter, this->spheres[shell]->getLocalBandwidth() ) );
878  this->rotSphericalHarmonics[shell][arrPos][0] += retR;
879  this->rotSphericalHarmonics[shell][arrPos][1] += retI;
880  }
881  }
882  }
883  }
884 
885  //================================================ Done
886  return ;
887 
888 }
889 
902 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 )
903 {
904  //================================================ Initialise internal variables
905  proshade_unsign oneDim = shBand * 2;
906 
907  //================================================ Allocate memory
908  sigR = new double [(oneDim * oneDim * 4)];
909  sigI = sigR + (oneDim * oneDim * 2);
910  rcoeffs = new double [(oneDim * oneDim * 2)];
911  icoeffs = rcoeffs + (oneDim * oneDim);
912  weights = new double [4 * shBand];
913  workspace = new double [( 20 * shBand * shBand ) + ( 24 * shBand )];
914 
915  //================================================ Check memory allocation
916  ProSHADE_internal_misc::checkMemoryAllocation ( sigR, __FILE__, __LINE__, __func__ );
917  ProSHADE_internal_misc::checkMemoryAllocation ( rcoeffs, __FILE__, __LINE__, __func__ );
918  ProSHADE_internal_misc::checkMemoryAllocation ( weights, __FILE__, __LINE__, __func__ );
919  ProSHADE_internal_misc::checkMemoryAllocation ( workspace, __FILE__, __LINE__, __func__ );
920 
921  //================================================ Create the cosine/sine transform plan
922  idctPlan = fftw_plan_r2r_1d ( oneDim, weights, workspace, FFTW_REDFT01, FFTW_ESTIMATE );
923 
924  //================================================ Set up the discrete Fourier transform
925  int rank, howmany_rank;
926  fftw_iodim dims[1], howmany_dims[1];
927 
928  rank = 1;
929  dims[0].n = 2 * shBand;
930  dims[0].is = 2 * shBand;
931  dims[0].os = 1;
932  howmany_rank = 1;
933  howmany_dims[0].n = 2 * shBand;
934  howmany_dims[0].is = 1;
935  howmany_dims[0].os = 2 * shBand;
936 
937  //================================================ Create the discrete Fourier transform
938  ifftPlan = fftw_plan_guru_split_dft ( rank, dims, howmany_rank, howmany_dims, sigR, sigI, rcoeffs, icoeffs, FFTW_ESTIMATE );
939 
940  //================================================ Done
941  return ;
942 
943 }
944 
949 {
950  //================================================ Initialise local variables
951  double *sigR = NULL, *sigI = NULL, *rcoeffs = NULL, *icoeffs = NULL, *weights = NULL, *workspace = NULL;
952  fftw_plan idctPlan, ifftPlan;
953 
954  //================================================ For each shell
955  for ( int shell = 0; shell < static_cast<int> ( this->noSpheres ); shell++ )
956  {
957  //=========================================== Initialise internal variables
958  proshade_unsign oneDim = this->spheres[shell]->getLocalBandwidth() * 2;
959 
960  //=========================================== Allocate memory
961  ProSHADE_internal_overlay::initialiseInverseSHComputation ( this->spheres[shell]->getLocalBandwidth(), sigR, sigI, rcoeffs, icoeffs, weights, workspace, idctPlan, ifftPlan );
962 
963  //=========================================== Compute weights for the transform using the appropriate shell related band
964  makeweights ( this->spheres[shell]->getLocalBandwidth(), weights );
965 
966  //============================================ Allocate rotated shell mapped data memory
967  this->spheres[shell]->allocateRotatedMap ( );
968 
969  //============================================ Load SH coeffs to arrays
970  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( oneDim * oneDim ); iter++ )
971  {
972  rcoeffs[iter] = this->rotSphericalHarmonics[shell][iter][0];
973  icoeffs[iter] = this->rotSphericalHarmonics[shell][iter][1];
974  sigR[iter] = 0.0;
975  sigI[iter] = 0.0;
976  }
977 
978  //============================================ Get inverse spherical harmonics transform for the shell
979  InvFST_semi_fly ( rcoeffs,
980  icoeffs,
981  sigR,
982  sigI,
983  this->spheres[shell]->getLocalBandwidth(),
984  workspace,
985  0,
986  this->spheres[shell]->getLocalBandwidth(),
987  &idctPlan,
988  &ifftPlan );
989 
990  //=========================================== Copy the results to the rotated shells array
991  for ( unsigned int iter = 0; iter < static_cast<unsigned int> ( oneDim * oneDim ); iter++ )
992  {
993  this->spheres[shell]->setRotatedMappedData ( iter, sigR[iter] );
994  }
995 
996  //=========================================== Release the plans
997  fftw_destroy_plan ( idctPlan );
998  fftw_destroy_plan ( ifftPlan );
999 
1000  //=========================================== Release the memory
1001  delete[] sigR;
1002  delete[] rcoeffs;
1003  delete[] weights;
1004  delete[] workspace;
1005  }
1006 
1007  //================================================ Done
1008  return ;
1009 
1010 }
1011 
1018 void ProSHADE_internal_overlay::computeAngularThreshold ( std::vector<proshade_double>* lonCO, std::vector<proshade_double>* latCO, proshade_unsign angRes )
1019 {
1020  //================================================ Longitude angular thresholds
1021  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1022  {
1023  ProSHADE_internal_misc::addToDoubleVector ( lonCO, static_cast<proshade_double> ( iter ) * ( ( static_cast<proshade_double> ( M_PI ) * 2.0 ) / static_cast<proshade_double> ( angRes ) ) - ( static_cast<double> ( M_PI ) ) );
1024  }
1025 
1026  //================================================ Lattitude angular thresholds
1027  for ( proshade_unsign iter = 0; iter <= angRes; iter++ )
1028  {
1029  ProSHADE_internal_misc::addToDoubleVector ( latCO, static_cast<proshade_double> ( iter ) * ( static_cast<proshade_double> ( M_PI ) / static_cast<proshade_double> ( angRes ) ) - ( static_cast<proshade_double> ( M_PI ) / 2.0 ) );
1030  }
1031 
1032  //================================================ Done
1033  return ;
1034 
1035 }
1036 
1042 void ProSHADE_internal_data::ProSHADE_data::interpolateMapFromSpheres ( ProSHADE_settings* settings, proshade_double*& densityMapRotated )
1043 {
1044  //================================================ Initialise variables
1045  proshade_double rad = 0.0, lon = 0.0, lat = 0.0, newU = 0.0, newV = 0.0, newW = 0.0;
1046  proshade_unsign lowerLonL = 0, upperLonL = 0, lowerLonU = 0, upperLonU = 0, lowerLatL = 0, upperLatL = 0, lowerLatU = 0, upperLatU = 0, lowerShell = 0, upperShell = 0;
1047  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;
1048  proshade_double lowerShellValue = 0.0, upperShellValue = 0.0;
1049  proshade_double xSamplingRate = static_cast<proshade_double> ( this->xDimSize ) / this->xDimIndices;
1050  proshade_double ySamplingRate = static_cast<proshade_double> ( this->yDimSize ) / this->yDimIndices;
1051  proshade_double zSamplingRate = static_cast<proshade_double> ( this->zDimSize ) / this->zDimIndices;
1052  proshade_signed arrPos;
1053  std::vector<proshade_double> lonCOU, latCOU, lonCOL, latCOL;
1054 
1055  for ( proshade_signed uIt = 0; uIt < static_cast<proshade_signed> (this->xDimIndices); uIt++ )
1056  {
1057  for ( proshade_signed vIt = 0; vIt < static_cast<proshade_signed> (this->yDimIndices); vIt++ )
1058  {
1059  for ( proshade_signed wIt = 0; wIt < static_cast<proshade_signed> (this->zDimIndices); wIt++ )
1060  {
1061  //==================================== Convert to centered coords
1062  newU = static_cast<proshade_double> ( uIt - ( static_cast<proshade_signed> (this->xDimIndices) / 2 ) );
1063  newV = static_cast<proshade_double> ( vIt - ( static_cast<proshade_signed> (this->yDimIndices) / 2 ) );
1064  newW = static_cast<proshade_double> ( wIt - ( static_cast<proshade_signed> (this->zDimIndices) / 2 ) );
1065 
1066  //==================================== Deal with 0 ; 0 ; 0
1067  if ( ( newU == 0.0 ) && ( newV == 0.0 ) && ( newW == 0.0 ) )
1068  {
1069  arrPos = wIt + this->zDimIndices * ( vIt + this->yDimIndices * uIt );
1070  densityMapRotated[arrPos] = this->internalMap[arrPos];
1071  continue;
1072  }
1073 
1074  //==================================== Convert to spherical coords
1075  rad = sqrt ( pow( ( newU * xSamplingRate ), 2.0 ) +
1076  pow( ( newV * ySamplingRate ), 2.0 ) +
1077  pow( ( newW * zSamplingRate ), 2.0 ) );
1078  lon = atan2 ( ( newV * ySamplingRate ), ( newU * xSamplingRate ) );
1079  lat = asin ( ( newW * zSamplingRate ) / rad );
1080 
1081  //==================================== Deal with nan's
1082  if ( rad != rad ) { rad = 0.0; }
1083  if ( lon != lon ) { lon = 0.0; }
1084  if ( lat != lat ) { lat = 0.0; }
1085 
1086  //==================================== Find shells above and below
1087  lowerShell = 0;
1088  upperShell = 0;
1089  for ( proshade_unsign iter = 0; iter < (this->noSpheres-1); iter++ )
1090  {
1091  if ( ( this->spherePos.at(iter) <= rad ) && ( this->spherePos.at(iter+1) > rad ) )
1092  {
1093  lowerShell = iter;
1094  upperShell = iter+1;
1095  break;
1096  }
1097  }
1098 
1099  if ( upperShell == 0 )
1100  {
1101  arrPos = wIt + this->zDimIndices * ( vIt + this->yDimIndices * uIt );
1102  densityMapRotated[arrPos] = 0.0;
1103  continue;
1104  }
1105 
1106  //==================================== Get the longitude and lattitude cut-offs for this shell resolution
1107  lonCOL.clear(); latCOL.clear(); lonCOU.clear(); latCOU.clear();
1108  ProSHADE_internal_overlay::computeAngularThreshold ( &lonCOL, &latCOL, this->spheres[lowerShell]->getLocalAngRes() );
1109  ProSHADE_internal_overlay::computeAngularThreshold ( &lonCOU, &latCOU, this->spheres[upperShell]->getLocalAngRes() );
1110 
1111  //==================================== Find the angle cutoffs around the point
1112  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( lonCOL.size() ); iter++ )
1113  {
1114  if ( iter == ( static_cast<proshade_unsign> ( lonCOL.size() ) - 1 ) )
1115  {
1116  lowerLonL = 0;
1117  upperLonL = 1;
1118  break;
1119  }
1120  if ( ( std::floor(10000. * lonCOL.at(iter)) <= std::floor(10000. * lon) ) && ( std::floor(10000. * lonCOL.at(iter+1)) > std::floor(10000. * lon) ) )
1121  {
1122  lowerLonL = iter;
1123  upperLonL = iter+1;
1124  break;
1125  }
1126  }
1127  if ( upperLonL == this->spheres[lowerShell]->getLocalAngRes() ) { upperLonL = 0; }
1128 
1129  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( lonCOU.size() ); iter++ )
1130  {
1131  if ( iter == ( static_cast<proshade_unsign> ( lonCOU.size() ) - 1 ) )
1132  {
1133  lowerLonU = 0;
1134  upperLonU = 1;
1135  break;
1136  }
1137  if ( ( std::floor(10000. * lonCOU.at(iter)) <= std::floor(10000. * lon) ) && ( std::floor(10000. * lonCOU.at(iter+1)) > std::floor(10000. * lon) ) )
1138  {
1139  lowerLonU = iter;
1140  upperLonU = iter+1;
1141  break;
1142  }
1143  }
1144  if ( upperLonU == this->spheres[upperShell]->getLocalAngRes() ) { upperLonU = 0; }
1145 
1146  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( latCOL.size() ); iter++ )
1147  {
1148  if ( iter == ( static_cast<proshade_unsign> ( latCOL.size() ) - 1 ) )
1149  {
1150  lowerLatL = 0;
1151  upperLatL = 1;
1152  break;
1153  }
1154  if ( ( std::floor(10000. * latCOL.at(iter)) <= std::floor(10000. * lat) ) && ( std::floor(10000. * latCOL.at(iter+1)) > std::floor(10000. * lat) ) )
1155  {
1156  lowerLatL = iter;
1157  upperLatL = iter+1;
1158  break;
1159  }
1160  }
1161  if ( upperLatL == this->spheres[lowerShell]->getLocalAngRes() ) { upperLatL = 0; }
1162 
1163  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( latCOU.size() ); iter++ )
1164  {
1165  if ( iter == ( static_cast<proshade_unsign> ( latCOU.size() ) - 1 ) )
1166  {
1167  lowerLatU = 0;
1168  upperLatU = 1;
1169  break;
1170  }
1171  if ( ( std::floor(10000. * latCOU.at(iter)) <= std::floor(10000. * lat) ) && ( std::floor(10000. * latCOU.at(iter+1)) > std::floor(10000. * lat) ) )
1172  {
1173  lowerLatU = iter;
1174  upperLatU = iter+1;
1175  break;
1176  }
1177  }
1178  if ( upperLatU == this->spheres[upperShell]->getLocalAngRes() ) { upperLatU = 0; }
1179 
1180  //==================================== Interpolate lower shell
1181  x00 = this->spheres[lowerShell]->getRotatedMappedData ( lowerLatL * this->spheres[lowerShell]->getLocalAngRes() + lowerLonL );
1182  x01 = this->spheres[lowerShell]->getRotatedMappedData ( lowerLatL * this->spheres[lowerShell]->getLocalAngRes() + upperLonL );
1183  x10 = this->spheres[lowerShell]->getRotatedMappedData ( upperLatL * this->spheres[lowerShell]->getLocalAngRes() + lowerLonL );
1184  x11 = this->spheres[lowerShell]->getRotatedMappedData ( upperLatL * this->spheres[lowerShell]->getLocalAngRes() + upperLonL );
1185 
1186  distLLon = std::abs ( lon - lonCOL.at(lowerLonL) ) / ( std::abs( lon - lonCOL.at(lowerLonL) ) + std::abs( lon - lonCOL.at(upperLonL) ) );
1187  valLLon = ( ( 1.0 - distLLon ) * x00 ) + ( distLLon * x01 );
1188  valULon = ( ( 1.0 - distLLon ) * x10 ) + ( distLLon * x11 );
1189 
1190  distLLat = std::abs ( lat - latCOL.at(lowerLatL) ) / ( std::abs( lat - latCOL.at(lowerLatL) ) + std::abs( lat - latCOL.at(upperLatL) ) );
1191  lowerShellValue = ( ( 1.0 - distLLat ) * valLLon ) + ( distLLat * valULon );
1192 
1193  //==================================== Interpolate upper shell
1194  x00 = this->spheres[upperShell]->getRotatedMappedData ( lowerLatU * this->spheres[upperShell]->getLocalAngRes() + lowerLonU );
1195  x01 = this->spheres[upperShell]->getRotatedMappedData ( lowerLatU * this->spheres[upperShell]->getLocalAngRes() + upperLonU );
1196  x10 = this->spheres[upperShell]->getRotatedMappedData ( upperLatU * this->spheres[upperShell]->getLocalAngRes() + lowerLonU );
1197  x11 = this->spheres[upperShell]->getRotatedMappedData ( upperLatU * this->spheres[upperShell]->getLocalAngRes() + upperLonU );
1198 
1199  distLLon = std::abs ( lon - lonCOU.at(lowerLonU) ) / ( std::abs( lon - lonCOU.at(lowerLonU) ) + std::abs( lon - lonCOU.at(upperLonU) ) );
1200  valLLon = ( ( 1.0 - distLLon ) * x00 ) + ( distLLon * x01 );
1201  valULon = ( ( 1.0 - distLLon ) * x10 ) + ( distLLon * x11 );
1202 
1203  distLLat = std::abs ( lat - latCOU.at(lowerLatU) ) / ( std::abs( lat - latCOU.at(lowerLatU) ) + std::abs( lat - latCOU.at(upperLatU) ) );
1204  upperShellValue = ( ( 1.0 - distLLat ) * valLLon ) + ( distLLat * valULon );
1205 
1206  //==================================== Interpolate between shells
1207  distLRad = std::abs ( rad - this->spherePos.at(lowerShell) ) / ( std::abs( rad - this->spherePos.at(lowerShell) ) +
1208  std::abs( rad - this->spherePos.at(upperShell) ) );
1209 
1210  arrPos = wIt + this->zDimIndices * ( vIt + this->yDimIndices * uIt );
1211  densityMapRotated[arrPos] = ( ( 1.0 - distLRad ) * lowerShellValue ) + ( distLRad * upperShellValue );
1212  }
1213 
1214  }
1215 
1216  }
1217 
1218  //================================================ Done
1219  return ;
1220 
1221 }
1222 
1229 {
1230  //================================================ Initialise local variables
1231  std::vector < proshade_double > ret;
1232  proshade_double eulA, eulB, eulG;
1233 
1234  //================================================ Get inverse SO(3) map top peak Euler angle values
1236  this->getMaxBand() * 2,
1237  &eulA, &eulB, &eulG, settings );
1238 
1239  //================================================ Re-format to vector
1243 
1244  //================================================ Done
1245  return ( ret );
1246 
1247 }
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:577
ProSHADE_settings::maxBandwidth
proshade_unsign maxBandwidth
The bandwidth of spherical harmonics decomposition for the largest sphere.
Definition: ProSHADE_settings.hpp:109
ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics
void computeSphericalHarmonics(ProSHADE_settings *settings)
This function computes the spherical harmonics decomposition for the whole structure.
Definition: ProSHADE_data.cpp:1608
ProSHADE_internal_data::ProSHADE_data::mapMovFromsChangeX
proshade_double mapMovFromsChangeX
When the map is translated, the xFrom and xTo values are changed. This variable holds how much they h...
Definition: ProSHADE_data.hpp:94
ProSHADE_internal_data::ProSHADE_data::getOverlayRotationFunction
void getOverlayRotationFunction(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *obj2)
This function computes the overlay rotation function (i.e. the correlation function in SO(3) space).
Definition: ProSHADE_overlay.cpp:35
ProSHADE_internal_data::ProSHADE_data::zFrom
proshade_signed zFrom
This is the starting index along the z axis.
Definition: ProSHADE_data.hpp:112
ProSHADE_internal_data::ProSHADE_data::getZDimSize
proshade_single getZDimSize(void)
This function allows access to the map size in angstroms along the Z axis.
Definition: ProSHADE_data.cpp:3303
ProSHADE_internal_data::ProSHADE_data::interpolateMapFromSpheres
void interpolateMapFromSpheres(ProSHADE_settings *settings, proshade_double *&densityMapRotated)
This function interpolates the density map from the sphere mapped data.
Definition: ProSHADE_overlay.cpp:1042
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:308
ProSHADE_internal_data::ProSHADE_data::computeOptimalTranslation
void computeOptimalTranslation(proshade_double eulA, proshade_double eulB, proshade_double eulG, proshade_double trsX, proshade_double trsY, proshade_double trsZ)
This function computes and saves the optimal translation vector from the already determined translati...
Definition: ProSHADE_overlay.cpp:109
ProSHADE_internal_data::ProSHADE_data::getXAxisOrigin
proshade_signed * getXAxisOrigin(void)
This function allows access to the map X axis origin value.
Definition: ProSHADE_data.cpp:3403
ProSHADE_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:3283
ProSHADE_internal_data::ProSHADE_data::getInternalMap
proshade_double *& getInternalMap(void)
This function allows access to the first map array value address.
Definition: ProSHADE_data.cpp:3433
ProSHADE_internal_data::ProSHADE_data
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:3293
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeX
proshade_double mapCOMProcessChangeX
The change in X axis between the creation of the structure (originalMapXCom) and just before rotation...
Definition: ProSHADE_data.hpp:97
ProSHADE_internal_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:1228
ProSHADE_internal_mapManip::moveMapByFourier
void moveMapByFourier(proshade_double *&map, proshade_single xMov, proshade_single yMov, proshade_single zMov, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
Function for moving map back to original PDB location by using Fourier transformation.
Definition: ProSHADE_mapManip.cpp:792
ProSHADE_internal_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:657
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:471
ProSHADE_internal_data::ProSHADE_data::getXFromPtr
proshade_signed * getXFromPtr(void)
This function allows access to the map start along the X axis.
Definition: ProSHADE_data.cpp:3343
ProSHADE_settings::verbose
proshade_signed verbose
Should the software report on the progress, or just be quiet? Value between -1 (nothing) and 4 (loud)
Definition: ProSHADE_settings.hpp:188
ProSHADE_internal_misc::addToDoubleVector
void addToDoubleVector(std::vector< proshade_double > *vecToAddTo, proshade_double elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:77
ProSHADE_internal_data::ProSHADE_data::getXDim
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
Definition: ProSHADE_data.cpp:3313
ProSHADE_internal_data::ProSHADE_data::rotateMap
void rotateMap(ProSHADE_settings *settings, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function rotates a map based on the given Euler angles.
Definition: ProSHADE_overlay.cpp:736
ProSHADE_internal_data::ProSHADE_data::mapToSpheres
void mapToSpheres(ProSHADE_settings *settings)
This function converts the internal map onto a set of concentric spheres.
Definition: ProSHADE_data.cpp:1564
ProSHADE_internal_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:260
ProSHADE_overlay.hpp
This header file declares the functions required for the structure overlay computation.
ProSHADE_internal_data::ProSHADE_data::mapMovFromsChangeY
proshade_double mapMovFromsChangeY
When the map is translated, the yFrom and yTo values are changed. This variable holds how much they h...
Definition: ProSHADE_data.hpp:95
ProSHADE_internal_peakSearch::getBestPeakEulerAngsNaive
void getBestPeakEulerAngsNaive(proshade_complex *map, proshade_unsign dim, proshade_double *eulA, proshade_double *eulB, proshade_double *eulG, ProSHADE_settings *settings)
This function finds the highest peaks optimised Euler angles using the "naive" approach.
Definition: ProSHADE_peakSearch.cpp:351
ProSHADE_internal_data::ProSHADE_data::computeTranslationMap
void computeTranslationMap(ProSHADE_internal_data::ProSHADE_data *obj1)
This function does the computation of the translation map and saves results internally.
Definition: ProSHADE_overlay.cpp:389
ProSHADE_internal_data::ProSHADE_data::xFrom
proshade_signed xFrom
This is the starting index along the x axis.
Definition: ProSHADE_data.hpp:110
ProSHADE_internal_data::ProSHADE_data::getMaxBand
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
Definition: ProSHADE_data.cpp:3004
ProSHADE_internal_data::ProSHADE_data::invertSHCoefficients
void invertSHCoefficients(void)
This function computes the shell mapped data from inverting the Spherical Harmonics coefficients.
Definition: ProSHADE_overlay.cpp:948
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:3323
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:89
ProSHADE_internal_data::ProSHADE_data::processInternalMap
void processInternalMap(ProSHADE_settings *settings)
This function simply clusters several other functions which should be called together.
Definition: ProSHADE_data.cpp:1464
ProSHADE_internal_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:515
ProSHADE_internal_data::ProSHADE_data::getTranslationFnPointer
proshade_complex * getTranslationFnPointer(void)
This function allows access to the translation function through a pointer.
Definition: ProSHADE_data.cpp:3443
ProSHADE_internal_data::ProSHADE_data::getZDim
proshade_unsign getZDim(void)
This function allows access to the map size in indices along the Z axis.
Definition: ProSHADE_data.cpp:3333
ProSHADE_internal_overlay::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:1018
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:902
ProSHADE_internal_overlay::findHighestValueInMap
void findHighestValueInMap(fftw_complex *resIn, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ, proshade_double *mapPeak)
This function simply finds the highest value in fftw_complex map and returns its position and value.
Definition: ProSHADE_overlay.cpp:558
ProSHADE_internal_data::ProSHADE_data::mapMovFromsChangeZ
proshade_double mapMovFromsChangeZ
When the map is translated, the zFrom and zTo values are changed. This variable holds how much they h...
Definition: ProSHADE_data.hpp:96
ProSHADE_internal_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:859
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:2994
ProSHADE_internal_data::ProSHADE_data::getYToPtr
proshade_signed * getYToPtr(void)
This function allows access to the map last position along the Y axis.
Definition: ProSHADE_data.cpp:3383
ProSHADE_internal_data::ProSHADE_data::getZAxisOrigin
proshade_signed * getZAxisOrigin(void)
This function allows access to the map Z axis origin value.
Definition: ProSHADE_data.cpp:3423
ProSHADE_internal_data::ProSHADE_data::translateMap
void translateMap(ProSHADE_settings *settings, 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:795
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:65
ProSHADE_internal_data::ProSHADE_data::getXToPtr
proshade_signed * getXToPtr(void)
This function allows access to the map last position along the X axis.
Definition: ProSHADE_data.cpp:3373
ProSHADE_internal_data::ProSHADE_data::getYFromPtr
proshade_signed * getYFromPtr(void)
This function allows access to the map start along the Y axis.
Definition: ProSHADE_data.cpp:3353
ProSHADE_internal_data::ProSHADE_data::mapCOMProcessChangeY
proshade_double mapCOMProcessChangeY
The change in Y axis between the creation of the structure (originalMapYCom) and just before rotation...
Definition: ProSHADE_data.hpp:98
ProSHADE_internal_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:432
ProSHADE_internal_data::ProSHADE_data::getInvSO3Coeffs
proshade_complex * getInvSO3Coeffs(void)
This function allows access to the inverse SO(3) coefficients array.
Definition: ProSHADE_data.cpp:3231
ProSHADE_internal_data::ProSHADE_data::yFrom
proshade_signed yFrom
This is the starting index along the y axis.
Definition: ProSHADE_data.hpp:111
ProSHADE_settings::inputFiles
std::vector< std::string > inputFiles
This vector contains the filenames of all input structure files.
Definition: ProSHADE_settings.hpp:95
ProSHADE_internal_maths::complexMultiplicationConjug
void complexMultiplicationConjug(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2, proshade_double *retReal, proshade_double *retImag)
Function to multiply two complex numbers by using the second number's conjugate.
Definition: ProSHADE_maths.cpp:62
ProSHADE_internal_data::ProSHADE_data::getYAxisOrigin
proshade_signed * getYAxisOrigin(void)
This function allows access to the map Y axis origin value.
Definition: ProSHADE_data.cpp:3413
ProSHADE_internal_data::ProSHADE_data::computeRotatedSH
void computeRotatedSH(ProSHADE_settings *settings)
This function multiplies the objects spherical harmonics with the Wigner D matrices,...
Definition: ProSHADE_overlay.cpp:848
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:686
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:497
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:706
ProSHADE_internal_messages::printProgressMessage
void printProgressMessage(proshade_signed verbose, proshade_signed messageLevel, std::string message)
General stdout message printing.
Definition: ProSHADE_messages.cpp:70
ProSHADE_internal_mapManip::moveMapByIndices
void moveMapByIndices(proshade_single *xMov, proshade_single *yMov, proshade_single *zMov, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed *xFrom, proshade_signed *xTo, proshade_signed *yFrom, proshade_signed *yTo, proshade_signed *zFrom, proshade_signed *zTo, proshade_signed *xOrigin, proshade_signed *yOrigin, proshade_signed *zOrigin)
Function for moving map back to original PDB location by changing the indices.
Definition: ProSHADE_mapManip.cpp:744
ProSHADE_internal_data::ProSHADE_data::zeroPaddToDims
void zeroPaddToDims(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim)
This function changes the size of a structure to fit the supplied new limits.
Definition: ProSHADE_overlay.cpp:598
ProSHADE_internal_data::ProSHADE_data::readInStructure
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
Definition: ProSHADE_data.cpp:447
ProSHADE_internal_data::ProSHADE_data::getZFromPtr
proshade_signed * getZFromPtr(void)
This function allows access to the map start along the Z axis.
Definition: ProSHADE_data.cpp:3363
ProSHADE_internal_data::ProSHADE_data::getZToPtr
proshade_signed * getZToPtr(void)
This function allows access to the map last position along the Z axis.
Definition: ProSHADE_data.cpp:3393
ProSHADE_internal_data::ProSHADE_data::allocateRotatedSHMemory
void allocateRotatedSHMemory(ProSHADE_settings *settings)
This function allocates the memory required for storing the rotated Spherical Harmonics coefficients.
Definition: ProSHADE_overlay.cpp:818