ProSHADE  0.7.6.6 (JUL 2022)
Protein Shape Detection
ProSHADE_maths.hpp
Go to the documentation of this file.
1 
23 //==================================================== ProSHADE
24 #include "ProSHADE_misc.hpp"
25 
26 //==================================================== almostEqual for comparing floating point numbers (BSD License, works on Windows as well as linux)
27 #include <almostEqual.hpp>
28 
29 //==================================================== Overinclusion protection
30 #ifndef PROSHADE_MATHS
31 #define PROSHADE_MATHS
32 
33 //==================================================== Declare LAPACK functions
34 extern "C"
35 {
36  // ... The real matrix singular value decomposition function from LAPACK
37  extern void dgesdd_ ( char* jobz, int* m, int* n, double* a, int* lda, double* s, double* u, int* ldu, double* vt, int* ldvt, double* work, int* lwork, double* rwork, int* iwork, int* info );
38  // ... The complex matrix singular value decomposition function from LAPACK
39  extern void zgesdd_ ( char* jobz, int* m, int* n, std::complex<double>* a, int* lda, double* s, std::complex<double>* u, int* ldu, std::complex<double>* vt, int* ldvt, std::complex<double>* work, int* lwork, double* rwork, int* iwork, int* info );
40  // ... The eigenvalue/eigenvector solver
41  extern void dgeev_ ( char* jobvl, char* jobvr, int* n, double* a, int* lda, double* wr, double* wi, double* vl, int* ldvl, double* vr, int* ldvr, double* work, int* lwork, int* info );
42 }
43 
44 //==================================================== ProSHADE_internal_spheres Namespace
45 // ... NOTE: The full namespace is defined in ProSHADE_spheres.h ; however, as forward declaration
46 // ... would not work here, the ProSHADE_rotFun_sphere class needs to be defined here, making the
47 // ... organisation a bit messier.
49 {
58  {
59  private:
60  proshade_double radius;
61  proshade_double radiusMax;
62  proshade_double radiusMin;
63  proshade_unsign angularDim;
64  proshade_unsign rotFunDim;
65  proshade_double representedAngle;
66  proshade_unsign sphereNumber;
67 
68  proshade_double* axesValues;
69  std::vector< std::pair< proshade_unsign,proshade_unsign > > peaks;
70  public:
71  ProSHADE_rotFun_sphere ( proshade_double rad, proshade_double radRange, proshade_unsign dim, proshade_unsign rfDim, proshade_double repAng, proshade_unsign sphNo );
72  ~ProSHADE_rotFun_sphere ( void );
73 
74  public:
75  proshade_double getRadius ( void );
76  proshade_double getMaxRadius ( void );
77  proshade_double getMinRadius ( void );
78  proshade_unsign getAngularDim ( void );
79  proshade_unsign getRotFunDim ( void );
80  proshade_double getRepresentedAngle ( void );
81  proshade_unsign getSphereNumber ( void );
82  std::vector<std::pair<proshade_unsign,proshade_unsign>> getPeaks ( void );
83  proshade_double getSphereLatLonPosition ( proshade_unsign lattitude, proshade_unsign longitude );
84  proshade_double getSphereLatLonLinearInterpolationPos ( proshade_double lattitude, proshade_double longitude );
85  std::vector< std::vector< proshade_double > > getCopyOfValues ( void );
86 
87 
88  public:
89  void interpolateSphereValues ( proshade_complex* rotFun );
90  void findAllPeaks ( proshade_signed noSmNeighbours, std::vector< proshade_double >* allHeights );
91  void removeSmallPeaks ( proshade_double peakThres );
92  };
93 }
94 
95 //==================================================== ProSHADE_internal_io Namespace
102 namespace ProSHADE_internal_maths
103 {
112  {
113  private:
114  proshade_double a00, a01, a02, a03, a10, a11, a12, a13, a20, a21, a22, a23, a30, a31, a32, a33;
115  proshade_double xStartIndex, xRange, yStartIndex, yRange;
116 
117  public:
125  BicubicInterpolator ( proshade_double** areaToInterpolate, proshade_double xStart, proshade_double yStart )
126  {
127  //======================================== Save the original non-unit square positions
128  this->xStartIndex = xStart;
129  this->yStartIndex = yStart;
130 
131  //======================================== Prepare variables for converting from original indices to unit square
132  this->xRange = 1.0;
133  this->yRange = 1.0;
134 
135  //======================================== Precompute interpolators
136  this->a00 = areaToInterpolate[1][1];
137  this->a01 = - ( 0.5 * areaToInterpolate[1][0] ) +
138  ( 0.5 * areaToInterpolate[1][2] );
139  this->a02 = areaToInterpolate[1][0] -
140  ( 2.5 * areaToInterpolate[1][1] ) +
141  ( 2.0 * areaToInterpolate[1][2] ) -
142  ( 0.5 * areaToInterpolate[1][3] );
143  this->a03 = - ( 0.5 * areaToInterpolate[1][0] ) +
144  ( 1.5 * areaToInterpolate[1][1] ) -
145  ( 1.5 * areaToInterpolate[1][2] ) +
146  ( 0.5 * areaToInterpolate[1][3] );
147  this->a10 = - ( 0.5 * areaToInterpolate[0][1] ) +
148  ( 0.5 * areaToInterpolate[2][1] );
149  this->a11 = ( 0.25 * areaToInterpolate[0][0] ) -
150  ( 0.25 * areaToInterpolate[0][2] ) -
151  ( 0.25 * areaToInterpolate[2][0] ) +
152  ( 0.25 * areaToInterpolate[2][2] );
153  this->a12 = - ( 0.5 * areaToInterpolate[0][0] ) +
154  ( 1.25 * areaToInterpolate[0][1] ) -
155  areaToInterpolate[0][2] +
156  ( 0.25 * areaToInterpolate[0][3] ) +
157  ( 0.5 * areaToInterpolate[2][0] ) -
158  ( 1.25 * areaToInterpolate[2][1] ) +
159  areaToInterpolate[2][2] -
160  ( 0.25 * areaToInterpolate[2][3] );
161  this->a13 = ( 0.25 * areaToInterpolate[0][0] ) -
162  ( 0.75 * areaToInterpolate[0][1] ) +
163  ( 0.75 * areaToInterpolate[0][2] ) -
164  ( 0.25 * areaToInterpolate[0][3] ) -
165  ( 0.25 * areaToInterpolate[2][0] ) +
166  ( 0.75 * areaToInterpolate[2][1] ) -
167  ( 0.75 * areaToInterpolate[2][2] ) +
168  ( 0.25 * areaToInterpolate[2][3] );
169  this->a20 = areaToInterpolate[0][1] -
170  ( 2.5 * areaToInterpolate[1][1] ) +
171  ( 2.0 * areaToInterpolate[2][1] ) -
172  ( 0.5 * areaToInterpolate[3][1] );
173  this->a21 = - ( 0.5 * areaToInterpolate[0][0] ) +
174  ( 0.5 * areaToInterpolate[0][2] ) +
175  ( 1.25 * areaToInterpolate[1][0] ) -
176  ( 1.25 * areaToInterpolate[1][2] ) -
177  areaToInterpolate[2][0] + areaToInterpolate[2][2] +
178  ( 0.25 * areaToInterpolate[3][0] ) -
179  ( 0.25 * areaToInterpolate[3][2] );
180  this->a22 = areaToInterpolate[0][0] -
181  ( 2.5 * areaToInterpolate[0][1] ) +
182  ( 2.0 * areaToInterpolate[0][2] ) -
183  ( 0.5 * areaToInterpolate[0][3] ) -
184  ( 2.5 * areaToInterpolate[1][0] ) +
185  ( 6.25 * areaToInterpolate[1][1] ) -
186  ( 5.0 * areaToInterpolate[1][2] ) +
187  ( 1.25 * areaToInterpolate[1][3] ) +
188  ( 2.0 * areaToInterpolate[2][0] ) -
189  ( 5.0 * areaToInterpolate[2][1] ) +
190  ( 4.0 * areaToInterpolate[2][2] ) -
191  areaToInterpolate[2][3] -
192  ( 0.5 * areaToInterpolate[3][0] ) +
193  ( 1.25 * areaToInterpolate[3][1] ) -
194  areaToInterpolate[3][2] +
195  ( 0.25 * areaToInterpolate[3][3] );
196  this->a23 = - ( 0.5 * areaToInterpolate[0][0] ) +
197  ( 1.5 * areaToInterpolate[0][1] ) -
198  ( 1.5 * areaToInterpolate[0][2] ) +
199  ( 0.5 * areaToInterpolate[0][3] ) +
200  ( 1.25 * areaToInterpolate[1][0] ) -
201  ( 3.75 * areaToInterpolate[1][1] ) +
202  ( 3.75 * areaToInterpolate[1][2] ) -
203  ( 1.25 * areaToInterpolate[1][3] ) -
204  areaToInterpolate[2][0] +
205  ( 3.0 * areaToInterpolate[2][1] ) -
206  ( 3.0 * areaToInterpolate[2][2] ) +
207  areaToInterpolate[2][3] +
208  ( 0.25 * areaToInterpolate[3][0] ) -
209  ( 0.75 * areaToInterpolate[3][1] ) +
210  ( 0.75 * areaToInterpolate[3][2] ) -
211  ( 0.25 * areaToInterpolate[3][3] );
212  this->a30 = - ( 0.5 * areaToInterpolate[0][1] ) +
213  ( 1.5 * areaToInterpolate[1][1] ) -
214  ( 1.5 * areaToInterpolate[2][1] ) +
215  ( 0.5*areaToInterpolate[3][1] );
216  this->a31 = ( 0.25 * areaToInterpolate[0][0] ) -
217  ( 0.25 * areaToInterpolate[0][2] ) -
218  ( 0.75 * areaToInterpolate[1][0] ) +
219  ( 0.75 * areaToInterpolate[1][2] ) +
220  ( 0.75 * areaToInterpolate[2][0] ) -
221  ( 0.75 * areaToInterpolate[2][2] ) -
222  ( 0.25 * areaToInterpolate[3][0] ) +
223  ( 0.25 * areaToInterpolate[3][2] );
224  this->a32 = - ( 0.5 * areaToInterpolate[0][0] ) +
225  ( 1.25 * areaToInterpolate[0][1] ) -
226  areaToInterpolate[0][2] +
227  ( 0.25 * areaToInterpolate[0][3] ) +
228  ( 1.5 * areaToInterpolate[1][0] ) -
229  ( 3.75 * areaToInterpolate[1][1] ) +
230  ( 3.0 * areaToInterpolate[1][2] ) -
231  ( 0.75 * areaToInterpolate[1][3] ) -
232  ( 1.5 * areaToInterpolate[2][0] ) +
233  ( 3.75 * areaToInterpolate[2][1] ) -
234  ( 3.0 * areaToInterpolate[2][2] ) +
235  ( 0.75 * areaToInterpolate[2][3] ) +
236  ( 0.5 * areaToInterpolate[3][0] ) -
237  ( 1.25 * areaToInterpolate[3][1] ) +
238  areaToInterpolate[3][2] -
239  ( 0.25 * areaToInterpolate[3][3] );
240  this->a33 = ( 0.25 * areaToInterpolate[0][0] ) -
241  ( 0.75 * areaToInterpolate[0][1] ) +
242  ( 0.75 * areaToInterpolate[0][2] ) -
243  ( 0.25 * areaToInterpolate[0][3] ) -
244  ( 0.75 * areaToInterpolate[1][0] ) +
245  ( 2.25 * areaToInterpolate[1][1] ) -
246  ( 2.25 * areaToInterpolate[1][2] ) +
247  ( 0.75 * areaToInterpolate[1][3] ) +
248  ( 0.75 * areaToInterpolate[2][0] ) -
249  ( 2.25 * areaToInterpolate[2][1] ) +
250  ( 2.25 * areaToInterpolate[2][2] ) -
251  ( 0.75 * areaToInterpolate[2][3] ) -
252  ( 0.25 * areaToInterpolate[3][0] ) +
253  ( 0.75 * areaToInterpolate[3][1] ) -
254  ( 0.75 * areaToInterpolate[3][2] ) +
255  ( 0.25 * areaToInterpolate[3][3] );
256  }
257 
262  ~BicubicInterpolator ( void ) { ; }
263 
273  proshade_double getValue ( proshade_double x, proshade_double y )
274  {
275  //======================================== Sanity check
276  if ( ( ( x < this->xStartIndex ) || ( x > ( this->xStartIndex + this->xRange ) ) ) ||
277  ( ( y < this->yStartIndex ) || ( y > ( this->yStartIndex + this->yRange ) ) ) )
278  {
279  if ( ( x < this->xStartIndex ) || ( x > ( this->xStartIndex + this->xRange ) ) ) { std::cout << "PROBLEM WITH LAT" << std::endl; }
280  if ( ( y < this->yStartIndex ) || ( y > ( this->yStartIndex + this->yRange ) ) ) { std::cout << "PROBLEM WITH LON" << std::endl; }
281 
282  throw ProSHADE_exception ( "Requested bicubic interpolation outside of pre-computed\n : square.", "ES00064", __FILE__, __LINE__, __func__, "The supplied x or y value(s) is outside of the range of\n : the bi-cubic interpolator's pre-computed square. Please\n : make sure the start values were correctly supplied when\n : the constructor was called or create a new interpolator\n : for these values." );
283  }
284 
285  //======================================== Convert x and y to unit square
286  proshade_double unitSquareX = ( x - this->xStartIndex ) / this->xRange;
287  proshade_double unitSquareY = ( y - this->yStartIndex ) / this->yRange;
288 
289  //======================================== Precompute powers
290  proshade_double x2 = std::pow ( unitSquareX, 2.0 );
291  proshade_double x3 = std::pow ( unitSquareX, 3.0 );
292  proshade_double y2 = std::pow ( unitSquareY, 2.0 );
293  proshade_double y3 = std::pow ( unitSquareY, 3.0 );
294 
295  //======================================== Done
296  return ( ( this->a00 + this->a01 * unitSquareY + this->a02 * y2 + this->a03 * y3 ) +
297  ( this->a10 + this->a11 * unitSquareY + this->a12 * y2 + this->a13 * y3 ) * unitSquareX +
298  ( this->a20 + this->a21 * unitSquareY + this->a22 * y2 + this->a23 * y3 ) * x2 +
299  ( this->a30 + this->a31 * unitSquareY + this->a32 * y2 + this->a33 * y3 ) * x3 );
300  }
301  };
302 
303  void complexMultiplication ( proshade_double* r1, proshade_double* i1, proshade_double* r2, proshade_double* i2,
304  proshade_double* retReal, proshade_double* retImag );
305  void complexMultiplicationConjug ( proshade_double* r1, proshade_double* i1, proshade_double* r2, proshade_double* i2,
306  proshade_double* retReal, proshade_double* retImag );
307  proshade_double complexMultiplicationRealOnly ( proshade_double* r1, proshade_double* i1, proshade_double* r2, proshade_double* i2 );
308  proshade_double complexMultiplicationConjugRealOnly ( proshade_double* r1, proshade_double* i1, proshade_double* r2, proshade_double* i2 );
309  void vectorMeanAndSD ( std::vector<proshade_double>* vec, proshade_double*& ret );
310  void vectorMedianAndIQR ( std::vector<proshade_double>* vec, proshade_double*& ret );
311  void arrayMedianAndIQR ( proshade_double* vec, proshade_unsign vecSize, proshade_double*& ret );
312  proshade_double pearsonCorrCoeff ( proshade_double* valSet1, proshade_double* valSet2, proshade_unsign length );
313  void getLegendreAbscAndWeights ( proshade_unsign order, proshade_double* abscissas, proshade_double* weights,
314  proshade_unsign noSteps );
315  void getGLPolyAtZero ( proshade_unsign order, proshade_double *polyValue, proshade_double *deriValue );
316  void getGLFirstRealRoot ( proshade_double polyAtZero, proshade_unsign order, proshade_double *abscAtZero,
317  proshade_double *weighAtZero, proshade_unsign taylorSeriesCap );
318  proshade_double evaluateGLPolynomial ( proshade_double *series, proshade_double target, proshade_unsign terms );
319  proshade_double advanceGLPolyValue ( proshade_double from, proshade_double to, proshade_double valAtFrom, proshade_unsign order, proshade_unsign noSteps );
320  void completeAbscissasAndWeights ( proshade_unsign order, proshade_double* abscissa, proshade_double* weights,
321  proshade_unsign noSteps );
322  proshade_double gaussLegendreIntegrationReal ( proshade_double* vals, proshade_unsign valsSize, proshade_unsign order,
323  proshade_double* abscissas, proshade_double* weights, proshade_double integralOverRange,
324  proshade_double maxSphereDists );
325  void gaussLegendreIntegration ( proshade_complex* vals, proshade_unsign valsSize, proshade_unsign order,
326  proshade_double* abscissas, proshade_double* weights, proshade_double integralOverRange,
327  proshade_double maxSphereDists, proshade_double* retReal, proshade_double* retImag );
328  void complexMatrixSVDSigmasOnly ( proshade_complex** mat, int dim, double*& singularValues );
329  void realMatrixSVDUandVOnly ( proshade_double* mat, int dim, proshade_double* uAndV, bool fail = true );
330  void getEulerZYZFromSOFTPosition ( proshade_signed band, proshade_signed x, proshade_signed y, proshade_signed z, proshade_double* eulerAlpha,
331  proshade_double* eulerBeta, proshade_double* eulerGamma );
332  void getSOFTPositionFromEulerZYZ ( proshade_signed band, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma,
333  proshade_double* x, proshade_double* y, proshade_double* z );
334  void getRotationMatrixFromEulerZYZAngles ( proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma, proshade_double* matrix );
335  void getRotationMatrixFromEulerZYZAngles ( proshade_single eulerAlpha, proshade_single eulerBeta, proshade_single eulerGamma, proshade_single* matrix );
336  void getAxisAngleFromRotationMatrix ( proshade_double* rotMat, proshade_double* x, proshade_double* y, proshade_double* z, proshade_double* ang, proshade_signed verbose = 1 );
337  void getAxisAngleFromRotationMatrix ( std::vector< proshade_double >* rotMat, proshade_double* x, proshade_double* y, proshade_double* z, proshade_double* ang, proshade_signed verbose = 1 );
338  void getRotationMatrixFromAngleAxis ( proshade_double* rotMat, proshade_double x, proshade_double y, proshade_double z, proshade_double ang );
339  void getRotationMatrixFromAngleAxis ( proshade_single* rotMat, proshade_double x, proshade_double y, proshade_double z, proshade_double ang );
340  void getEulerZYZFromRotMatrix ( proshade_double* rotMat, proshade_double* eA, proshade_double* eB, proshade_double* eG );
341  void getEulerZYZFromAngleAxis ( proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, proshade_double* eA,
342  proshade_double* eB, proshade_double* eG );
343  void multiplyTwoSquareMatrices ( proshade_double* A, proshade_double* B, proshade_double* res, proshade_unsign dim = 3 );
344  std::vector < proshade_signed > primeFactorsDecomp ( proshade_signed number );
345  proshade_double normalDistributionValue ( proshade_double mean, proshade_double standardDev, proshade_double value );
346  proshade_double computeDotProduct ( proshade_double* x1, proshade_double* y1, proshade_double* z1, proshade_double* x2, proshade_double* y2,
347  proshade_double* z2 );
348  proshade_double computeDotProduct ( proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2,
349  proshade_double z2 );
350  proshade_double* computeCrossProduct ( proshade_double* x1, proshade_double* y1, proshade_double* z1, proshade_double* x2, proshade_double* y2,
351  proshade_double* z2 );
352  proshade_double* computeCrossProduct ( proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2,
353  proshade_double z2 );
354  proshade_double* compute3x3MatrixMultiplication ( proshade_double* mat1, proshade_double* mat2 );
355  proshade_double* compute3x3MatrixVectorMultiplication ( proshade_double* mat, proshade_double x, proshade_double y, proshade_double z );
356  proshade_single* compute3x3MatrixVectorMultiplication ( proshade_single* mat, proshade_single x, proshade_single y, proshade_single z );
357  proshade_double* compute3x3MatrixInverse ( proshade_double* mat );
358  void transpose3x3MatrixInPlace ( proshade_single* mat );
359  void transpose3x3MatrixInPlace ( proshade_double* mat );
360  proshade_double* build3x3MatrixFromDiag ( proshade_double* diag );
361  proshade_double* build3x3MatrixFromXYZRotations ( proshade_double xRot, proshade_double yRot, proshade_double zRot );
362  proshade_double* findRotMatMatchingVectors ( proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2,
363  proshade_double z2 );
364  proshade_double* compute3x3MoorePenrosePseudoInverseOfIMinusMat ( std::vector < proshade_double >* rMat, proshade_signed verbose );
365  std::vector < proshade_double > findVectorFromTwoVAndTwoD ( proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2,
366  proshade_double z2, proshade_double dot1, proshade_double dot2 );
367  std::vector < proshade_double > findVectorFromThreeVAndThreeD ( proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2,
368  proshade_double z2, proshade_double x3, proshade_double y3, proshade_double z3, proshade_double dot1,
369  proshade_double dot2, proshade_double dot3 );
370  std::vector< proshade_double > multiplyGroupElementMatrices ( std::vector< proshade_double >* el1, std::vector< proshade_double >* el2 );
371  bool rotationMatrixSimilarity ( std::vector< proshade_double >* mat1, std::vector< proshade_double >* mat2, proshade_double tolerance = 0.1 );
372  bool rotationMatrixSimilarity ( proshade_double* mat1, proshade_double* mat2, proshade_double tolerance = 0.1 );
373  proshade_double rotationMatrixSimilarityValue ( std::vector< proshade_double >* mat1, std::vector< proshade_double >* mat2 );
374  proshade_double rotationMatrixSimilarityValue ( proshade_double* mat1, proshade_double* mat2 );
375  bool vectorOrientationSimilarity ( proshade_double a1, proshade_double a2, proshade_double a3, proshade_double b1, proshade_double b2,
376  proshade_double b3, proshade_double tolerance = 0.1 );
377  bool vectorOrientationSimilaritySameDirection ( proshade_double a1, proshade_double a2, proshade_double a3, proshade_double b1, proshade_double b2,
378  proshade_double b3, proshade_double tolerance = 0.1 );
379  void optimiseAxisBiCubicInterpolation ( proshade_double* bestLattitude, proshade_double* bestLongitude, proshade_double* bestSum, std::vector<proshade_unsign>* sphereList,
380  std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*>* sphereMappedRotFun, proshade_double step = 0.05 );
381  void prepareBiCubicInterpolatorsMinusMinus ( proshade_double bestLattitude, proshade_double bestLongitude, std::vector<proshade_unsign>* sphereList,
382  std::vector<ProSHADE_internal_maths::BicubicInterpolator*>* interpols, std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*>* sphereMappedRotFun );
383  void prepareBiCubicInterpolatorsMinusPlus ( proshade_double bestLattitude, proshade_double bestLongitude, std::vector<proshade_unsign>* sphereList,
384  std::vector<ProSHADE_internal_maths::BicubicInterpolator*>* interpols, std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*>* sphereMappedRotFun );
385  void prepareBiCubicInterpolatorsPlusMinus ( proshade_double bestLattitude, proshade_double bestLongitude, std::vector<proshade_unsign>* sphereList,
386  std::vector<ProSHADE_internal_maths::BicubicInterpolator*>* interpols, std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*>* sphereMappedRotFun );
387  void prepareBiCubicInterpolatorsPlusPlus ( proshade_double bestLattitude, proshade_double bestLongitude, std::vector<proshade_unsign>* sphereList,
388  std::vector<ProSHADE_internal_maths::BicubicInterpolator*>* interpols, std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*>* sphereMappedRotFun );
389  bool isAxisUnique ( std::vector< proshade_double* >* CSymList, proshade_double* axis, proshade_double tolerance = 0.1, bool improve = false );
390  bool isAxisUnique ( std::vector< proshade_double* >* CSymList, proshade_double X, proshade_double Y, proshade_double Z, proshade_double fold, proshade_double tolerance );
391  proshade_signed whichAxisUnique ( std::vector< proshade_double* >* CSymList, proshade_double* axis, proshade_double tolerance );
392  proshade_signed whichAxisUnique ( std::vector< proshade_double* >* CSymList, proshade_double X, proshade_double Y, proshade_double Z, proshade_double fold, proshade_double tolerance );
393  std::vector< proshade_unsign > findAllPrimes ( proshade_unsign upTo );
394  bool isPrime ( proshade_unsign toCheck );
395  proshade_double computeGaussian ( proshade_double val, proshade_double sigma );
396  std::vector < proshade_double > smoothen1D ( proshade_double step, proshade_signed windowSize, proshade_double sigma, std::vector< proshade_double > data, proshade_signed decRound = 2 );
397  proshade_single getResolutionOfReflection ( proshade_single h, proshade_single k, proshade_single l, proshade_single xDim, proshade_single yDim, proshade_single zDim );
398  void binReciprocalSpaceReflections ( proshade_unsign xInds, proshade_unsign yInds, proshade_unsign zInds, proshade_single xSize, proshade_single ySize, proshade_single zSize,
399  proshade_signed* noBin, proshade_signed*& binIndexing, std::vector< proshade_single >*& resArray );
400  void cutIndicesToResolution ( proshade_signed xInds, proshade_signed yInds, proshade_signed zInds, proshade_single resolution, proshade_signed* binIndexing,
401  std::vector< proshade_single >* resArray, proshade_signed* cutXDim, proshade_signed* cutYDim,
402  proshade_signed* cutZDim, proshade_signed*& cutBinIndices, proshade_signed*& noBins );
403  void cutArrayToResolution ( proshade_signed xInds, proshade_signed yInds, proshade_signed zInds, proshade_signed noBins, fftw_complex* inputMap, fftw_complex*& cutMap );
404  proshade_double computeFSC ( fftw_complex *fCoeffs1, fftw_complex *fCoeffs2, proshade_signed xInds, proshade_signed yInds, proshade_signed zInds,
405  proshade_signed noBins, proshade_signed* binIndexing, proshade_double**& binData, proshade_signed*& binCounts, proshade_double*& fscByBin,
406  bool averageByBinSize = false );
407  void computeFSCWeightByBin ( proshade_double*& weights1, proshade_double*& weights2, proshade_signed* binIndexing, proshade_double* fscByBin, proshade_signed noBins,
408  proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
409  proshade_double computeTheFValue ( proshade_complex* fCoeffs, proshade_double* weights, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim );
410  void computeTrFunDerivatives ( proshade_complex* fCoeffs, proshade_double* weights1, proshade_double* weights2, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim,
411  proshade_double*& firstDers, proshade_double*& secondDers );
412  proshade_double* computeTrFunStep ( proshade_double* firstDers, proshade_double* secondDers );
413  std::vector< proshade_signed > findPeaks1D ( std::vector< proshade_double > data );
414  proshade_double findTopGroupSmooth ( std::vector< proshade_double* >* CSym, size_t peakPos, proshade_double step, proshade_double sigma, proshade_signed windowSize, proshade_double maxLim = 1.0 );
415  proshade_double findTopGroupSmooth ( std::vector< std::vector< proshade_double > >* CSym, size_t peakPos, proshade_double step, proshade_double sigma, proshade_signed windowSize,
416  proshade_double maxLim = 1.0 );
417  void combineFourierForTranslation ( fftw_complex* tmpOut1, fftw_complex* tmpOut2, fftw_complex*& resOut, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD );
418  void findHighestValueInMap ( fftw_complex* resIn, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD, proshade_double* trsX,
419  proshade_double* trsY, proshade_double* trsZ, proshade_double* mapPeak );
420 }
421 
422 #endif
ProSHADE_internal_maths::getSOFTPositionFromEulerZYZ
void getSOFTPositionFromEulerZYZ(proshade_signed band, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma, proshade_double *x, proshade_double *y, proshade_double *z)
Function to find the index position in the inverse SOFT map from given Euler angles (ZYZ convention).
Definition: ProSHADE_maths.cpp:991
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_internal_maths::findVectorFromThreeVAndThreeD
std::vector< proshade_double > findVectorFromThreeVAndThreeD(proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2, proshade_double z2, proshade_double x3, proshade_double y3, proshade_double z3, proshade_double dot1, proshade_double dot2, proshade_double dot3)
Function for finding a vector which would have a given three dot products to three other vectors.
Definition: ProSHADE_maths.cpp:2479
ProSHADE_internal_maths::getEulerZYZFromSOFTPosition
void getEulerZYZFromSOFTPosition(proshade_signed band, proshade_signed x, proshade_signed y, proshade_signed z, proshade_double *eulerAlpha, proshade_double *eulerBeta, proshade_double *eulerGamma)
Function to find Euler angles (ZYZ convention) from index position in the inverse SOFT map.
Definition: ProSHADE_maths.cpp:966
ProSHADE_internal_maths::computeTrFunStep
proshade_double * computeTrFunStep(proshade_double *firstDers, proshade_double *secondDers)
This function computes the step sizes for translation function optimisation from the first and second...
Definition: ProSHADE_maths.cpp:4019
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_internal_maths::vectorMeanAndSD
void vectorMeanAndSD(std::vector< proshade_double > *vec, proshade_double *&ret)
Function to get vector mean and standard deviation.
Definition: ProSHADE_maths.cpp:121
ProSHADE_internal_maths::prepareBiCubicInterpolatorsMinusPlus
void prepareBiCubicInterpolatorsMinusPlus(proshade_double bestLattitude, proshade_double bestLongitude, std::vector< proshade_unsign > *sphereList, std::vector< ProSHADE_internal_maths::BicubicInterpolator * > *interpols, std::vector< ProSHADE_internal_spheres::ProSHADE_rotFun_sphere * > *sphereMappedRotFun)
This function prepares the interpolation objects for the bi-cubic interpolation.
Definition: ProSHADE_maths.cpp:2899
ProSHADE_internal_maths::isPrime
bool isPrime(proshade_unsign toCheck)
This function check is the supplied number is prime or not.
Definition: ProSHADE_maths.cpp:3259
ProSHADE_internal_maths::whichAxisUnique
proshade_signed whichAxisUnique(std::vector< proshade_double * > *CSymList, proshade_double *axis, proshade_double tolerance)
This function checks if new axis is unique, or already detected and returns the position of match or ...
Definition: ProSHADE_maths.cpp:3150
ProSHADE_internal_maths::computeTheFValue
proshade_double computeTheFValue(proshade_complex *fCoeffs, proshade_double *weights, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function computes the real part of the sum of all coefficients except where the weight is less t...
Definition: ProSHADE_maths.cpp:3892
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getMinRadius
proshade_double getMinRadius(void)
Accessor function for the private variable minimal radius.
Definition: ProSHADE_spheres.cpp:730
ProSHADE_internal_maths::rotationMatrixSimilarityValue
proshade_double rotationMatrixSimilarityValue(std::vector< proshade_double > *mat1, std::vector< proshade_double > *mat2)
This function computes the distance between two rotation matrices and returns it.
Definition: ProSHADE_maths.cpp:2624
ProSHADE_internal_maths::computeCrossProduct
proshade_double * computeCrossProduct(proshade_double *x1, proshade_double *y1, proshade_double *z1, proshade_double *x2, proshade_double *y2, proshade_double *z2)
Simple 3D vector cross product computation.
Definition: ProSHADE_maths.cpp:1817
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getSphereNumber
proshade_unsign getSphereNumber(void)
Accessor function for the private variable sphere number.
Definition: ProSHADE_spheres.cpp:752
ProSHADE_internal_maths::computeFSC
proshade_double computeFSC(fftw_complex *fCoeffs1, fftw_complex *fCoeffs2, proshade_signed xInds, proshade_signed yInds, proshade_signed zInds, proshade_signed noBins, proshade_signed *binIndexing, proshade_double **&binData, proshade_signed *&binCounts, proshade_double *&fscByBin, bool averageByBinSize=false)
This function computes the FSC.
Definition: ProSHADE_maths.cpp:3693
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere
This class contains all inputed data for the rotation function angle-axis converted spheres.
Definition: ProSHADE_maths.hpp:58
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
ProSHADE_internal_maths::normalDistributionValue
proshade_double normalDistributionValue(proshade_double mean, proshade_double standardDev, proshade_double value)
Function to the heiht of normal distribution given by mean and standard deviation for a given value.
Definition: ProSHADE_maths.cpp:1768
ProSHADE_internal_maths::prepareBiCubicInterpolatorsPlusMinus
void prepareBiCubicInterpolatorsPlusMinus(proshade_double bestLattitude, proshade_double bestLongitude, std::vector< proshade_unsign > *sphereList, std::vector< ProSHADE_internal_maths::BicubicInterpolator * > *interpols, std::vector< ProSHADE_internal_spheres::ProSHADE_rotFun_sphere * > *sphereMappedRotFun)
This function prepares the interpolation objects for the bi-cubic interpolation.
Definition: ProSHADE_maths.cpp:2954
ProSHADE_internal_maths::getGLFirstRealRoot
void getGLFirstRealRoot(proshade_double polyAtZero, proshade_unsign order, proshade_double *abscAtZero, proshade_double *weighAtZero, proshade_unsign taylorSeriesCap)
This function finds the first root for Legendre polynomials of odd order.
Definition: ProSHADE_maths.cpp:386
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::~ProSHADE_rotFun_sphere
~ProSHADE_rotFun_sphere(void)
Destructor for releasing memory from the ProSHADE_rotFun_sphere class.
Definition: ProSHADE_spheres.cpp:674
ProSHADE_internal_maths::getRotationMatrixFromEulerZYZAngles
void getRotationMatrixFromEulerZYZAngles(proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma, proshade_double *matrix)
Function to find the rotation matrix from Euler angles (ZYZ convention).
Definition: ProSHADE_maths.cpp:1019
ProSHADE_internal_maths::build3x3MatrixFromDiag
proshade_double * build3x3MatrixFromDiag(proshade_double *diag)
Function for building a 3x3 matrix from diagonal (and assuming zero padding).
Definition: ProSHADE_maths.cpp:2028
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getCopyOfValues
std::vector< std::vector< proshade_double > > getCopyOfValues(void)
Function for obtaining a copy of all sphere values.
Definition: ProSHADE_spheres.cpp:936
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getRepresentedAngle
proshade_double getRepresentedAngle(void)
Accessor function for the private variable represented angle.
Definition: ProSHADE_spheres.cpp:741
ProSHADE_internal_maths::BicubicInterpolator::getValue
proshade_double getValue(proshade_double x, proshade_double y)
This function allows accessing the interpolated value for a given position x and y.
Definition: ProSHADE_maths.hpp:273
ProSHADE_misc.hpp
This header file declares all miscellaneous functions.
ProSHADE_internal_maths::advanceGLPolyValue
proshade_double advanceGLPolyValue(proshade_double from, proshade_double to, proshade_double valAtFrom, proshade_unsign order, proshade_unsign noSteps)
This function finds the next value of a Legendre polynomial using the Runge-Kutta method.
Definition: ProSHADE_maths.cpp:480
ProSHADE_internal_maths::complexMultiplication
void complexMultiplication(proshade_double *r1, proshade_double *i1, proshade_double *r2, proshade_double *i2, proshade_double *retReal, proshade_double *retImag)
Function to multiply two complex numbers.
Definition: ProSHADE_maths.cpp:38
ProSHADE_internal_maths::BicubicInterpolator::~BicubicInterpolator
~BicubicInterpolator(void)
This is the destructor for the BicubicInterpolator class.
Definition: ProSHADE_maths.hpp:262
ProSHADE_internal_maths::optimiseAxisBiCubicInterpolation
void optimiseAxisBiCubicInterpolation(proshade_double *bestLattitude, proshade_double *bestLongitude, proshade_double *bestSum, std::vector< proshade_unsign > *sphereList, std::vector< ProSHADE_internal_spheres::ProSHADE_rotFun_sphere * > *sphereMappedRotFun, proshade_double step=0.05)
This function provides axis optimisation given starting lattitude and longitude indices.
Definition: ProSHADE_maths.cpp:2744
ProSHADE_internal_maths::vectorOrientationSimilarity
bool vectorOrientationSimilarity(proshade_double a1, proshade_double a2, proshade_double a3, proshade_double b1, proshade_double b2, proshade_double b3, proshade_double tolerance=0.1)
This function compares two vectors using cosine distance and decides if they are similar using tolera...
Definition: ProSHADE_maths.cpp:2680
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_internal_maths::prepareBiCubicInterpolatorsMinusMinus
void prepareBiCubicInterpolatorsMinusMinus(proshade_double bestLattitude, proshade_double bestLongitude, std::vector< proshade_unsign > *sphereList, std::vector< ProSHADE_internal_maths::BicubicInterpolator * > *interpols, std::vector< ProSHADE_internal_spheres::ProSHADE_rotFun_sphere * > *sphereMappedRotFun)
This function prepares the interpolation objects for the bi-cubic interpolation.
Definition: ProSHADE_maths.cpp:2844
ProSHADE_internal_maths::vectorOrientationSimilaritySameDirection
bool vectorOrientationSimilaritySameDirection(proshade_double a1, proshade_double a2, proshade_double a3, proshade_double b1, proshade_double b2, proshade_double b3, proshade_double tolerance=0.1)
This function compares two vectors using cosine distance and decides if they are similar using tolera...
Definition: ProSHADE_maths.cpp:2714
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getMaxRadius
proshade_double getMaxRadius(void)
Accessor function for the private variable maximum radius.
Definition: ProSHADE_spheres.cpp:697
ProSHADE_internal_maths::computeGaussian
proshade_double computeGaussian(proshade_double val, proshade_double sigma)
This function computes a Gaussian (normal) distribution value given distance from mean and sigma.
Definition: ProSHADE_maths.cpp:3286
ProSHADE_internal_maths::getEulerZYZFromRotMatrix
void getEulerZYZFromRotMatrix(proshade_double *rotMat, proshade_double *eA, proshade_double *eB, proshade_double *eG)
This function converts rotation matrix to the Euler ZYZ angles representation.
Definition: ProSHADE_maths.cpp:1561
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::ProSHADE_rotFun_sphere
ProSHADE_rotFun_sphere(proshade_double rad, proshade_double radRange, proshade_unsign dim, proshade_unsign rfDim, proshade_double repAng, proshade_unsign sphNo)
Constructor for getting empty ProSHADE_rotFun_sphere class.
Definition: ProSHADE_spheres.cpp:650
ProSHADE_internal_maths::cutIndicesToResolution
void cutIndicesToResolution(proshade_signed xInds, proshade_signed yInds, proshade_signed zInds, proshade_single resolution, proshade_signed *binIndexing, std::vector< proshade_single > *resArray, proshade_signed *cutXDim, proshade_signed *cutYDim, proshade_signed *cutZDim, proshade_signed *&cutBinIndices, proshade_signed *&noBins)
This function cuts the bin assignment array into a smaller array containing all bins up to a given re...
Definition: ProSHADE_maths.cpp:3563
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_maths::BicubicInterpolator::BicubicInterpolator
BicubicInterpolator(proshade_double **areaToInterpolate, proshade_double xStart, proshade_double yStart)
This is the constructor for the BicubicInterpolator class.
Definition: ProSHADE_maths.hpp:125
ProSHADE_internal_maths::computeDotProduct
proshade_double computeDotProduct(proshade_double *x1, proshade_double *y1, proshade_double *z1, proshade_double *x2, proshade_double *y2, proshade_double *z2)
Simple 3D vector dot product computation.
Definition: ProSHADE_maths.cpp:1785
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::findAllPeaks
void findAllPeaks(proshade_signed noSmNeighbours, std::vector< proshade_double > *allHeights)
Function for finding all peaks in the sampling grid.
Definition: ProSHADE_spheres.cpp:965
ProSHADE_internal_maths::findPeaks1D
std::vector< proshade_signed > findPeaks1D(std::vector< proshade_double > data)
This function simply finds all the peaks in a 1D data array.
Definition: ProSHADE_maths.cpp:4047
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::removeSmallPeaks
void removeSmallPeaks(proshade_double peakThres)
Function for removing peaks with too small height.
Definition: ProSHADE_spheres.cpp:1027
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_maths::getRotationMatrixFromAngleAxis
void getRotationMatrixFromAngleAxis(proshade_double *rotMat, proshade_double x, proshade_double y, proshade_double z, proshade_double ang)
This function converts the axis-angle representation to the rotation matrix representation.
Definition: ProSHADE_maths.cpp:1458
ProSHADE_internal_maths::findHighestValueInMap
void findHighestValueInMap(fftw_complex *resIn, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD, proshade_double *trsX, proshade_double *trsY, proshade_double *trsZ, proshade_double *mapPeak)
This function simply finds the highest value in fftw_complex map and returns its position and value.
Definition: ProSHADE_maths.cpp:4289
ProSHADE_internal_maths::getAxisAngleFromRotationMatrix
void getAxisAngleFromRotationMatrix(proshade_double *rotMat, proshade_double *x, proshade_double *y, proshade_double *z, proshade_double *ang, proshade_signed verbose=1)
This function converts rotation matrix to the axis-angle representation.
Definition: ProSHADE_maths.cpp:1083
ProSHADE_internal_maths::build3x3MatrixFromXYZRotations
proshade_double * build3x3MatrixFromXYZRotations(proshade_double xRot, proshade_double yRot, proshade_double zRot)
Function for building a 3x3 rotation matrix from the x, y and z rotations in degrees.
Definition: ProSHADE_maths.cpp:2057
ProSHADE_internal_maths::getEulerZYZFromAngleAxis
void getEulerZYZFromAngleAxis(proshade_double axX, proshade_double axY, proshade_double axZ, proshade_double axAng, proshade_double *eA, proshade_double *eB, proshade_double *eG)
This function converts angle-axis representation to the Euler ZYZ angles representation.
Definition: ProSHADE_maths.cpp:1612
ProSHADE_internal_maths::realMatrixSVDUandVOnly
void realMatrixSVDUandVOnly(proshade_double *mat, int dim, proshade_double *uAndV, bool fail=true)
Function to compute the real matrix SVD and return the U and V matrices.
Definition: ProSHADE_maths.cpp:874
ProSHADE_internal_maths
This namespace contains the internal functions for common mathematical operations.
ProSHADE_internal_maths::combineFourierForTranslation
void combineFourierForTranslation(fftw_complex *tmpOut1, fftw_complex *tmpOut2, fftw_complex *&resOut, proshade_unsign xD, proshade_unsign yD, proshade_unsign zD)
This function combines Fourier coefficients of two structures in a way, so that inverse Fourier of th...
Definition: ProSHADE_maths.cpp:4242
ProSHADE_internal_maths::BicubicInterpolator
Definition: ProSHADE_maths.hpp:112
ProSHADE_internal_maths::primeFactorsDecomp
std::vector< proshade_signed > primeFactorsDecomp(proshade_signed number)
Function to find prime factors of an integer.
Definition: ProSHADE_maths.cpp:1719
ProSHADE_internal_maths::getGLPolyAtZero
void getGLPolyAtZero(proshade_unsign order, proshade_double *polyValue, proshade_double *deriValue)
This function obtains the Legendre polynomial values and its derivative at zero for any positive inte...
Definition: ProSHADE_maths.cpp:349
ProSHADE_internal_maths::isAxisUnique
bool isAxisUnique(std::vector< proshade_double * > *CSymList, proshade_double *axis, proshade_double tolerance=0.1, bool improve=false)
This function checks if new axis is unique, or already detected.
Definition: ProSHADE_maths.cpp:3064
ProSHADE_internal_maths::completeAbscissasAndWeights
void completeAbscissasAndWeights(proshade_unsign order, proshade_double *abscissa, proshade_double *weights, proshade_unsign noSteps)
This function completes the abscissas and weights series from the first roots computed beforehand.
Definition: ProSHADE_maths.cpp:524
ProSHADE_internal_maths::findRotMatMatchingVectors
proshade_double * findRotMatMatchingVectors(proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2, proshade_double z2)
Computation of rotation matrix rotating one vector onto the other.
Definition: ProSHADE_maths.cpp:2119
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getSphereLatLonPosition
proshade_double getSphereLatLonPosition(proshade_unsign lattitude, proshade_unsign longitude)
Accessor function for specific lattitude and longitude point of the sphere sampling grid.
Definition: ProSHADE_spheres.cpp:877
ProSHADE_internal_maths::compute3x3MatrixVectorMultiplication
proshade_double * compute3x3MatrixVectorMultiplication(proshade_double *mat, proshade_double x, proshade_double y, proshade_double z)
Function for computing a 3x3 matrix to 3x1 vector multiplication.
Definition: ProSHADE_maths.cpp:1895
ProSHADE_internal_spheres
This namespace contains the structure and functions required for storing internal map projections ont...
Definition: ProSHADE_maths.hpp:49
ProSHADE_internal_maths::compute3x3MoorePenrosePseudoInverseOfIMinusMat
proshade_double * compute3x3MoorePenrosePseudoInverseOfIMinusMat(std::vector< proshade_double > *rMat, proshade_signed verbose)
This function computes the Moore-Penrose pseudo-inverse of equation I - input matrix.
Definition: ProSHADE_maths.cpp:2186
ProSHADE_internal_maths::rotationMatrixSimilarity
bool rotationMatrixSimilarity(std::vector< proshade_double > *mat1, std::vector< proshade_double > *mat2, proshade_double tolerance=0.1)
This function compares the distance between two rotation matrices and decides if they are similar usi...
Definition: ProSHADE_maths.cpp:2576
ProSHADE_internal_maths::vectorMedianAndIQR
void vectorMedianAndIQR(std::vector< proshade_double > *vec, proshade_double *&ret)
Function to get vector median and inter-quartile range.
Definition: ProSHADE_maths.cpp:149
ProSHADE_internal_maths::computeTrFunDerivatives
void computeTrFunDerivatives(proshade_complex *fCoeffs, proshade_double *weights1, proshade_double *weights2, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_double *&firstDers, proshade_double *&secondDers)
This function computes the first and second derivatives of the translation function at coefficient [0...
Definition: ProSHADE_maths.cpp:3930
ProSHADE_internal_maths::cutArrayToResolution
void cutArrayToResolution(proshade_signed xInds, proshade_signed yInds, proshade_signed zInds, proshade_signed noBins, fftw_complex *inputMap, fftw_complex *&cutMap)
This function re-sizes data array to contain only values up to a particular resolution bin.
Definition: ProSHADE_maths.cpp:3628
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_maths::smoothen1D
std::vector< proshade_double > smoothen1D(proshade_double step, proshade_signed windowSize, proshade_double sigma, std::vector< proshade_double > data, proshade_signed decRound=2)
This function takes a 1D vector and computes smoothened version based on the parameters.
Definition: ProSHADE_maths.cpp:3313
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_internal_spheres::ProSHADE_rotFun_sphere::getRadius
proshade_double getRadius(void)
Accessor function for the private variable radius.
Definition: ProSHADE_spheres.cpp:687
ProSHADE_internal_maths::findVectorFromTwoVAndTwoD
std::vector< proshade_double > findVectorFromTwoVAndTwoD(proshade_double x1, proshade_double y1, proshade_double z1, proshade_double x2, proshade_double y2, proshade_double z2, proshade_double dot1, proshade_double dot2)
Function for finding a vector which would have a given two dot products to two other vectors.
Definition: ProSHADE_maths.cpp:2333
ProSHADE_internal_maths::prepareBiCubicInterpolatorsPlusPlus
void prepareBiCubicInterpolatorsPlusPlus(proshade_double bestLattitude, proshade_double bestLongitude, std::vector< proshade_unsign > *sphereList, std::vector< ProSHADE_internal_maths::BicubicInterpolator * > *interpols, std::vector< ProSHADE_internal_spheres::ProSHADE_rotFun_sphere * > *sphereMappedRotFun)
This function prepares the interpolation objects for the bi-cubic interpolation.
Definition: ProSHADE_maths.cpp:3009
ProSHADE_internal_maths::multiplyTwoSquareMatrices
void multiplyTwoSquareMatrices(proshade_double *A, proshade_double *B, proshade_double *res, proshade_unsign dim=3)
Function to compute matrix multiplication.
Definition: ProSHADE_maths.cpp:1693
ProSHADE_internal_maths::getResolutionOfReflection
proshade_single getResolutionOfReflection(proshade_single h, proshade_single k, proshade_single l, proshade_single xDim, proshade_single yDim, proshade_single zDim)
This function computes the resolution of a particular reflection.
Definition: ProSHADE_maths.cpp:3359
ProSHADE_internal_maths::computeFSCWeightByBin
void computeFSCWeightByBin(proshade_double *&weights1, proshade_double *&weights2, proshade_signed *binIndexing, proshade_double *fscByBin, proshade_signed noBins, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim)
This function computes the weights for each reflection using its bin belonging.
Definition: ProSHADE_maths.cpp:3800
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getSphereLatLonLinearInterpolationPos
proshade_double getSphereLatLonLinearInterpolationPos(proshade_double lattitude, proshade_double longitude)
Function for obtaining sphere values outside of the grid points.
Definition: ProSHADE_spheres.cpp:891
ProSHADE_internal_maths::findTopGroupSmooth
proshade_double findTopGroupSmooth(std::vector< proshade_double * > *CSym, size_t peakPos, proshade_double step, proshade_double sigma, proshade_signed windowSize, proshade_double maxLim=1.0)
This function finds a subgroup of axes with distinctly higher correlation value.
Definition: ProSHADE_maths.cpp:4113
ProSHADE_internal_maths::findAllPrimes
std::vector< proshade_unsign > findAllPrimes(proshade_unsign upTo)
This function finds all prime numbers up to the supplied limit.
Definition: ProSHADE_maths.cpp:3218
ProSHADE_internal_maths::compute3x3MatrixInverse
proshade_double * compute3x3MatrixInverse(proshade_double *mat)
Function for computing a 3x3 matrix inverse.
Definition: ProSHADE_maths.cpp:1940
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getRotFunDim
proshade_unsign getRotFunDim(void)
Accessor function for the private variable rotastion function dim.
Definition: ProSHADE_spheres.cpp:708
ProSHADE_internal_maths::arrayMedianAndIQR
void arrayMedianAndIQR(proshade_double *vec, proshade_unsign vecSize, proshade_double *&ret)
Function to get array median and inter-quartile range.
Definition: ProSHADE_maths.cpp:200
ProSHADE_internal_maths::binReciprocalSpaceReflections
void binReciprocalSpaceReflections(proshade_unsign xInds, proshade_unsign yInds, proshade_unsign zInds, proshade_single xSize, proshade_single ySize, proshade_single zSize, proshade_signed *noBin, proshade_signed *&binIndexing, std::vector< proshade_single > *&resArray)
This function does binning of the reciprocal space reflections.
Definition: ProSHADE_maths.cpp:3397
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getAngularDim
proshade_unsign getAngularDim(void)
Accessor function for the private variable angular dim.
Definition: ProSHADE_spheres.cpp:719
ProSHADE_internal_maths::multiplyGroupElementMatrices
std::vector< proshade_double > multiplyGroupElementMatrices(std::vector< proshade_double > *el1, std::vector< proshade_double > *el2)
This function computes matrix multiplication using the ProSHADE group element matrix format as input ...
Definition: ProSHADE_maths.cpp:2523
ProSHADE_internal_maths::evaluateGLPolynomial
proshade_double evaluateGLPolynomial(proshade_double *series, proshade_double target, proshade_unsign terms)
This function evaluates the decomposed Legendre polynomial at a given position.
Definition: ProSHADE_maths.cpp:450
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::interpolateSphereValues
void interpolateSphereValues(proshade_complex *rotFun)
Function for interpolating the sphere grid values from angle-axis converted rotation function.
Definition: ProSHADE_spheres.cpp:779
ProSHADE_internal_spheres::ProSHADE_rotFun_sphere::getPeaks
std::vector< std::pair< proshade_unsign, proshade_unsign > > getPeaks(void)
Accessor function for the private variable containing all detected peaks.
Definition: ProSHADE_spheres.cpp:763
ProSHADE_internal_maths::transpose3x3MatrixInPlace
void transpose3x3MatrixInPlace(proshade_single *mat)
Transposes 3x3 matrix in place.
Definition: ProSHADE_maths.cpp:1973
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_internal_maths::compute3x3MatrixMultiplication
proshade_double * compute3x3MatrixMultiplication(proshade_double *mat1, proshade_double *mat2)
Function for computing a 3x3 matrix multiplication.
Definition: ProSHADE_maths.cpp:1865