ProSHADE  0.7.6.6 (JUL 2022)
Protein Shape Detection
ProSHADE_distances.cpp
Go to the documentation of this file.
1 
23 //==================================================== ProSHADE
24 #include "ProSHADE_distances.hpp"
25 
32 {
33  //================================================ Allocate the required memory
34  this->rrpMatrices = new proshade_double** [this->maxShellBand];
35  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices, __FILE__, __LINE__, __func__ );
36 
37  for ( proshade_unsign bwIt = 0; bwIt < this->maxShellBand; bwIt++ )
38  {
39  //============================================ For rach sphere
40  this->rrpMatrices[bwIt] = new proshade_double* [this->noSpheres];
41  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices[bwIt], __FILE__, __LINE__, __func__ );
42 
43  for ( proshade_unsign shIt = 0; shIt < this->noSpheres; shIt++ )
44  {
45  this->rrpMatrices[bwIt][shIt] = new double [this->noSpheres];
46  ProSHADE_internal_misc::checkMemoryAllocation ( this->rrpMatrices[bwIt][shIt], __FILE__, __LINE__, __func__ );
47  }
48  }
49 }
50 
60 {
61  //================================================ Report progress
62  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Computing RRP matrices for structure " + this->fileName, settings->messageShift );
63 
64  //================================================ Allocate the memory
65  this->allocateRRPMemory ( );
66 
67  //================================================ Start computation: For each band (l)
68  proshade_double descValR = 0.0;
69  proshade_unsign arrPos1, arrPos2;
70  for ( proshade_unsign band = 0; band < this->maxShellBand; band++ )
71  {
72  //============================================ For each unique shell couple
73  for ( proshade_unsign shell1 = 0; shell1 < this->noSpheres; shell1++ )
74  {
75  //======================================== Does the band exist for this shell1?
76  if ( !ProSHADE_internal_distances::isBandWithinShell ( band, shell1, this->spheres ) )
77  {
78  for ( proshade_unsign shell2 = 0; shell2 < this->noSpheres; shell2++ )
79  {
80  this->rrpMatrices[band][shell1][shell2] = 0.0;
81  this->rrpMatrices[band][shell2][shell1] = 0.0;
82  }
83  continue;
84  }
85 
86  for ( proshade_unsign shell2 = 0; shell2 < this->noSpheres; shell2++ )
87  {
88  //==================================== Compute each values only once
89  if ( shell1 > shell2 ) { continue; }
90 
91  //==================================== Check if band exists for this shell2?
92  if ( !ProSHADE_internal_distances::isBandWithinShell ( band, shell2, this->spheres ) )
93  {
94  this->rrpMatrices[band][shell1][shell2] = 0.0;
95  this->rrpMatrices[band][shell2][shell1] = 0.0;
96  continue;
97  }
98 
99  //==================================== Initialise
100  descValR = 0.0;
101 
102  //==================================== Sum over order (m)
103  for ( proshade_unsign order = 0; order < static_cast< proshade_unsign > ( ( 2 * band ) + 1 ); order++ )
104  {
105  arrPos1 = static_cast< proshade_unsign > ( seanindex ( static_cast< int > ( order ) - static_cast<int > ( band ),
106  static_cast< int > ( band ), static_cast< int > ( this->spheres[shell1]->getLocalBandwidth() ) ) );
107  arrPos2 = static_cast< proshade_unsign > ( seanindex ( static_cast< int > ( order ) - static_cast< int > ( band ),
108  static_cast< int > ( band ), static_cast< int > ( this->spheres[shell2]->getLocalBandwidth() ) ) );
109  descValR += ProSHADE_internal_maths::complexMultiplicationConjugRealOnly ( &this->sphericalHarmonics[shell1][arrPos1][0],
110  &this->sphericalHarmonics[shell1][arrPos1][1],
111  &this->sphericalHarmonics[shell2][arrPos2][0],
112  &this->sphericalHarmonics[shell2][arrPos2][1] );
113  }
114 
115  //==================================== Save the matrices
116  this->rrpMatrices[band][shell1][shell2] = descValR;
117  this->rrpMatrices[band][shell2][shell1] = descValR;
118  }
119  }
120  }
121 
122  //================================================ Report progress
123  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "RRP matrices successfully computed.", settings->messageShift );
124 
125  //================================================ Done
126  return ;
127 
128 }
129 
139 bool ProSHADE_internal_distances::isBandWithinShell ( proshade_unsign bandInQuestion, proshade_unsign shellInQuestion, ProSHADE_internal_spheres::ProSHADE_sphere** spheres )
140 {
141  if ( bandInQuestion < spheres[shellInQuestion]->getLocalBandwidth() )
142  {
143  return ( true );
144  }
145  else
146  {
147  return ( false );
148  }
149 }
150 
162 {
163  //================================================ Report starting the task
164  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting energy levels distance computation.", settings->messageShift );
165 
166  //================================================ Initialise local variables
167  proshade_double ret = 0.0;
168  std::vector<proshade_double> bandDists;
169 
170  //================================================ Sanity check
171  if ( !settings->computeEnergyLevelsDesc )
172  {
173  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!" );
174  }
175 
176  //================================================ Get the RRP matrices for both objects
177  obj1->computeRRPMatrices ( settings );
178  obj2->computeRRPMatrices ( settings );
179 
180  //================================================ Find the minimium comparable shells and bands
181  proshade_unsign minCommonShells = std::min ( obj1->getMaxSpheres(), obj2->getMaxSpheres() );
182  proshade_unsign minCommonBands = std::min ( obj1->getMaxBand(), obj2->getMaxBand() );
183 
184  //================================================ Get the Pearson's coefficients for each common band
185  computeRRPPearsonCoefficients ( obj1, obj2, settings, minCommonBands, minCommonShells, &bandDists );
186 
187  //================================================ Get distance (by averaging Patterson's coefficients)
188  ret = static_cast<proshade_double> ( std::accumulate ( bandDists.begin(), bandDists.end(), 0.0 ) ) /
189  static_cast<proshade_double> ( bandDists.size() );
190 
191  //================================================ Report completion
192  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Energy levels distance computation complete.", settings->messageShift );
193 
194  //================================================ Done
195  return ( ret );
196 
197 }
198 
211 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 )
212 {
213  //================================================ Report completion
214  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Correlating RRP matrices.", settings->messageShift );
215 
216  //================================================ Initialise local variables
217  proshade_double *str1Vals = new proshade_double[minCommonShells * minCommonShells];
218  proshade_double *str2Vals = new proshade_double[minCommonShells * minCommonShells];
219  ProSHADE_internal_misc::checkMemoryAllocation ( str1Vals, __FILE__, __LINE__, __func__ );
220  ProSHADE_internal_misc::checkMemoryAllocation ( str2Vals, __FILE__, __LINE__, __func__ );
221  proshade_unsign arrIter = 0;
222 
223  //================================================ Start computation: For each band (l)
224  for ( proshade_unsign band = 0; band < minCommonBands; band++ )
225  {
226  //============================================ Reset local counter
227  arrIter = 0;
228 
229  //============================================ For each shell pair
230  for ( proshade_unsign shell1 = 0; shell1 < minCommonShells; shell1++ )
231  {
232  //======================================== Check if band exists (progressive only)
233  if ( settings->progressiveSphereMapping ) { if ( !obj1->shellBandExists( shell1, band ) || !obj2->shellBandExists( shell1, band ) ) { continue; } }
234 
235  for ( proshade_unsign shell2 = 0; shell2 < minCommonShells; shell2++ )
236  {
237  //============================ Check the other shell as well
238  if ( !obj1->shellBandExists( shell2, band ) || !obj2->shellBandExists( shell2, band ) ) { continue; }
239 
240  //==================================== Set values between which the Person's correlation coefficient should be computed
241  str1Vals[arrIter] = obj1->getRRPValue ( band, shell1, shell2 ) *
242  pow ( static_cast<proshade_double> ( shell1 ), settings->enLevMatrixPowerWeight ) *
243  pow ( static_cast<proshade_double> ( shell2 ), settings->enLevMatrixPowerWeight );
244  str2Vals[arrIter] = obj2->getRRPValue ( band, shell1, shell2 ) *
245  pow ( static_cast<proshade_double> ( shell1 ), settings->enLevMatrixPowerWeight ) *
246  pow ( static_cast<proshade_double> ( shell2 ), settings->enLevMatrixPowerWeight );
247 
248  arrIter += 1;
249  }
250  }
251 
252  //============================================ Get Pearson's Correlation Coefficient
253  ProSHADE_internal_misc::addToDoubleVector ( bandDists, ProSHADE_internal_maths::pearsonCorrCoeff ( str1Vals, str2Vals, arrIter ) );
254  }
255 
256  //================================================ Clean up
257  delete[] str1Vals;
258  delete[] str2Vals;
259 
260  //================================================ Report completion
261  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "RRP matrices correlation computed.", settings->messageShift );
262 
263  //================================================ Done
264  return ;
265 }
266 
276 void ProSHADE_internal_data::ProSHADE_data::allocateEMatrices ( proshade_unsign band, proshade_single oversamplingRatio )
277 {
278  //================================================ Compute oversampled band
279  proshade_unsign oversampleEMatricesBy = static_cast< proshade_unsign > ( static_cast< proshade_single > ( band ) * ( 1.0f - oversamplingRatio ) );
280 
281  //================================================ Save the maximum band to the object
282  this->maxEMatDim = band + oversampleEMatricesBy;
283 
284  //================================================ Allocate the required memory
285  this->eMatrices = new proshade_complex** [this->maxEMatDim];
286  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices, __FILE__, __LINE__, __func__ );
287 
288  for ( proshade_unsign bandIter = 0; bandIter < this->maxEMatDim; bandIter++ )
289  {
290  //============================================ Allocate the data structure
291  this->eMatrices[bandIter] = new proshade_complex* [static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 )];
292  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices[bandIter], __FILE__, __LINE__, __func__ );
293 
294  for ( proshade_unsign band2Iter = 0; band2Iter < static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 ); band2Iter++ )
295  {
296  this->eMatrices[bandIter][band2Iter] = new proshade_complex [static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 )];
297  ProSHADE_internal_misc::checkMemoryAllocation ( this->eMatrices[bandIter][band2Iter], __FILE__, __LINE__, __func__ );
298  }
299  }
300 
301  //================================================ Set memory to zero
302  for ( proshade_unsign bandIter = 0; bandIter < this->maxEMatDim; bandIter++ )
303  {
304  for ( proshade_unsign band2Iter = 0; band2Iter < static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 ); band2Iter++ )
305  {
306  for ( proshade_unsign band3Iter = 0; band3Iter < static_cast<proshade_unsign> ( ( bandIter * 2 ) + 1 ); band3Iter++ )
307  {
308  this->eMatrices[bandIter][band2Iter][band3Iter][0] = 0.0;
309  this->eMatrices[bandIter][band2Iter][band3Iter][1] = 0.0;
310  }
311  }
312  }
313 
314  //================================================ Done
315  return ;
316 
317 }
318 
329 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 )
330 {
331  //================================================ Allocate the memory
332  obj1Vals = new proshade_double [minSpheres];
333  obj2Vals = new proshade_double [minSpheres];
334  radiiVals = new proshade_complex[minSpheres];
335  GLabscissas = new proshade_double [intOrder];
336  GLweights = new proshade_double [intOrder];
337 
338  //================================================ Check the memory allocation
339  ProSHADE_internal_misc::checkMemoryAllocation ( obj1Vals, __FILE__, __LINE__, __func__ );
340  ProSHADE_internal_misc::checkMemoryAllocation ( obj2Vals, __FILE__, __LINE__, __func__ );
341  ProSHADE_internal_misc::checkMemoryAllocation ( radiiVals, __FILE__, __LINE__, __func__ );
342  ProSHADE_internal_misc::checkMemoryAllocation ( GLabscissas, __FILE__, __LINE__, __func__ );
343  ProSHADE_internal_misc::checkMemoryAllocation ( GLweights, __FILE__, __LINE__, __func__ );
344 
345  //================================================ Done
346  return ;
347 
348 }
349 
358 void ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude ( ProSHADE_internal_data::ProSHADE_data* obj, int band, int order, proshade_unsign radius, proshade_double* result )
359 {
360  //================================================ Pre-compute values
361  int locBand = static_cast< int > ( obj->spheres[radius]->getLocalBandwidth() );
362  int objArrPos = seanindex ( order - band, band, locBand );
363 
364  //================================================ Find the magnitude
366  &obj->sphericalHarmonics[radius][objArrPos][1],
367  &obj->sphericalHarmonics[radius][objArrPos][0],
368  &obj->sphericalHarmonics[radius][objArrPos][1] );
369 
370  //================================================ Weight by radius^2 for the integration that will follow
371  *result *= pow ( static_cast<proshade_double> ( obj->getAnySphereRadius( radius ) ), 2.0 );
372 
373  //================================================ Done
374  return ;
375 
376 }
377 
391 void ProSHADE_internal_distances::computeEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, int bandIter, int orderIter, proshade_complex* radiiVals, int integOrder, proshade_double* abscissas, proshade_double* weights, proshade_double integRange, proshade_double sphereDist )
392 {
393  //================================================ Initialise local variables
394  proshade_unsign objCombValsIter = 0;
395  proshade_double hlpReal, hlpImag, rSquared;
396  proshade_complex arrVal;
397  int o1ArrPos, o2ArrPos, locBand;
398  proshade_unsign integOrderU = static_cast< proshade_unsign > ( integOrder );
399 
400  //================================================ For each combination of m and m' for E matrices
401  for ( int order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
402  {
403  //============================================ Reset loop
404  objCombValsIter = 0;
405 
406  //============================================ Find the c*conj(c) values for different radii
407  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
408  {
409  //======================================== Get only values where the shell has the band
410  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= static_cast< proshade_unsign > ( bandIter ) ) { continue; }
411 
412  //======================================== Pre-compute values
413  rSquared = pow ( ( static_cast<proshade_double> ( obj1->getAnySphereRadius( radiusIter ) ) ), 2.0 );
414 
415  //======================================== Multiply coeffs
416  locBand = static_cast< int > ( obj1->spheres[radiusIter]->getLocalBandwidth() );
417  o1ArrPos = seanindex ( orderIter - bandIter, bandIter, locBand );
418  o2ArrPos = seanindex ( order2Iter - bandIter, bandIter, locBand );
419 
421  &obj1->sphericalHarmonics[radiusIter][o1ArrPos][1],
422  &obj2->sphericalHarmonics[radiusIter][o2ArrPos][0],
423  &obj2->sphericalHarmonics[radiusIter][o2ArrPos][1],
424  &hlpReal, &hlpImag );
425 
426  //======================================== Apply r^2 integral weight
427  radiiVals[objCombValsIter][0] = hlpReal * rSquared;
428  radiiVals[objCombValsIter][1] = hlpImag * rSquared;
429 
430  objCombValsIter += 1;
431  }
432 
433  //============================================ Integrate over all radii using n-point Gauss-Legendre integration
434  ProSHADE_internal_maths::gaussLegendreIntegration ( radiiVals, objCombValsIter, integOrderU, abscissas, weights, integRange, sphereDist, &hlpReal, &hlpImag );
435 
436  //============================================ Save the result into E matrices
437  arrVal[0] = hlpReal;
438  arrVal[1] = hlpImag;
439  obj2->setEMatrixValue ( bandIter, orderIter, order2Iter, arrVal );
440  }
441 
442  //================================================ Done
443  return ;
444 
445 }
446 
462 proshade_double ProSHADE_internal_distances::computeWeightsForEMatricesForLM ( ProSHADE_internal_data::ProSHADE_data* obj1, ProSHADE_internal_data::ProSHADE_data* obj2, int bandIter, int orderIter, proshade_double* obj1Vals, proshade_double* obj2Vals, int integOrder, proshade_double* abscissas, proshade_double* weights, proshade_single sphereDist )
463 {
464  //================================================ Initialise local values
465  proshade_unsign obj1ValsIter = 0;
466  proshade_unsign obj2ValsIter = 0;
467  proshade_unsign integOrderU = static_cast< proshade_unsign > ( integOrder );
468 
469  //================================================ Set sphere counters
470  proshade_unsign minSphere = std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() );
471  proshade_unsign maxSphere = 0;
472 
473  //================================================ For each radius, deal with weights
474  for ( proshade_unsign radiusIter = 0; radiusIter < std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ); radiusIter++ )
475  {
476  //============================================ Get only values where the shell has the band
477  if ( std::min ( obj1->getShellBandwidth ( radiusIter ), obj2->getShellBandwidth ( radiusIter ) ) <= static_cast< proshade_unsign > ( bandIter ) ) { continue; }
478  minSphere = std::min ( radiusIter, minSphere );
479  maxSphere = std::max ( radiusIter, maxSphere );
480 
481  //============================================ Get the magnitudes for weighting
482  computeSphericalHarmonicsMagnitude ( obj1, bandIter, orderIter, radiusIter, &(obj1Vals[obj1ValsIter]) );
483  computeSphericalHarmonicsMagnitude ( obj2, bandIter, orderIter, radiusIter, &(obj2Vals[obj2ValsIter]) );
484  obj1ValsIter += 1;
485  obj2ValsIter += 1;
486  }
487 
488  //================================================ Integrate weights
489  proshade_single minSphereRad = obj1->getSpherePosValue ( minSphere ) - ( sphereDist * 0.5f );
490  proshade_single maxSphereRad = obj1->getSpherePosValue ( maxSphere ) + ( sphereDist * 0.5f );
491 
492  obj1->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj1Vals, obj1ValsIter, integOrderU, abscissas, weights, static_cast< proshade_double > ( maxSphereRad - minSphereRad ), static_cast< proshade_double > ( sphereDist ) ) );
493  obj2->setIntegrationWeightCumul ( ProSHADE_internal_maths::gaussLegendreIntegrationReal ( obj2Vals, obj2ValsIter, integOrderU, abscissas, weights, static_cast< proshade_double > ( maxSphereRad - minSphereRad ), static_cast< proshade_double > ( sphereDist ) ) );
494 
495  //================================================ Done
496  return ( static_cast< proshade_double > ( maxSphereRad - minSphereRad ) );
497 
498 }
499 
508 void ProSHADE_internal_distances::releaseTrSigmaWorkspace ( proshade_double*& obj1Vals, proshade_double*& obj2Vals, proshade_double*& GLabscissas, proshade_double*& GLweights, proshade_complex*& radiiVals )
509 {
510  //================================================ Release memory
511  delete[] obj1Vals;
512  delete[] obj2Vals;
513  delete[] radiiVals;
514  delete[] GLabscissas;
515  delete[] GLweights;
516 
517  //================================================ Set to NULL
518  obj1Vals = nullptr;
519  obj2Vals = nullptr;
520  radiiVals = nullptr;
521  GLabscissas = nullptr;
522  GLweights = nullptr;
523 
524  //================================================ Done
525  return ;
526 
527 }
528 
541 {
542  //================================================ Report progress
543  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Starting computation of E matrices.", settings->messageShift );
544 
545  //================================================ Allocatre memory for E matrices in the second object (first may be compared to more structures and therefore its data would be written over)
546  obj2->allocateEMatrices ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ), settings->resolutionOversampling );
547 
548  //================================================ Initialise local variables
549  proshade_double *obj1Vals, *obj2Vals, *GLAbscissas, *GLWeights;
550  proshade_complex* radiiVals;
551  proshade_double integRange;
552 
553  //================================================ Allocate workspace memory
554  allocateTrSigmaWorkspace ( std::min( obj1->getMaxSpheres(), obj2->getMaxSpheres() ), settings->integOrder, obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals);
555 
556  //================================================ For each band (l), compute the E matrix integrals
557  for ( int bandIter = 0; bandIter < static_cast< int > ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) ); bandIter++ )
558  {
559  //============================================ Initialise abscissas and weights for integration
560  int localIntegOrder = std::min ( static_cast< int > ( std::max ( bandIter, 5 ) ), static_cast< int > ( settings->integOrder ) );
561  if ( settings->noIntegrationSpeedup ) { localIntegOrder = static_cast< int > ( settings->integOrder ); }
562  ProSHADE_internal_maths::getLegendreAbscAndWeights ( static_cast< proshade_unsign > ( localIntegOrder ), GLAbscissas, GLWeights, settings->integApproxSteps );
563 
564  //============================================ For each order (m)
565  for ( int orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
566  {
567  //======================================== Get weights for the required band(l) and order (m)
568  integRange = computeWeightsForEMatricesForLM ( obj1, obj2, bandIter, orderIter, obj1Vals, obj2Vals, localIntegOrder, GLAbscissas, GLWeights, settings->maxSphereDists );
569 
570  //======================================== Compute E matrices value for given band (l) and order(m)
571  computeEMatricesForLM ( obj1, obj2, bandIter, orderIter, radiiVals, localIntegOrder, GLAbscissas, GLWeights, integRange, static_cast< proshade_double > ( settings->maxSphereDists ) );
572  }
573 
574  //============================================ Report progress
575  if ( settings->verbose > 3 )
576  {
577  std::stringstream hlpSS;
578  hlpSS << "E matrices computed for band " << bandIter;
579  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, hlpSS.str(), settings->messageShift );
580  }
581  }
582 
583  //================================================ Release the workspace memory
584  releaseTrSigmaWorkspace ( obj1Vals, obj2Vals, GLAbscissas, GLWeights, radiiVals );
585 
586  //================================================ Report progress
587  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices computed.", settings->messageShift );
588 
589  //================================================ Done
590  return ;
591 
592 }
593 
605 {
606  //================================================ Report progress
607  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Starting E matrices normalisation.", settings->messageShift );
608 
609  //================================================ Normalise by the Pearson's c.c. like formula
610  proshade_double eMatNormFactor = std::sqrt ( obj1->getIntegrationWeight() * obj2->getIntegrationWeight() );
611 
612  //================================================ If this is self-correlation (i.e. obj1 == obj2), then divide normalisation factor by 2 as the weight was applied cumulatively!
613  if ( obj1->inputOrder == obj2->inputOrder ) { eMatNormFactor /= 2.0; }
614 
615  for ( proshade_unsign bandIter = 0; bandIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); bandIter++ )
616  {
617  //============================================ For each combination of m and m' for E matrices
618  for ( proshade_unsign orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
619  {
620  for ( proshade_unsign order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
621  {
622  obj2->normaliseEMatrixValue ( bandIter, orderIter, order2Iter, eMatNormFactor );
623  }
624  }
625  }
626 
627  //================================================ Report progress
628  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 4, "E matrices normalised.", settings->messageShift );
629 
630  //================================================ Done
631  return ;
632 
633 }
634 
649 {
650  //================================================ Report starting the task
651  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting trace sigma distance computation.", settings->messageShift );
652 
653  //================================================ Initialise return variable
654  proshade_double ret = 0.0;
655 
656  //================================================ Sanity check
657  if ( !settings->computeTraceSigmaDesc )
658  {
659  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!" );
660  }
661 
662  //================================================ Empty the cumulative weights back to 0.0 for each structure
663  obj1->setIntegrationWeight ( 0.0 );
664  obj1->setIntegrationWeight ( 0.0 );
665 
666  //================================================ Compute un-weighted E matrices and their weights
667  computeEMatrices ( obj1, obj2, settings );
668 
669  //================================================ Normalise E matrices by the magnitudes
670  normaliseEMatrices ( obj1, obj2, settings );
671 
672  //================================================ Allocate the required memory
673  double* singularValues = new double[( ( std::min ( obj1->getMaxBand(), obj2->getMaxBand() ) * 2 ) + 1 )];
674  ProSHADE_internal_misc::checkMemoryAllocation ( singularValues, __FILE__, __LINE__, __func__ );
675 
676  //================================================ Compute the distance
677  for ( proshade_unsign lIter = 0; lIter < std::min ( obj1->getMaxBand(), obj2->getMaxBand() ); lIter++ )
678  {
679  //============================================ Find the complex matrix SVD singular values
680  ProSHADE_internal_maths::complexMatrixSVDSigmasOnly ( obj2->getEMatrixByBand ( lIter ), static_cast<int> ( ( lIter * 2 ) + 1 ), singularValues );
681 
682  //============================================ Now sum the trace
683  for ( proshade_unsign iter = 0; iter < ( ( lIter * 2 ) + 1 ); iter++ )
684  {
685  ret += singularValues[iter];
686  }
687  }
688 
689  //================================================ Report completion
690  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "E matrices decomposed to singular values.", settings->messageShift );
691 
692  //================================================ Release the memory
693  delete[] singularValues;
694 
695  //================================================ Report completion
696  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Trace sigma distance computation complete.", settings->messageShift );
697 
698  //================================================ Done
699  return ( ret );
700 
701 }
702 
708 {
709  //================================================ Allocate the memory
710  this->so3Coeffs = reinterpret_cast< fftw_complex* > ( fftw_malloc ( sizeof ( fftw_complex ) * static_cast<proshade_unsign>( ( 4 * pow( static_cast<proshade_double> ( band ), 3.0 ) - static_cast<proshade_double> ( band ) ) / 3.0 ) ) );
711  this->so3CoeffsInverse = reinterpret_cast< fftw_complex* > ( fftw_malloc ( sizeof ( fftw_complex ) * static_cast<proshade_unsign>( pow( static_cast<proshade_double> ( band ) * 2.0, 3.0 ) ) ) );
712 
713  //================================================ Check memory allocation
714  ProSHADE_internal_misc::checkMemoryAllocation ( this->so3Coeffs, __FILE__, __LINE__, __func__ );
715  ProSHADE_internal_misc::checkMemoryAllocation ( this->so3CoeffsInverse, __FILE__, __LINE__, __func__ );
716 
717  //================================================ Done
718  return ;
719 
720 }
721 
733 {
734  //================================================ Report progress
735  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Converting E matrices to SO(3) coefficients.", settings->messageShift );
736 
737  //================================================ Allocate memory for the coefficients
738  obj2->allocateSO3CoeffsSpace ( obj2->getEMatDim ( ) );
739 
740  //================================================ Initialise local variables
741  proshade_double wigNorm, hlpValReal, hlpValImag;
742  proshade_double signValue = 1.0;
743  proshade_unsign indexO;
744  proshade_complex hlpVal;
745 
746  //================================================ For each band (l)
747  for ( proshade_signed bandIter = 0; bandIter < static_cast<proshade_signed> ( obj2->getEMatDim ( ) ); bandIter++ )
748  {
749  //============================================ Get wigner normalisation factor
750  wigNorm = 2.0 * M_PI * sqrt ( 2.0 / (2.0 * static_cast< proshade_double > ( bandIter ) + 1.0 ) );
751 
752  //============================================ For each order (m)
753  for ( proshade_signed orderIter = 0; orderIter < ( ( bandIter * 2 ) + 1 ); orderIter++ )
754  {
755  //======================================== Set the sign
756  if ( ( orderIter - bandIter + bandIter ) % 2 ) { signValue = -1.0 ; }
757  else { signValue = 1.0 ; }
758 
759  //======================================== For each order2 (m')
760  for ( proshade_signed order2Iter = 0; order2Iter < ( ( bandIter * 2 ) + 1 ); order2Iter++ )
761  {
762  //==================================== Find output index
763  indexO = static_cast< proshade_unsign > ( so3CoefLoc ( static_cast< int > ( orderIter - bandIter ), static_cast< int > ( order2Iter - bandIter ), static_cast< int > ( bandIter ), static_cast< int > ( obj2->getEMatDim ( ) ) ) );
764 
765  //==================================== Compute and save the SO(3) coefficients
766  obj2->getEMatrixValue ( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( orderIter ), static_cast< proshade_unsign > ( order2Iter ), &hlpValReal, &hlpValImag );
767  hlpVal[0] = hlpValReal * wigNorm * signValue;
768  hlpVal[1] = hlpValImag * wigNorm * signValue;
769  obj2->setSO3CoeffValue ( indexO, hlpVal );
770 
771  //==================================== Switch the sign value
772  signValue *= -1.0;
773  }
774  }
775  }
776 
777  //================================================ Report progress
778  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "SO(3) coefficients obtained.", settings->messageShift );
779 
780  //================================================ Done
781  return ;
782 
783 }
784 
792 void ProSHADE_internal_distances::allocateInvSOFTWorkspaces ( proshade_complex*& work1, proshade_complex*& work2, proshade_double*& work3, proshade_unsign band )
793 {
794  //================================================ Allocate memory
795  work1 = new proshade_complex[8 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 3.0 ) )];
796  work2 = new proshade_complex[14 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (48 * band)];
797  work3 = new proshade_double [2 * static_cast<proshade_unsign> ( pow( static_cast<double> ( band ), 2.0 ) ) + (24 * band)];
798 
799  //================================================ Check the memory allocation
800  ProSHADE_internal_misc::checkMemoryAllocation ( work1, __FILE__, __LINE__, __func__ );
801  ProSHADE_internal_misc::checkMemoryAllocation ( work2, __FILE__, __LINE__, __func__ );
802  ProSHADE_internal_misc::checkMemoryAllocation ( work3, __FILE__, __LINE__, __func__ );
803 
804  //================================================ Done
805  return ;
806 
807 }
808 
816 void ProSHADE_internal_distances::prepareInvSOFTPlan ( fftw_plan* inverseSO3, int band, fftw_complex* work1, proshade_complex* invCoeffs )
817 {
818  //================================================ Prepare the plan describing variables
819  int howmany = 4 * band * band;
820  int idist = 2 * band;
821  int odist = 2 * band;
822  int rank = 2;
823 
824  int inembed[2], onembed[2];
825  inembed[0] = 2 * band;
826  inembed[1] = 4 * band * band;
827  onembed[0] = 2 * band;
828  onembed[1] = 4 * band * band;
829 
830  int istride = 1;
831  int ostride = 1;
832 
833  int na[2];
834  na[0] = 1;
835  na[1] = 2 * band;
836 
837  //================================================ Create the plan
838  *inverseSO3 = fftw_plan_many_dft ( rank,
839  na,
840  howmany,
841  work1,
842  inembed,
843  istride,
844  idist,
845  invCoeffs,
846  onembed,
847  ostride,
848  odist,
849  FFTW_FORWARD,
850  FFTW_ESTIMATE );
851 
852  //================================================ Done
853  return ;
854 
855 }
856 
863 void ProSHADE_internal_distances::releaseInvSOFTMemory ( proshade_complex*& work1, proshade_complex*& work2, proshade_double*& work3 )
864 {
865  //================================================ Release memory
866  delete[] work1;
867  delete[] work2;
868  delete[] work3;
869 
870  //================================================ Done
871  return ;
872 }
873 
886 {
887  //================================================ Report progress
888  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Computing inverse SO(3) Fourier transform.", settings->messageShift );
889 
890  //================================================ Initialise local variables
891  proshade_complex *workspace1, *workspace2;
892  proshade_double *workspace3;
893  fftw_plan inverseSO3;
894 
895  //================================================ Allocate memory for the workspaces
896  allocateInvSOFTWorkspaces ( workspace1, workspace2, workspace3, obj2->getEMatDim ( ) );
897 
898  //================================================ Prepare the FFTW plan
899  prepareInvSOFTPlan ( &inverseSO3, static_cast< int > ( obj2->getEMatDim ( ) ), workspace1, obj2->getInvSO3Coeffs ( ) );
900 
901  //================================================ Compute the transform
902  Inverse_SO3_Naive_fftw ( static_cast< int > ( obj2->getEMatDim ( ) ),
903  obj2->getSO3Coeffs ( ),
904  obj2->getInvSO3Coeffs ( ),
905  workspace1,
906  workspace2,
907  workspace3,
908  &inverseSO3,
909  0 );
910 
911  //================================================ Release memory
912  releaseInvSOFTMemory ( workspace1, workspace2, workspace3 );
913  fftw_destroy_plan ( inverseSO3 );
914 
915  //================================================ Report progress
916  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 3, "Inverse SO(3) Fourier transform computed.", settings->messageShift );
917 
918 
919 
920  //================================================ Done
921  return ;
922 
923 }
924 
940 {
941  //================================================ Report starting the task
942  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 1, "Starting rotation function distance computation.", settings->messageShift );
943 
944  //================================================ Initialise return variable
945  proshade_double ret = 0.0;
946  proshade_double eulA, eulB, eulG, EMatR, EMatI, WigDR, WigDI;
947 
948  //================================================ Sanity check
949  if ( !settings->computeRotationFuncDesc )
950  {
951  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!" );
952  }
953 
954  //================================================ Compute weighted E matrices if not already present
955  if ( !settings->computeTraceSigmaDesc )
956  {
957  computeEMatrices ( obj1, obj2, settings );
958  normaliseEMatrices ( obj1, obj2, settings );
959  }
960 
961  //================================================ Generate SO(3) coefficients
962  generateSO3CoeffsFromEMatrices ( obj2, settings );
963 
964  //================================================ Compute the inverse SO(3) Fourier Transform (SOFT) on the newly computed coefficients
965  computeInverseSOFTTransform ( obj2, settings );
966 
967  //================================================ Get inverse SO(3) map top peak Euler angle values
969  obj2->getEMatDim ( ) * 2,
970  &eulA, &eulB, &eulG, settings );
971 
972  //================================================ Compute the Wigner D matrices for the Euler angles
973  ProSHADE_internal_wigner::computeWignerMatricesForRotation ( settings, obj2, eulA, eulB, eulG );
974 
975  //================================================ Compute the distance
976  for ( proshade_unsign bandIter = 0; bandIter < obj2->getEMatDim(); bandIter++ )
977  {
978  //============================================ For each order1
979  for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
980  {
981  //======================================== For each order2
982  for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
983  {
984  //==================================== Multiply D_{l} * E_{l} and get sum over l of traces (i.e. just sum all together)
985  obj2->getEMatrixValue ( bandIter, order1, order2, &EMatR, &EMatI );
986  obj2->getWignerMatrixValue ( bandIter, order2, order1, &WigDR, &WigDI );
987  ret += ProSHADE_internal_maths::complexMultiplicationRealOnly ( &WigDR, &WigDI, &EMatR, &EMatI );
988  }
989  }
990  }
991 
992  //================================================ Report completion
993  ProSHADE_internal_messages::printProgressMessage ( settings->verbose, 2, "Rotation function distance computation complete.", settings->messageShift );
994 
995  //================================================ Done
996  return ( ret );
997 
998 }
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:604
ProSHADE_settings::integOrder
proshade_unsign integOrder
The order required for full Gauss-Legendre integration between the spheres.
Definition: ProSHADE_settings.hpp:70
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:3962
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:714
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:119
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:624
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:120
ProSHADE_internal_spheres::ProSHADE_sphere::getLocalBandwidth
proshade_unsign getLocalBandwidth(void)
This function returns the local bandwidth.
Definition: ProSHADE_spheres.cpp:391
ProSHADE_internal_data::ProSHADE_data::sphericalHarmonics
proshade_complex ** sphericalHarmonics
A set of spherical harmonics values arrays for each sphere.
Definition: ProSHADE_data.hpp:121
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:329
ProSHADE_internal_distances::computeEMatricesForLM
void computeEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, int bandIter, int orderIter, proshade_complex *radiiVals, int 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:391
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:4305
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:4354
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
ProSHADE_internal_data::ProSHADE_data::allocateRRPMemory
void allocateRRPMemory()
This function allocates the required memory for the RRP matrices.
Definition: ProSHADE_distances.cpp:31
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:3977
ProSHADE_internal_distances::computeRotationFunctionDescriptor
proshade_double computeRotationFunctionDescriptor(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:939
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::prepareInvSOFTPlan
void prepareInvSOFTPlan(fftw_plan *inverseSO3, int band, fftw_complex *work1, proshade_complex *invCoeffs)
This function prepares the FFTW plan for the inverse SO(3) transform.
Definition: ProSHADE_distances.cpp:816
ProSHADE_internal_data::ProSHADE_data::inputOrder
proshade_unsign inputOrder
This value is the input order - it is useful to know for writing out files, so that they would not ov...
Definition: ProSHADE_data.hpp:140
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:863
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:3727
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:4370
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:139
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:59
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:3915
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:111
ProSHADE_internal_data::ProSHADE_data::getSO3Coeffs
proshade_complex * getSO3Coeffs(void)
This function allows access to the SO(3) coefficients array.
Definition: ProSHADE_data.cpp:4003
ProSHADE_settings::maxSphereDists
proshade_single maxSphereDists
The distance between spheres in spherical mapping for the largest sphere.
Definition: ProSHADE_settings.hpp:67
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::getIntegrationWeight
proshade_double getIntegrationWeight(void)
This function allows access to the integration weight for the object.
Definition: ProSHADE_data.cpp:3926
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:4319
ProSHADE_internal_spheres::ProSHADE_sphere
This class contains all inputed and derived data for a single sphere.
Definition: ProSHADE_spheres.hpp:49
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:152
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:807
ProSHADE_settings::messageShift
proshade_signed messageShift
This value allows shifting the messages to create more readable log for sub-processes.
Definition: ProSHADE_settings.hpp:153
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:4029
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:3938
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:211
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:246
ProSHADE_internal_distances::generateSO3CoeffsFromEMatrices
void generateSO3CoeffsFromEMatrices(ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function converts the E matrices to SO(3) coefficients.
Definition: ProSHADE_distances.cpp:732
ProSHADE_internal_wigner::computeWignerMatricesForRotation
void computeWignerMatricesForRotation(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *obj, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function computes the Wigner D matrices for a particular set of Euler angles.
Definition: ProSHADE_wignerMatrices.cpp:258
ProSHADE_internal_data::ProSHADE_data::rrpMatrices
proshade_double *** rrpMatrices
The energy levels descriptor shell correlation tables.
Definition: ProSHADE_data.hpp:127
ProSHADE_internal_peakSearch::getBestPeakEulerAngsNaive
void getBestPeakEulerAngsNaive(proshade_complex *map, proshade_unsign dim, proshade_double *eulA, proshade_double *eulB, proshade_double *eulG, ProSHADE_settings *settings)
This function finds the highest peaks optimised Euler angles using the "naive" approach.
Definition: ProSHADE_peakSearch.cpp:352
ProSHADE_internal_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:3774
ProSHADE_internal_data::ProSHADE_data::getMaxBand
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
Definition: ProSHADE_data.cpp:3748
ProSHADE_settings::resolutionOversampling
proshade_single resolutionOversampling
How much (%) should the requested resolution be over-sampled by map re-sampling?
Definition: ProSHADE_settings.hpp:53
ProSHADE_settings
This class stores all the settings and is passed to the executive classes instead of a multitude of p...
Definition: ProSHADE_settings.hpp:37
ProSHADE_settings::integApproxSteps
proshade_unsign integApproxSteps
The number of steps taken in the approximation of Legendre polynomial decomposition to terms.
Definition: ProSHADE_settings.hpp:71
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:540
ProSHADE_internal_data::ProSHADE_data::getSpherePosValue
proshade_single getSpherePosValue(proshade_unsign shell)
This function allows access to sphere positions.
Definition: ProSHADE_data.cpp:3950
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:3758
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:648
ProSHADE_internal_distances::computeInverseSOFTTransform
void computeInverseSOFTTransform(ProSHADE_internal_data::ProSHADE_data *obj2, ProSHADE_settings *settings)
This function computes the inverse SO(3) transform.
Definition: ProSHADE_distances.cpp:885
ProSHADE_internal_distances::computeWeightsForEMatricesForLM
proshade_double computeWeightsForEMatricesForLM(ProSHADE_internal_data::ProSHADE_data *obj1, ProSHADE_internal_data::ProSHADE_data *obj2, int bandIter, int orderIter, proshade_double *obj1Vals, proshade_double *obj2Vals, int integOrder, proshade_double *abscissas, proshade_double *weights, proshade_single sphereDist)
This function computes the E matrix weight values for a given band and order and saves these into the...
Definition: ProSHADE_distances.cpp:462
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:73
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:508
ProSHADE_internal_distances::computeSphericalHarmonicsMagnitude
void computeSphericalHarmonicsMagnitude(ProSHADE_internal_data::ProSHADE_data *obj, int band, int 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:358
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:792
ProSHADE_internal_data::ProSHADE_data::setEMatrixValue
void setEMatrixValue(int band, int order1, int order2, proshade_complex val)
This function allows setting the E matrix value.
Definition: ProSHADE_data.cpp:4336
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:3992
ProSHADE_internal_data::ProSHADE_data::getEMatDim
proshade_unsign getEMatDim(void)
This function allows access to the maximum band for the E matrix.
Definition: ProSHADE_data.cpp:4014
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::allocateEMatrices
void allocateEMatrices(proshade_unsign band, proshade_single oversamplingRatio)
This function allocates the required memory for the E matrices.
Definition: ProSHADE_distances.cpp:276
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:161
ProSHADE_internal_maths::getLegendreAbscAndWeights
void getLegendreAbscAndWeights(proshade_unsign order, proshade_double *abscissas, proshade_double *weights, proshade_unsign noSteps)
Function to prepare abscissas and weights for Gauss-Legendre integration using the Glaser-Liu-Rokhlin...
Definition: ProSHADE_maths.cpp:289
ProSHADE_settings::noIntegrationSpeedup
bool noIntegrationSpeedup
This option turns off the integration speedup (only using abscissas and weights up to appropriate l f...
Definition: ProSHADE_settings.hpp:72
ProSHADE_internal_data::ProSHADE_data::spheres
ProSHADE_internal_spheres::ProSHADE_sphere ** spheres
The set of concentric spheres to which the intermal density map has been projected.
Definition: ProSHADE_data.hpp:120
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:707
ProSHADE_internal_messages::printProgressMessage
void printProgressMessage(proshade_signed verbose, proshade_signed messageLevel, std::string message, proshade_signed messageShift=0)
General stdout message printing.
Definition: ProSHADE_messages.cpp:71
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:118
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:117