ProSHADE  0.7.5.4 (MAR 2021)
Protein Shape Detection
ProSHADE_distances.cpp
Go to the documentation of this file.
1 
23 //==================================================== ProSHADE
24 #include "ProSHADE_distances.hpp"
25 
34 {
35  //================================================ Allocate the required memory
36  this->rrpMatrices = new proshade_double** [this->maxShellBand];
37  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices, __FILE__, __LINE__, __func__ );
38 
39  for ( proshade_unsign bwIt = 0; bwIt < this->maxShellBand; bwIt++ )
40  {
41  //============================================ For rach sphere
42  this->rrpMatrices[bwIt] = new proshade_double* [this->noSpheres];
43  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices[bwIt], __FILE__, __LINE__, __func__ );
44 
45  for ( proshade_unsign shIt = 0; shIt < this->noSpheres; shIt++ )
46  {
47  this->rrpMatrices[bwIt][shIt] = new double [this->noSpheres];
48  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices[bwIt][shIt], __FILE__, __LINE__, __func__ );
49  }
50  }
51 }
52 
62 {
63  //================================================ Report progress
64  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Computing RRP matrices for structure " + this->fileName );
65 
66  //================================================ Allocate the memory
67  this->allocateRRPMemory ( settings );
68 
69  //================================================ Start computation: For each band (l)
70  proshade_double descValR = 0.0;
71  proshade_unsign arrPos1, arrPos2;
72  for ( proshade_unsign band = 0; band < this->maxShellBand; band++ )
73  {
74  //============================================ For each unique shell couple
75  for ( proshade_unsign shell1 = 0; shell1 < this->noSpheres; shell1++ )
76  {
77  //======================================== Does the band exist for this shell1?
78  if ( !ProSHADE_internal_distances::isBandWithinShell ( band, shell1, this->spheres ) )
79  {
80  for ( proshade_unsign shell2 = 0; shell2 < this->noSpheres; shell2++ )
81  {
82  this->rrpMatrices[band][shell1][shell2] = 0.0;
83  this->rrpMatrices[band][shell2][shell1] = 0.0;
84  }
85  continue;
86  }
87 
88  for ( proshade_unsign shell2 = 0; shell2 < this->noSpheres; shell2++ )
89  {
90  //==================================== Compute each values only once
91  if ( shell1 > shell2 ) { continue; }
92 
93  //==================================== Check if band exists for this shell2?
94  if ( !ProSHADE_internal_distances::isBandWithinShell ( band, shell2, this->spheres ) )
95  {
96  this->rrpMatrices[band][shell1][shell2] = 0.0;
97  this->rrpMatrices[band][shell2][shell1] = 0.0;
98  continue;
99  }
100 
101  //==================================== Initialise
102  descValR = 0.0;
103 
104  //==================================== Sum over order (m)
105  for ( proshade_unsign order = 0; order < static_cast<proshade_unsign> ( ( 2 * band ) + 1 ); order++ )
106  {
107  arrPos1 = static_cast<proshade_unsign> ( seanindex ( static_cast<proshade_signed> ( order ) -
108  static_cast<proshade_signed> ( band ),
109  band, this->spheres[shell1]->getLocalBandwidth() ) );
110  arrPos2 = static_cast<proshade_unsign> ( seanindex ( static_cast<proshade_signed> ( order ) -
111  static_cast<proshade_signed> ( band ),
112  band, this->spheres[shell2]->getLocalBandwidth() ) );
113  descValR += ProSHADE_internal_maths::complexMultiplicationConjugRealOnly ( &this->sphericalHarmonics[shell1][arrPos1][0],
114  &this->sphericalHarmonics[shell1][arrPos1][1],
115  &this->sphericalHarmonics[shell2][arrPos2][0],
116  &this->sphericalHarmonics[shell2][arrPos2][1] );
117  }
118 
119  //==================================== Save the matrices
120  this->rrpMatrices[band][shell1][shell2] = descValR;
121  this->rrpMatrices[band][shell2][shell1] = descValR;
122  }
123  }
124  }
125 
126  //================================================ Report progress
127  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "RRP matrices successfully computed." );
128 
129  //================================================ Done
130  return ;
131 
132 }
133 
143 bool ProSHADE_internal_distances::isBandWithinShell ( proshade_unsign bandInQuestion, proshade_unsign shellInQuestion, ProSHADE_internal_spheres::ProSHADE_sphere** spheres )
144 {
145  if ( bandInQuestion < spheres[shellInQuestion]->getLocalBandwidth() )
146  {
147  return ( true );
148  }
149  else
150  {
151  return ( false );
152  }
153 }
154 
166 {
167  //================================================ Report starting the task
168  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting energy levels distance computation." );
169 
170  //================================================ Initialise local variables
171  proshade_double ret = 0.0;
172  std::vector<proshade_double> bandDists;
173 
174  //================================================ Sanity check
175  if ( !settings->computeEnergyLevelsDesc )
176  {
177  throw ProSHADE_exception ( "Attempted computing energy levels descriptors when it was not required.", "ED00017", __FILE__, __LINE__, __func__, "Attempted to pre-compute the RRP matrices, when the user\n : has specifically stated that these should not be computed.\n : Unless you manipulated the code, this error should never\n : occur; if you see this, I made a large blunder. Please let\n : me know!" );
178  }
179 
180  //================================================ Get the RRP matrices for both objects
181  obj1->computeRRPMatrices ( settings );
182  obj2->computeRRPMatrices ( settings );
183 
184  //================================================ Find the minimium comparable shells and bands
185  proshade_unsign minCommonShells = std::min ( obj1->getMaxSpheres(), obj2->getMaxSpheres() );
186  proshade_unsign minCommonBands = std::min ( obj1->getMaxBand(), obj2->getMaxBand() );
187 
188  //================================================ Get the Pearson's coefficients for each common band
189  computeRRPPearsonCoefficients ( obj1, obj2, settings, minCommonBands, minCommonShells, &bandDists );
190 
191  //================================================ Get distance (by averaging Patterson's coefficients)
192  ret = static_cast<proshade_double> ( std::accumulate ( bandDists.begin(), bandDists.end(), 0.0 ) ) /
193  static_cast<proshade_double> ( bandDists.size() );
194 
195  //================================================ Report completion
196  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Energy levels distance computation complete." );
197 
198  //================================================ Done
199  return ( ret );
200 
201 }
202 
215 void ProSHADE_internal_distances::computeRRPPearsonCoefficients ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, ProSHADE_settings* settings, proshade_unsign minCommonBands, proshade_unsign minCommonShells, std::vector<proshade_double>* bandDists )
216 {
217  //================================================ Report completion
218  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Correlating RRP matrices." );
219 
220  //================================================ Initialise local variables
221  proshade_double *str1Vals = new proshade_double[minCommonShells * minCommonShells];
222  proshade_double *str2Vals = new proshade_double[minCommonShells * minCommonShells];
223  ProSHADE_internal_misc::checkMemoryAllocation ( str1Vals, __FILE__, __LINE__, __func__ );
224  ProSHADE_internal_misc::checkMemoryAllocation ( str2Vals, __FILE__, __LINE__, __func__ );
225  proshade_unsign arrIter = 0;
226 
227  //================================================ Start computation: For each band (l)
228  for ( proshade_unsign band = 0; band < minCommonBands; band++ )
229  {
230  //============================================ Reset local counter
231  arrIter = 0;
232 
233  //============================================ For each shell pair
234  for ( proshade_unsign shell1 = 0; shell1 < minCommonShells; shell1++ )
235  {
236  //======================================== Check if band exists (progressive only)
237  if ( settings->progressiveSphereMapping ) { if ( !obj1->shellBandExists( shell1, band ) || !obj2->shellBandExists( shell1, band ) ) { continue; } }
238 
239  for ( proshade_unsign shell2 = 0; shell2 < minCommonShells; shell2++ )
240  {
241  //============================ Check the other shell as well
242  if ( !obj1->shellBandExists( shell2, band ) || !obj2->shellBandExists( shell2, band ) ) { continue; }
243 
244  //==================================== Set values between which the Person's correlation coefficient should be computed
245  str1Vals[arrIter] = obj1->getRRPValue ( band, shell1, shell2 ) *
246  pow ( static_cast<proshade_double> ( shell1 ), settings->enLevMatrixPowerWeight ) *
247  pow ( static_cast<proshade_double> ( shell2 ), settings->enLevMatrixPowerWeight );
248  str2Vals[arrIter] = obj2->getRRPValue ( band, shell1, shell2 ) *
249  pow ( static_cast<proshade_double> ( shell1 ), settings->enLevMatrixPowerWeight ) *
250  pow ( static_cast<proshade_double> ( shell2 ), settings->enLevMatrixPowerWeight );
251 
252  arrIter += 1;
253  }
254  }
255 
256  //============================================ Get Pearson's Correlation Coefficient
257  ProSHADE_internal_misc::addToDoubleVector ( bandDists, ProSHADE_internal_maths::pearsonCorrCoeff ( str1Vals, str2Vals, arrIter ) );
258  }
259 
260  //================================================ Clean up
261  delete[] str1Vals;
262  delete[] str2Vals;
263 
264  //================================================ Report completion
265  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "RRP matrices correlation computed." );
266 
267  //================================================ Done
268  return ;
269 }
270 
281 {
282  //================================================ Save the maximum band to the object
283  this->maxCompBand = band;
284 
285  //================================================ Allocate the required memory
286  this->eMatrices = new proshade_complex** [this->maxCompBand];
287  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices, __FILE__, __LINE__, __func__ );
288 
289  for ( proshade_unsign bandIter = 0; bandIter < this->maxCompBand; bandIter++ )
290  {
291  //============================================ Allocate the data structure
292  this->eMatrices[bandIter] = new proshade_complex* [static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 )];
293  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices[bandIter], __FILE__, __LINE__, __func__ );
294 
295  for ( proshade_unsign band2Iter = 0; band2Iter < static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 ); band2Iter++ )
296  {
297  this->eMatrices[bandIter][band2Iter] = new proshade_complex [static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 )];
298  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices[bandIter][band2Iter], __FILE__, __LINE__, __func__ );
299  }
300  }
301 
302  //================================================ Done
303  return ;
304 
305 }
306 
317 void ProSHADE_internal_distances::allocateTrSigmaWorkspace ( proshade_unsign minSpheres, proshade_unsign intOrder, proshade_double*& obj1Vals, proshade_double*& obj2Vals, proshade_double*& GLabscissas, proshade_double*& GLweights, proshade_complex*& radiiVals )
318 {
319  //================================================ Allocate the memory
320  obj1Vals = new proshade_double [minSpheres];
321  obj2Vals = new proshade_double [minSpheres];
322  radiiVals = new proshade_complex[minSpheres];
323  GLabscissas = new proshade_double [intOrder];
324  GLweights = new proshade_double [intOrder];
325 
326  //================================================ Check the memory allocation
327  ProSHADE_internal_misc::checkMemoryAllocation ( obj1Vals, __FILE__, __LINE__, __func__ );
328  ProSHADE_internal_misc::checkMemoryAllocation ( obj2Vals, __FILE__, __LINE__, __func__ );
329  ProSHADE_internal_misc::checkMemoryAllocation ( radiiVals, __FILE__, __LINE__, __func__ );
330  ProSHADE_internal_misc::checkMemoryAllocation ( GLabscissas, __FILE__, __LINE__, __func__ );
331  ProSHADE_internal_misc::checkMemoryAllocation ( GLweights, __FILE__, __LINE__, __func__ );
332 
333  //================================================ Done
334  return ;
335 
336 }
337 
346 void ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude ( ProSHADE_internal_data::ProSHADE_data* obj, proshade_unsign band, proshade_unsign order, proshade_unsign radius, proshade_double* result )
347 {
348  //================================================ Find the magnitude
350  obj->getImagSphHarmValue ( band, order, radius ),
351  obj->getRealSphHarmValue ( band, order, radius ),
352  obj->getImagSphHarmValue ( band, order, radius ) );
353 
354  //================================================ Weight by radius^2 for the integration that will follow
355  *result *= pow ( static_cast<proshade_double> ( obj->getAnySphereRadius( radius ) ), 2.0 );
356 
357  //================================================ Done
358  return ;
359 
360 }
361 
375 void ProSHADE_internal_distances::computeEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_complex* radiiVals, proshade_unsign integOrder, proshade_double* abscissas, proshade_double* weights, proshade_double integRange, proshade_double sphereDist )
376 {
377  //================================================ Initialise local variables
378  proshade_unsign objCombValsIter = 0;
379  proshade_double hlpReal, hlpImag;
380  proshade_complex arrVal;
381 
382  //================================================ For each combination of m and m' for E matrices
383  for ( proshade_unsign order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
384  {
385  //============================================ Reset loop
386  objCombValsIter = 0;
387 
388  //============================================ Find the c*conj(c) values for different radii
389  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
390  {
391 
392  //======================================== Get only values where the shell has the band
393  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= bandIter ) { continue; }
394 
395  //======================================== Multiply coeffs
396  ProSHADE_internal_maths::complexMultiplicationConjug ( obj1->getRealSphHarmValue ( bandIter, orderIter, radiusIter ),
397  obj1->getImagSphHarmValue ( bandIter, orderIter, radiusIter ),
398  obj2->getRealSphHarmValue ( bandIter, order2Iter, radiusIter ),
399  obj2->getImagSphHarmValue ( bandIter, order2Iter, radiusIter ),
400  &hlpReal, &hlpImag );
401 
402  //======================================== Apply r^2 integral weight
403  radiiVals[objCombValsIter][0] = hlpReal * pow ( ( static_cast<proshade_double> ( obj1->getAnySphereRadius( radiusIter ) ) ), 2.0 );
404  radiiVals[objCombValsIter][1] = hlpImag * pow ( ( static_cast<proshade_double> ( obj1->getAnySphereRadius( radiusIter ) ) ), 2.0 );
405 
406  objCombValsIter += 1;
407  }
408 
409  //============================================ Integrate over all radii using n-point Gauss-Legendre integration
410  ProSHADE_internal_maths::gaussLegendreIntegration ( radiiVals, objCombValsIter, integOrder, abscissas, weights, integRange, sphereDist, &hlpReal, &hlpImag );
411 
412  //============================================ Save the result into E matrices
413  arrVal[0] = hlpReal;
414  arrVal[1] = hlpImag;
415  obj2->setEMatrixValue ( bandIter, orderIter, order2Iter, arrVal );
416  }
417 
418  //================================================ Done
419  return ;
420 
421 }
422 
438 proshade_double ProSHADE_internal_distances::computeWeightsForEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_double* obj1Vals, proshade_double* obj2Vals, proshade_unsign integOrder, proshade_double* abscissas, proshade_double* weights, proshade_double sphereDist )
439 {
440  //================================================ Initialise local values
441  proshade_unsign obj1ValsIter = 0;
442  proshade_unsign obj2ValsIter = 0;
443 
444  //================================================ Set sphere counters
445  proshade_unsign minSphere = std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() );
446  proshade_unsign maxSphere = 0;
447 
448  //================================================ For each radius, deal with weights
449  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
450  {
451  //============================================ Get only values where the shell has the band
452  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= bandIter ) { continue; }
453  minSphere = std::min ( radiusIter, minSphere );
454  maxSphere = std::max ( radiusIter, maxSphere );
455 
456  //============================================ Get the magnitudes for weighting
457  computeSphericalHarmonicsMagnitude ( obj1, bandIter, orderIter, radiusIter, &(obj1Vals[obj1ValsIter]) );
458  computeSphericalHarmonicsMagnitude ( obj2, bandIter, orderIter, radiusIter, &(obj2Vals[obj2ValsIter]) );
459  obj1ValsIter += 1;
460  obj2ValsIter += 1;
461  }
462 
463  //================================================ Integrate weights
464  proshade_double minSphereRad = obj1->getSpherePosValue ( minSphere ) - ( sphereDist * 0.5 );
465  proshade_double maxSphereRad = obj1->getSpherePosValue ( maxSphere ) + ( sphereDist * 0.5 );
466 
467  obj1->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj1Vals, obj1ValsIter, integOrder, abscissas, weights, maxSphereRad - minSphereRad, sphereDist ) );
468  obj2->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj2Vals, obj2ValsIter, integOrder, abscissas, weights, maxSphereRad - minSphereRad, sphereDist ) );
469 
470  //================================================ Done
471  return ( maxSphereRad - minSphereRad );
472 
473 }
474 
483 void ProSHADE_internal_distances::releaseTrSigmaWorkspace ( proshade_double*& obj1Vals, proshade_double*& obj2Vals, proshade_double*& GLabscissas, proshade_double*& GLweights, proshade_complex*& radiiVals )
484 {
485  //================================================ Release memory
486  delete[] obj1Vals;
487  delete[] obj2Vals;
488  delete[] radiiVals;
489  delete[] GLabscissas;
490  delete[] GLweights;
491 
492  //================================================ Set to NULL
493  obj1Vals = NULL;
494  obj2Vals = NULL;
495  radiiVals = NULL;
496  GLabscissas = NULL;
497  GLweights = NULL;
498 
499  //================================================ Done
500  return ;
501 
502 }
503 
516 {
517  //================================================ Report progress
518  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Starting computation of E matrices." );
519 
520  //================================================ Allocatre memory for E matrices in the second object (first may be compared to more structures and therefore its data would be written over)
521  obj2->allocateEMatrices ( settings, std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) );
522 
523  //================================================ Initialise local variables
524  proshade_double *obj1Vals, *obj2Vals, *GLAbscissas, *GLWeights;
525  proshade_complex* radiiVals;
526  proshade_double integRange;
527 
528  //================================================ Allocate workspace memory
529  allocateTrSigmaWorkspace ( std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ), settings->integOrder, obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals);
530 
531  //================================================ Initialise abscissas and weights for integration
532  ProSHADE_internal_maths::getLegendreAbscAndWeights ( settings->integOrder, GLAbscissas, GLWeights, settings->taylorSeriesCap );
533 
534  //================================================ For each band (l), compute the E matrix integrals
535  for ( proshade_unsign bandIter = 0; bandIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); bandIter++ )
536  {
537  //============================================ For each order (m)
538  for ( proshade_unsign orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
539  {
540  //======================================== Get weights for the required band(l) and order (m)
541  integRange = computeWeightsForEMatricesForLM ( obj1, obj2, bandIter, orderIter, obj1Vals, obj2Vals, settings->integOrder, GLAbscissas, GLWeights, settings->maxSphereDists );
542 
543  //======================================== Compute E matrices value for given band (l) and order(m)
544  computeEMatricesForLM ( obj1, obj2, bandIter, orderIter, radiiVals, settings->integOrder, GLAbscissas, GLWeights, integRange, settings->maxSphereDists );
545  }
546 
547  //============================================ Report progress
548  if ( settings->verbose > 3 )
549  {
550  std::stringstream hlpSS;
551  hlpSS << "E matrices computed for band " << bandIter;
552  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, hlpSS.str() );
553  }
554  }
555 
556  //================================================ Release the workspace memory
557  releaseTrSigmaWorkspace ( obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals );
558 
559  //================================================ Report progress
560  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices computed." );
561 
562  //================================================ Done
563  return ;
564 
565 }
566 
578 {
579  //================================================ Report progress
580  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Starting E matrices normalisation." );
581 
582  //================================================ Normalise by the Pearson's c.c. like formula
583  proshade_double eMatNormFactor = std::sqrt ( obj1->getIntegrationWeight() * obj2->getIntegrationWeight() );
584 
585  //================================================ If this is self-correlation (i.e. obj1 == obj2), then divide normalisation factor by 2 as the weight was applied cumulatively!
586  if ( settings->task == Symmetry ) { eMatNormFactor /= 2.0; }
587 
588  for ( proshade_unsign bandIter = 0; bandIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); bandIter++ )
589  {
590  //============================================ For each combination of m and m' for E matrices
591  for ( proshade_unsign orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
592  {
593  for ( proshade_unsign order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
594  {
595  obj2->normaliseEMatrixValue ( bandIter, orderIter, order2Iter, eMatNormFactor );
596  }
597  }
598  }
599 
600  //================================================ Report progress
601  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, "E matrices normalised." );
602 
603  //================================================ Done
604  return ;
605 
606 }
607 
622 {
623  //================================================ Report starting the task
624  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting trace sigma distance computation." );
625 
626  //================================================ Initialise return variable
627  proshade_double ret = 0.0;
628 
629  //================================================ Sanity check
630  if ( !settings->computeTraceSigmaDesc )
631  {
632  throw ProSHADE_exception ( "Attempted computing trace sigma descriptors when it was\n : not required.", "ED00018", __FILE__, __LINE__, __func__, "Attempted to pre-compute the E matrices, when the user\n : has specifically stated that these should not be computed.\n : Unless you manipulated the code, this error should never\n : occur; if you see this, I made a large blunder. Please let\n : me know!" );
633  }
634 
635  //================================================ Empty the cumulative weights back to 0.0 for each structure
636  obj1->setIntegrationWeight ( 0.0 );
637  obj1->setIntegrationWeight ( 0.0 );
638 
639  //================================================ Compute un-weighted E matrices and their weights
640  computeEMatrices ( obj1, obj2, settings );
641 
642  //================================================ Normalise E matrices by the magnitudes
643  normaliseEMatrices ( obj1, obj2, settings );
644 
645  //================================================ Allocate the required memory
646  double* singularValues = new double[( ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) * 2 ) + 1 )];
647  ProSHADE_internal_misc::checkMemoryAllocation ( singularValues, __FILE__, __LINE__, __func__ );
648 
649  //================================================ Compute the distance
650  for ( proshade_unsign lIter = 0; lIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); lIter++ )
651  {
652  //============================================ Find the complex matrix SVD singular values
653  ProSHADE_internal_maths::complexMatrixSVDSigmasOnly ( obj2->getEMatrixByBand ( lIter ), static_cast<int> ( ( lIter * 2 ) + 1 ), singularValues );
654 
655  //============================================ Now sum the trace
656  for ( proshade_unsign iter = 0; iter < ( ( lIter * 2 ) + 1 ); iter++ )
657  {
658  ret += singularValues[iter];
659  }
660  }
661 
662  //================================================ Report completion
663  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices decomposed to singular values." );
664 
665  //================================================ Release the memory
666  delete[] singularValues;
667 
668  //================================================ Report completion
669  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Trace sigma distance computation complete." );
670 
671  //================================================ Done
672  return ( ret );
673 
674 }
675 
681 {
682  //================================================ Allocate the memory
683  this->so3Coeffs = new fftw_complex [static_cast<proshade_unsign>( ( 4 * pow( static_cast<proshade_double> ( band ), 3.0 ) - static_cast<proshade_double> ( band ) ) / 3.0 )];
684  this->so3CoeffsInverse = new fftw_complex [static_cast<proshade_unsign>( pow( static_cast<proshade_double> ( band ) * 2.0, 3.0 ) )];
685 
686  //================================================ Check memory allocation
687  ProSHADE_internal_misc::checkMemoryAllocation ( this->so3Coeffs, __FILE__, __LINE__, __func__ );
688  ProSHADE_internal_misc::checkMemoryAllocation ( this->so3CoeffsInverse, __FILE__, __LINE__, __func__ );
689 
690  //================================================ Done
691  return ;
692 
693 }
694 
707 {
708  //================================================ Report progress
709  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Converting E matrices to SO(3) coefficients." );
710 
711  //================================================ Allocate memory for the coefficients
712  obj2->allocateSO3CoeffsSpace ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) );
713 
714  //================================================ Initialise local variables
715  proshade_double wigNorm, hlpValReal, hlpValImag;
716  proshade_double signValue = 1.0;
717  proshade_unsign indexO;
718  proshade_complex hlpVal;
719 
720  //================================================ For each band (l)
721  for ( proshade_signed bandIter = 0; bandIter < static_cast<proshade_signed> ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ); bandIter++ )
722  {
723  //============================================ Get wigner normalisation factor
724  wigNorm = 2.0 * M_PI * sqrt ( 2.0 / (2.0 * bandIter + 1.0 ) ) ;
725 
726  //============================================ For each order (m)
727  for ( proshade_signed orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
728  {
729  //======================================== Set the sign
730  if ( ( orderIter - bandIter + bandIter ) % 2 ) { signValue = -1.0 ; }
731  else { signValue = 1.0 ; }
732 
733  //======================================== For each order2 (m')
734  for ( proshade_signed order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
735  {
736  //==================================== Find output index
737  indexO = static_cast<proshade_unsign> ( so3CoefLoc ( orderIter - bandIter, order2Iter - bandIter, bandIter, std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ) );
738 
739  //==================================== Compute and save the SO(3) coefficients
740  obj2->getEMatrixValue ( bandIter, orderIter, order2Iter, &hlpValReal, &hlpValImag );
741  hlpVal[0] = hlpValReal * wigNorm * signValue;
742  hlpVal[1] = hlpValImag * wigNorm * signValue;
743  obj2->setSO3CoeffValue ( indexO, hlpVal );
744 
745  //==================================== Switch the sign value
746  signValue *= -1.0;
747  }
748  }
749  }
750 
751  //================================================ Report progress
752  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "SO(3) coefficients obtained." );
753 
754  //================================================ Done
755  return ;
756 
757 }
758 
766 void ProSHADE_internal_distances::allocateInvSOFTWorkspaces ( proshade_complex*& work1, proshade_complex*& work2, proshade_double*& work3, proshade_unsign band )
767 {
768  //================================================ Allocate memory
769  work1 = new proshade_complex[8 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 3.0 ) )];
770  work2 = new proshade_complex[14 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (48 * band)];
771  work3 = new proshade_double [2 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (24 * band)];
772 
773  //================================================ Check the memory allocation
774  ProSHADE_internal_misc::checkMemoryAllocation ( work1, __FILE__, __LINE__, __func__ );
775  ProSHADE_internal_misc::checkMemoryAllocation ( work2, __FILE__, __LINE__, __func__ );
776  ProSHADE_internal_misc::checkMemoryAllocation ( work3, __FILE__, __LINE__, __func__ );
777 
778  //================================================ Done
779  return ;
780 
781 }
782 
790 void ProSHADE_internal_distances::prepareInvSOFTPlan ( fftw_plan* inverseSO3, proshade_unsign band, fftw_complex* work1, proshade_complex* invCoeffs )
791 {
792  //================================================ Prepare the plan describing variables
793  int howmany = 4 * band * band;
794  int idist = 2 * band;
795  int odist = 2 * band;
796  int rank = 2;
797 
798  int inembed[2], onembed[2];
799  inembed[0] = 2 * band;
800  inembed[1] = 4 * band * band;
801  onembed[0] = 2 * band;
802  onembed[1] = 4 * band * band;
803 
804  int istride = 1;
805  int ostride = 1;
806 
807  int na[2];
808  na[0] = 1;
809  na[1] = 2 * band;
810 
811  //================================================ Create the plan
812  *inverseSO3 = fftw_plan_many_dft ( rank,
813  na,
814  howmany,
815  work1,
816  inembed,
817  istride,
818  idist,
819  invCoeffs,
820  onembed,
821  ostride,
822  odist,
823  FFTW_FORWARD,
824  FFTW_ESTIMATE );
825 
826  //================================================ Done
827  return ;
828 
829 }
830 
837 void ProSHADE_internal_distances::releaseInvSOFTMemory ( proshade_complex*& work1, proshade_complex*& work2, proshade_double*& work3 )
838 {
839  //================================================ Release memory
840  delete[] work1;
841  delete[] work2;
842  delete[] work3;
843 
844  //================================================ Done
845  return ;
846 }
847 
860 {
861  //================================================ Report progress
862  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Computing inverse SO(3) Fourier transform." );
863 
864  //================================================ Initialise local variables
865  proshade_complex *workspace1, *workspace2;
866  proshade_double *workspace3;
867  fftw_plan inverseSO3;
868 
869  //================================================ Allocate memory for the workspaces
870  allocateInvSOFTWorkspaces ( workspace1, workspace2, workspace3, std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) );
871 
872  //================================================ Prepare the FFTW plan
873  prepareInvSOFTPlan ( &inverseSO3, std::min ( obj1->getMaxBand(), obj2->getMaxBand() ), workspace1, obj2->getInvSO3Coeffs ( ) );
874 
875  //================================================ Compute the transform
876  Inverse_SO3_Naive_fftw ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ),
877  obj2->getSO3Coeffs ( ),
878  obj2->getInvSO3Coeffs ( ),
879  workspace1,
880  workspace2,
881  workspace3,
882  &inverseSO3,
883  0 );
884 
885  //================================================ Release memory
886  releaseInvSOFTMemory ( workspace1, workspace2, workspace3 );
887  fftw_destroy_plan ( inverseSO3 );
888 
889  //================================================ Report progress
890  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Inverse SO(3) Fourier transform computed." );
891 
892  //================================================ Done
893  return ;
894 
895 }
896 
912 {
913  //================================================ Report starting the task
914  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting rotation function distance computation." );
915 
916  //================================================ Initialise return variable
917  proshade_double ret = 0.0;
918  proshade_double eulA, eulB, eulG, EMatR, EMatI, WigDR, WigDI;
919 
920  //================================================ Sanity check
921  if ( !settings->computeRotationFuncDesc )
922  {
923  throw ProSHADE_exception ( "Attempted computing rotation function descriptors when it\n : was not required.", "ED00023", __FILE__, __LINE__, __func__, "Attempted to compute the SO(3) transform and the rotation \n : function descriptor when the user did not request this. \n : Unless you manipulated the code, this error should never \n : occur; if you see this, I made a large blunder. \n : Please let me know!" );
924  }
925 
926  //================================================ Compute weighted E matrices if not already present
927  if ( !settings->computeTraceSigmaDesc )
928  {
929  computeEMatrices ( obj1, obj2, settings );
930  normaliseEMatrices ( obj1, obj2, settings );
931  }
932 
933  //================================================ Generate SO(3) coefficients
934  generateSO3CoeffsFromEMatrices ( obj1, obj2, settings );
935 
936  //================================================ Compute the inverse SO(3) Fourier Transform (SOFT) on the newly computed coefficients
937  computeInverseSOFTTransform ( obj1, obj2, settings );
938 
939  //================================================ Get inverse SO(3) map top peak Euler angle values
941  std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) * 2,
942  &eulA, &eulB, &eulG, settings );
943 
944  //================================================ Compute the Wigner D matrices for the Euler angles
945  ProSHADE_internal_wigner::computeWignerMatricesForRotation ( settings, obj2, eulA, eulB, eulG );
946 
947  //================================================ Compute the distance
948  for ( proshade_unsign bandIter = 0; bandIter < obj2->getComparisonBand(); bandIter++ )
949  {
950  //============================================ For each order1
951  for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
952  {
953  //======================================== For each order2
954  for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
955  {
956  //==================================== Multiply D_{l} * E_{l} and get sum over l of traces (i.e. just sum all together)
957  obj2->getEMatrixValue ( bandIter, order1, order2, &EMatR, &EMatI );
958  obj2->getWignerMatrixValue ( bandIter, order2, order1, &WigDR, &WigDI );
959  ret += ProSHADE_internal_maths::complexMultiplicationRealOnly ( &WigDR, &WigDI, &EMatR, &EMatI );
960  }
961  }
962  }
963 
964  //================================================ Report completion
965  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Rotation function distance computation complete." );
966 
967  //================================================ Done
968  return ( ret );
969 
970 }
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::integOrder
proshade_unsign integOrder
The order required for full Gauss-Legendre integration between the spheres.
Definition: ProSHADE_settings.hpp:119
ProSHADE_internal_data::ProSHADE_data::getEMatrixByBand
proshade_complex ** getEMatrixByBand(proshade_unsign band)
This function allows access to E matrix for a particular band.
Definition: ProSHADE_data.cpp:3201
ProSHADE_internal_maths::gaussLegendreIntegration
void gaussLegendreIntegration(proshade_complex *vals, proshade_unsign valsSize, proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_double integralOverRange, proshade_double maxSphereDists, proshade_double *retReal, proshade_double *retImag)
Function to compute the complete complex Gauss-Legendre integration over spherical harmonic values in...
Definition: ProSHADE_maths.cpp:709
ProSHADE_settings::computeTraceSigmaDesc
bool computeTraceSigmaDesc
If true, the trace sigma descriptor will be computed, otherwise all its computations will be omitted.
Definition: ProSHADE_settings.hpp:160
ProSHADE_distances.hpp
This is the header file containing declarations of functions required for computation of shape distan...
ProSHADE_internal_maths::gaussLegendreIntegrationReal
proshade_double gaussLegendreIntegrationReal(proshade_double *vals, proshade_unsign valsSize, proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_double integralOverRange, proshade_double maxSphereDists)
Function to compute real part of the Gauss-Legendre integration over spherical harmonic values in dif...
Definition: ProSHADE_maths.cpp:619
ProSHADE_settings::computeRotationFuncDesc
bool computeRotationFuncDesc
If true, the rotation function descriptor will be computed, otherwise all its computations will be om...
Definition: ProSHADE_settings.hpp:161
ProSHADE_settings::taylorSeriesCap
proshade_unsign taylorSeriesCap
The max limit on the Taylor series expansion done for the abscissas of the Gauss-Legendre integration...
Definition: ProSHADE_settings.hpp:120
ProSHADE_internal_data::ProSHADE_data::getComparisonBand
proshade_unsign getComparisonBand(void)
This function allows access to the maximum band for the comparison.
Definition: ProSHADE_data.cpp:3253
ProSHADE_internal_distances::allocateTrSigmaWorkspace
void allocateTrSigmaWorkspace(proshade_unsign minSpheres, proshade_unsign intOrder, proshade_double *&obj1Vals, proshade_double *&obj2Vals, proshade_double *&GLabscissas, proshade_double *&glWeights, proshade_complex *&radiiVals)
This helper function is responsible for allocating the workspace memory required for trace sigma desc...
Definition: ProSHADE_distances.cpp:317
ProSHADE_internal_data::ProSHADE_data::setIntegrationWeight
void setIntegrationWeight(proshade_double intW)
This function allows setting the integration weight for the object.
Definition: ProSHADE_data.cpp:3453
ProSHADE_internal_distances::computeWeightsForEMatricesForLM
proshade_double computeWeightsForEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_double *obj1Vals, proshade_double *obj2Vals, proshade_unsign integOrder, proshade_double *abscissas, proshade_double *weights, proshade_double sphereDist)
This function computes the E matrix weight values for a given band and order and saves these into the...
Definition: ProSHADE_distances.cpp:438
ProSHADE_internal_data::ProSHADE_data::normaliseEMatrixValue
void normaliseEMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double normF)
This function allows normalising the E matrix value.
Definition: ProSHADE_data.cpp:3502
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
ProSHADE_internal_data::ProSHADE_data::getEMatrixValue
void getEMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double *valueReal, proshade_double *valueImag)
This function allows access to E matrix by knowing the band, order1 and order2 indices.
Definition: ProSHADE_data.cpp:3216
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::noSpheres
proshade_unsign noSpheres
The number of spheres with map projected onto them.
Definition: ProSHADE_data.hpp:119
ProSHADE_internal_distances::releaseInvSOFTMemory
void releaseInvSOFTMemory(proshade_complex *&work1, proshade_complex *&work2, proshade_double *&work3)
This function releases the memory used for computation of the inverse SOFT transform.
Definition: ProSHADE_distances.cpp:837
ProSHADE_internal_data::ProSHADE_data::getMaxSpheres
proshade_unsign getMaxSpheres(void)
This function returns the number of spheres which contain the whole object.
Definition: ProSHADE_data.cpp:2983
ProSHADE_internal_data::ProSHADE_data::setSO3CoeffValue
void setSO3CoeffValue(proshade_unsign position, proshade_complex val)
This function allows setting the SOFT coefficient values using array position and value.
Definition: ProSHADE_data.cpp:3518
ProSHADE_internal_distances::isBandWithinShell
bool isBandWithinShell(proshade_unsign bandInQuestion, proshade_unsign shellInQuestion, ProSHADE_internal_spheres::ProSHADE_sphere **spheres)
This function checks if a band is available for a given shell.
Definition: ProSHADE_distances.cpp:143
ProSHADE_internal_data::ProSHADE_data::computeRRPMatrices
void computeRRPMatrices(ProSHADE_settings *settings)
This function pre-computes the RRP matrices for a data object.
Definition: ProSHADE_distances.cpp:61
ProSHADE_internal_data::ProSHADE_data::getAnySphereRadius
proshade_double getAnySphereRadius(proshade_unsign shell)
This function allows access to the radius of any particular sphere.
Definition: ProSHADE_data.cpp:3154
ProSHADE_settings::progressiveSphereMapping
bool progressiveSphereMapping
If true, each shell will have its own angular resolution dependent on the actual number of map points...
Definition: ProSHADE_settings.hpp:152
ProSHADE_internal_data::ProSHADE_data::getSO3Coeffs
proshade_complex * getSO3Coeffs(void)
This function allows access to the SO(3) coefficients array.
Definition: ProSHADE_data.cpp:3242
ProSHADE_settings::maxSphereDists
proshade_single maxSphereDists
The distance between spheres in spherical mapping for the largest sphere.
Definition: ProSHADE_settings.hpp:116
ProSHADE_internal_data::ProSHADE_data::maxShellBand
proshade_unsign maxShellBand
The maximum band for any shell of the object.
Definition: ProSHADE_data.hpp:123
ProSHADE_internal_data::ProSHADE_data::getSpherePosValue
proshade_double getSpherePosValue(proshade_unsign shell)
This function allows access to sphere positions.
Definition: ProSHADE_data.cpp:3189
ProSHADE_internal_data::ProSHADE_data::getIntegrationWeight
proshade_double getIntegrationWeight(void)
This function allows access to the integration weight for the object.
Definition: ProSHADE_data.cpp:3165
ProSHADE_internal_data::ProSHADE_data::setIntegrationWeightCumul
void setIntegrationWeightCumul(proshade_double intW)
This function allows setting the cumulative integration weight for the object.
Definition: ProSHADE_data.cpp:3467
ProSHADE_internal_spheres::ProSHADE_sphere
This class contains all inputed and derived data for a single sphere.
Definition: ProSHADE_spheres.hpp:49
ProSHADE_internal_distances::computeRotationunctionDescriptor
proshade_double computeRotationunctionDescriptor(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the rotation function descriptor value between two objects.
Definition: ProSHADE_distances.cpp:911
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_maths::complexMatrixSVDSigmasOnly
void complexMatrixSVDSigmasOnly(proshade_complex **mat, int dim, double *&singularValues)
Function to compute the complete complex matrix SVD and return only the sigmas.
Definition: ProSHADE_maths.cpp:802
ProSHADE_internal_data::ProSHADE_data::getWignerMatrixValue
void getWignerMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_double *valueReal, proshade_double *valueImag)
This function allows access to the Wigner D matrix by knowing the band, order1 and order2 indices.
Definition: ProSHADE_data.cpp:3268
ProSHADE_internal_data::ProSHADE_data::getShellBandwidth
proshade_unsign getShellBandwidth(proshade_unsign shell)
This function allows access to the bandwidth of a particular shell.
Definition: ProSHADE_data.cpp:3177
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_distances::computeRRPPearsonCoefficients
void computeRRPPearsonCoefficients(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings, proshade_unsign minCommonBands, proshade_unsign minCommonShells, std::vector< proshade_double > *bandDists)
This function gets the Pearson's coefficients or all bands between two objects.
Definition: ProSHADE_distances.cpp:215
ProSHADE_internal_distances::prepareInvSOFTPlan
void prepareInvSOFTPlan(fftw_plan *inverseSO3, proshade_unsign band, fftw_complex *work1, proshade_complex *invCoeffs)
This function prepares the FFTW plan for the inverse SO(3) transform.
Definition: ProSHADE_distances.cpp:790
ProSHADE_internal_maths::pearsonCorrCoeff
proshade_double pearsonCorrCoeff(proshade_double *valSet1, proshade_double *valSet2, proshade_unsign length)
Function for computing the Pearson's correlation coefficient.
Definition: ProSHADE_maths.cpp:244
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_internal_data::ProSHADE_data::rrpMatrices
proshade_double *** rrpMatrices
The energy levels descriptor shell correlation tables.
Definition: ProSHADE_data.hpp:126
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_settings::task
ProSHADE_Task task
This custom type variable determines which task to perfom (i.e. symmetry detection,...
Definition: ProSHADE_settings.hpp:92
ProSHADE_internal_data::ProSHADE_data::getRealSphHarmValue
proshade_double * getRealSphHarmValue(proshade_unsign band, proshade_unsign order, proshade_unsign shell)
This function allows access to the private internal real spherical harmonics values.
Definition: ProSHADE_data.cpp:3128
ProSHADE_internal_maths::complexMultiplicationConjugRealOnly
proshade_double complexMultiplicationConjugRealOnly(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2)
Function to conjuggate multiply two complex numbers and return the real part only.
Definition: ProSHADE_maths.cpp:103
ProSHADE_internal_data::ProSHADE_data::shellBandExists
bool shellBandExists(proshade_unsign shell, proshade_unsign bandVal)
This function checks if particular shell has a particular band.
Definition: ProSHADE_data.cpp:3030
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_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::getImagSphHarmValue
proshade_double * getImagSphHarmValue(proshade_unsign band, proshade_unsign order, proshade_unsign shell)
This function allows access to the private internal imaginary spherical harmonics values.
Definition: ProSHADE_data.cpp:3141
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::getRRPValue
proshade_double getRRPValue(proshade_unsign band, proshade_unsign sh1, proshade_unsign sh2)
This function allows access to the priva internal RRP matrices.
Definition: ProSHADE_data.cpp:3014
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_distances::computeEMatricesForLM
void computeEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, proshade_unsign bandIter, proshade_unsign orderIter, proshade_complex *radiiVals, proshade_unsign integOrder, proshade_double *abscissas, proshade_double *weights, proshade_double integRange, proshade_double sphereDist)
This function computes the E matrix un-weighted values for a given band and order and saves these int...
Definition: ProSHADE_distances.cpp:375
ProSHADE_internal_maths::getLegendreAbscAndWeights
void getLegendreAbscAndWeights(proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_unsign taylorSeriesCap)
Function to prepare abscissas and weights for Gauss-Legendre integration.
Definition: ProSHADE_maths.cpp:287
ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude
void computeSphericalHarmonicsMagnitude(ProSHADE_internal_data::ProSHADE_data *obj, proshade_unsign band, proshade_unsign order, proshade_unsign radius, proshade_double *result)
This function computes the magnitude of a particular spherical harmonics position for a given object,...
Definition: ProSHADE_distances.cpp:346
ProSHADE_internal_distances::computeTraceSigmaDescriptor
proshade_double computeTraceSigmaDescriptor(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the trace sigma descriptor value between two objects.
Definition: ProSHADE_distances.cpp:621
ProSHADE_internal_data::ProSHADE_data::allocateEMatrices
void allocateEMatrices(ProSHADE_settings *settings, proshade_unsign band)
This function allocates the required memory for the E matrices.
Definition: ProSHADE_distances.cpp:280
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_distances::releaseTrSigmaWorkspace
void releaseTrSigmaWorkspace(proshade_double *&obj1Vals, proshade_double *&obj2Vals, proshade_double *&GLabscissas, proshade_double *&glWeights, proshade_complex *&radiiVals)
This helper function is responsible for deleting the workspace memory required for trace sigma descri...
Definition: ProSHADE_distances.cpp:483
ProSHADE_internal_distances::allocateInvSOFTWorkspaces
void allocateInvSOFTWorkspaces(proshade_complex *&work1, proshade_complex *&work2, proshade_double *&work3, proshade_unsign band)
This function allocates the workspaces required to compute the inverse SOFT transform.
Definition: ProSHADE_distances.cpp:766
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_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::allocateRRPMemory
void allocateRRPMemory(ProSHADE_settings *settings)
This function allocates the required memory for the RRP matrices.
Definition: ProSHADE_distances.cpp:33
ProSHADE_internal_distances::computeEnergyLevelsDescriptor
proshade_double computeEnergyLevelsDescriptor(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the energy levels descriptor value between two objects.
Definition: ProSHADE_distances.cpp:165
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_data::ProSHADE_data::setEMatrixValue
void setEMatrixValue(proshade_unsign band, proshade_unsign order1, proshade_unsign order2, proshade_complex val)
This function allows setting the E matrix value.
Definition: ProSHADE_data.cpp:3484
ProSHADE_internal_data::ProSHADE_data::allocateSO3CoeffsSpace
void allocateSO3CoeffsSpace(proshade_unsign band)
This function allocates the memory for the SO(3) coefficients and the inverse for that calling object...
Definition: ProSHADE_distances.cpp:680
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_settings::enLevMatrixPowerWeight
proshade_double enLevMatrixPowerWeight
If RRP matrices shell position is to be weighted by putting the position as an exponent,...
Definition: ProSHADE_settings.hpp:159
ProSHADE_internal_maths::complexMultiplicationRealOnly
proshade_double complexMultiplicationRealOnly(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2)
Function to multiply two complex numbers and return the real part only.
Definition: ProSHADE_maths.cpp:83
ProSHADE_settings::computeEnergyLevelsDesc
bool computeEnergyLevelsDesc
If true, the energy levels descriptor will be computed, otherwise all its computations will be omitte...
Definition: ProSHADE_settings.hpp:158