ProSHADE  0.7.6.0 (JUL 2021)
Protein Shape Detection
ProSHADE_internal_sphericalHarmonics Namespace Reference

This namespace contains the internal functions for computing spherical harmonics and their related computations. More...

Functions

void allocateComputationMemory (proshade_unsign band, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&outputReal, proshade_double *&outputImag, double *&shWeights, double *&tableSpaceHelper, fftw_complex *&workspace)
 This function determines the integration order for the between spheres integration. More...
 
void placeWithinWorkspacePointers (fftw_complex *&workspace, proshade_unsign oDim, proshade_double *&rres, proshade_double *&ires, proshade_double *&fltres, proshade_double *&scratchpad)
 This function takes the workspace pointer and correctly places the other internal pointers. More...
 
void initialiseFFTWPlans (proshade_unsign band, fftw_plan &fftPlan, fftw_plan &dctPlan, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&rres, proshade_double *&ires, proshade_double *&scratchpad)
 This function initialises the FFTW plans. More...
 
void releaseSphericalMemory (proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&outputReal, proshade_double *&outputImag, double *&tableSpaceHelper, double **&tableSpace, double *&shWeights, fftw_complex *&workspace, fftw_plan &fftPlan, fftw_plan &dctPlan)
 This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer. More...
 
void initialiseAllMemory (proshade_unsign band, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&outputReal, proshade_double *&outputImag, double *&shWeights, double **&tableSpace, double *&tableSpaceHelper, fftw_complex *&workspace, proshade_double *&rres, proshade_double *&ires, proshade_double *&fltres, proshade_double *&scratchpad, fftw_plan &fftPlan, fftw_plan &dctPlan)
 This function initialises all the memory required for spherical harmonics computation. More...
 
void initialSplitDiscreteTransform (proshade_unsign oneDim, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&rres, proshade_double *&ires, proshade_double *mappedData, fftw_plan &fftPlan, proshade_double normCoeff)
 This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer. More...
 
void computeSphericalTransformCoeffs (proshade_unsign band, proshade_double *&rdataptr, proshade_double *&idataptr, proshade_double *&outputReal, proshade_double *&outputImag, proshade_double *&rres, proshade_double *&ires, proshade_double *&fltres, proshade_double *&scratchpad, double **&tablePml, double *&shWeights, fftw_plan &dctPlan)
 This function takes the split discrete transform and proceeds to complete the spherical harmonics decomposition. More...
 
void applyCondonShortleyPhase (proshade_unsign band, proshade_double *outputReal, proshade_double *outputImag, proshade_complex *&shArray)
 This is the final step in computing the full spherical harmonics decomposition of the input data. More...
 
void computeSphericalHarmonics (proshade_unsign band, proshade_double *sphereMappedData, proshade_complex *&shArray)
 This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer. More...
 

Detailed Description

This namespace contains the internal functions for computing spherical harmonics and their related computations.

The ProSHADE_internal_sphericalHarmonics namespace contains helper functions for spherical harmonics computation for each sphere as created in the sphere object as well as the further processing computations. None of these functions should be used directly be the user.

Function Documentation

◆ allocateComputationMemory()

void ProSHADE_internal_sphericalHarmonics::allocateComputationMemory ( proshade_unsign  band,
proshade_double *&  inputReal,
proshade_double *&  inputImag,
proshade_double *&  outputReal,
proshade_double *&  outputImag,
double *&  shWeights,
double *&  tableSpaceHelper,
fftw_complex *&  workspace 
)

This function determines the integration order for the between spheres integration.

This function simply takes all pointer variables required for the spherical harmonics computation and allocates the required amount of memory for them. It also does the memory checks in case memory allocation fails.

Parameters
[in]bandThe bandwidth to which the computation will be done.
[in]inputRealThe real part of the input will be copied here.
[in]inputImagThe immaginary part of the input will be copied here.
[in]outputRealThe real part of the output will be saved here.
[in]outputImagThe immaginary part of the output will be saved here.
[in]shWeightsThe weights for spherical harmonics computation will be stored here.
[in]tableSpaceHelperThis space is required by SOFT for pre-computing values into this table.
[in]workspaceThe space where multiple minor results are saved by SOFT.

Definition at line 39 of file ProSHADE_sphericalHarmonics.cpp.

40 {
41  //================================================ Initialise local variables
42  proshade_unsign oneDimmension = 2 * band;
43 
44  //================================================ Allocate Input Memory
45  inputReal = new proshade_double [oneDimmension * oneDimmension];
46  inputImag = new proshade_double [oneDimmension * oneDimmension];
47 
48  //================================================ Allocate Output Memory
49  outputReal = new proshade_double [oneDimmension * oneDimmension];
50  outputImag = new proshade_double [oneDimmension * oneDimmension];
51 
52  //================================================ Allocate Working Memory
53  shWeights = new proshade_double [band * 4];
54  workspace = new fftw_complex [( 8 * band * band ) + ( 10 * band )];
55 
56  //================================================ Allocate table
57  tableSpaceHelper = new proshade_double [static_cast<proshade_unsign> ( Reduced_Naive_TableSize ( static_cast< int > ( band ), static_cast< int > ( band ) ) +
58  Reduced_SpharmonicTableSize ( static_cast< int > ( band ), static_cast< int > ( band ) ) )];
59 
60  //================================================ Check memory allocation success
61  ProSHADE_internal_misc::checkMemoryAllocation ( inputReal, __FILE__, __LINE__, __func__ );
62  ProSHADE_internal_misc::checkMemoryAllocation ( inputImag, __FILE__, __LINE__, __func__ );
63  ProSHADE_internal_misc::checkMemoryAllocation ( outputReal, __FILE__, __LINE__, __func__ );
64  ProSHADE_internal_misc::checkMemoryAllocation ( outputImag, __FILE__, __LINE__, __func__ );
65  ProSHADE_internal_misc::checkMemoryAllocation ( shWeights, __FILE__, __LINE__, __func__ );
66  ProSHADE_internal_misc::checkMemoryAllocation ( tableSpaceHelper, __FILE__, __LINE__, __func__ );
67  ProSHADE_internal_misc::checkMemoryAllocation ( workspace, __FILE__, __LINE__, __func__ );
68 
69  //================================================ Done
70  return ;
71 
72 }

◆ applyCondonShortleyPhase()

void ProSHADE_internal_sphericalHarmonics::applyCondonShortleyPhase ( proshade_unsign  band,
proshade_double *  outputReal,
proshade_double *  outputImag,
proshade_complex *&  shArray 
)

This is the final step in computing the full spherical harmonics decomposition of the input data.

This is the final function that is needed for the complete spherical harmonics decomposition. It applied the Condon-Shortley phase as well as computing the negative orders, then saving the output into the final results array for further processing.

Parameters
[in]bandThe bandwidth to which the computation will be done.
[in]outputRealThe real results of the complete transform as done by the initialSplitDiscreteTransform() and computeSphericalTransformCoeffs() functions.
[in]outputImagThe imaginary results of the complete transform as done by the initialSplitDiscreteTransform() and computeSphericalTransformCoeffs() functions.
[in]shArrayAn array of complex numbers to which the results of the spherical harmonics decomposition are to be saved.

Definition at line 352 of file ProSHADE_sphericalHarmonics.cpp.

353 {
354  //================================================ Copy the results into the final holder
355  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( (band * 2) * (band * 2) ); iter++ )
356  {
357  shArray[iter][0] = outputReal[iter];
358  shArray[iter][1] = outputImag[iter];
359  }
360 
361  //================================================ Apply the Condon-Shortley phase sign
362  proshade_double powerOne = 1.0;
363  proshade_unsign hlp1 = 0;
364  proshade_unsign hlp2 = 0;
365  for ( proshade_signed order = 1; order < static_cast<proshade_signed> ( band ); order++)
366  {
367  powerOne *= -1.0;
368  for ( proshade_signed bandIter = order; bandIter < static_cast<proshade_signed> ( band ); bandIter++)
369  {
370  hlp1 = static_cast< proshade_unsign > ( seanindex ( static_cast< int > ( order ), static_cast< int > ( bandIter ), static_cast< int > ( band ) ) );
371  hlp2 = static_cast< proshade_unsign > ( seanindex ( static_cast< int > ( -order ), static_cast< int > ( bandIter ), static_cast< int > ( band ) ) );
372 
373  shArray[hlp2][0] = powerOne * static_cast<proshade_double> ( outputReal[hlp1] );
374  shArray[hlp2][1] = -powerOne * static_cast<proshade_double> ( outputImag[hlp1] );
375  }
376  }
377 
378  //================================================ DONE
379  return ;
380 
381 }

◆ computeSphericalHarmonics()

void ProSHADE_internal_sphericalHarmonics::computeSphericalHarmonics ( proshade_unsign  band,
proshade_double *  sphereMappedData,
proshade_complex *&  shArray 
)

This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer.

This function does all the spherical harmonics computations for a single shell, including the memory allocation and releasing and the FFTW transforms. Because the shells can have different resolutions, the memory management is left until here, but possible speed-up could be gained from having the same resolution on all shells as in the older ProSHADE versions.

Parameters
[in]bandThe bandwidth to which the computation will be done.
[in]sphereMappedDataAn array of doubles containing the mapped data onto a sphere for the sphere to be decomposed.
[in]shArrayAn array of complex numbers to which the results of the spherical harmonics decomposition are to be saved.

Definition at line 394 of file ProSHADE_sphericalHarmonics.cpp.

395 {
396  //================================================ Initialise local variables
397  proshade_double *inputReal = nullptr, *inputImag = nullptr, *outputReal = nullptr, *outputImag = nullptr;
398  double *shWeights = nullptr, *tableSpaceHelper = nullptr;
399  double** tablePml = nullptr;
400  fftw_complex* workspace = nullptr;
401  proshade_unsign oneDim = static_cast<proshade_unsign> ( band * 2 );
402  proshade_double normCoeff = ( 1.0 / ( static_cast<proshade_double> ( band * 2 ) ) ) * sqrt( 2.0 * M_PI );
403 
404  //================================================ Set output to zeroes (so that all unfilled data are not random)
405  for ( proshade_unsign i = 0; i < ( 2 * band * 2 * band); i++ )
406  {
407  shArray[i][0] = 0.0;
408  shArray[i][1] = 0.0;
409  }
410 
411  //================================================ Within workspace pointers
412  proshade_double *rres = nullptr, *ires = nullptr, *fltres = nullptr, *scratchpad = nullptr, *rdataptr = nullptr, *idataptr = nullptr;
413 
414  //================================================ FFTW Plans
415  fftw_plan fftPlan = nullptr;
416  fftw_plan dctPlan = nullptr;
417 
418  //================================================ Initialise all memory
419  initialiseAllMemory ( band, inputReal, inputImag, outputReal, outputImag, shWeights, tablePml, tableSpaceHelper, workspace,
420  rres, ires, fltres, scratchpad, fftPlan, dctPlan );
421 
422  //================================================ Do the initial discrete split transform
423  initialSplitDiscreteTransform ( oneDim, inputReal, inputImag, rres, ires, sphereMappedData, fftPlan, normCoeff );
424 
425  //================================================ Complete the spherical harmonics transform
426  computeSphericalTransformCoeffs ( band, rdataptr, idataptr, outputReal, outputImag, rres, ires, fltres, scratchpad, tablePml, shWeights, dctPlan );
427 
428  //================================================ Apply the Condon-Shortley phase and save result to the final array
429  applyCondonShortleyPhase ( band, outputReal, outputImag, shArray );
430 
431  //================================================ Free memory
432  releaseSphericalMemory ( inputReal, inputImag, outputReal, outputImag, tableSpaceHelper, tablePml, shWeights, workspace, fftPlan, dctPlan );
433 
434  //================================================ Done
435  return ;
436 
437 }

◆ computeSphericalTransformCoeffs()

void ProSHADE_internal_sphericalHarmonics::computeSphericalTransformCoeffs ( proshade_unsign  band,
proshade_double *&  rdataptr,
proshade_double *&  idataptr,
proshade_double *&  outputReal,
proshade_double *&  outputImag,
proshade_double *&  rres,
proshade_double *&  ires,
proshade_double *&  fltres,
proshade_double *&  scratchpad,
double **&  tablePml,
double *&  shWeights,
fftw_plan &  dctPlan 
)

This function takes the split discrete transform and proceeds to complete the spherical harmonics decomposition.

This function takes the results of the initial split discrete transform and proceeds to compute the spherical harmonics coefficients for all applicable bands using all the pre-computed valus (i.e. the Legendre polynomials table and the weights). To do this, the SOFT2.0 function is called and for more details, see SOFT2.0.

Parameters
[in]bandThe bandwidth to which the computation will be done.
[in]rdataptrPointer to be used as iterative locator in the large results array for real results.
[in]idataptrPointer to be used as iterative locator in the large results array for imaginary results.
[in]outputRealAn array for storing the real results of the completed transform.
[in]outputImagAn array for storing the imaginary results of the completed transform.
[in]rresArray containing the real results of the initial split discrete transform.
[in]iresArray containing the imaginary results of the initial split discrete transform.
[in]fltresHelper array for transform computations.
[in]scratchpadArray for keeping temporary results of the transform computations.
[in]tablePmlPre-computed array of the Legendre polynomials as done by SOFT2.0.
[in]shWeightsThe weights for the spherical harmonics.
[in]dctPlanThe FFTW plan for the final spherical harmonics transform.

Definition at line 300 of file ProSHADE_sphericalHarmonics.cpp.

301 {
302  //================================================ Calculate the coefficients for each band
303  rdataptr = outputReal;
304  idataptr = outputImag;
305  for ( proshade_unsign bandIter = 0; bandIter < band; bandIter++ )
306  {
307  //============================================ Real part calculation
308  SemiNaiveReduced ( rres + ( bandIter * ( band * 2 ) ),
309  static_cast< int > ( band ),
310  static_cast< int > ( bandIter ),
311  fltres,
312  scratchpad,
313  tablePml[bandIter],
314  shWeights,
315  &dctPlan);
316 
317  //============================================ Save the real results to temporary holder
318  memcpy ( rdataptr, fltres, sizeof(proshade_double) * ( band - bandIter ) );
319  rdataptr += band - bandIter;
320 
321  //============================================ Imaginary part calculation
322  SemiNaiveReduced ( ires + ( bandIter * ( band * 2 ) ),
323  static_cast< int > ( band ),
324  static_cast< int > ( bandIter ),
325  fltres,
326  scratchpad,
327  tablePml[bandIter],
328  shWeights,
329  &dctPlan );
330 
331  //============================================ Save the imaginary results
332  memcpy ( idataptr, fltres, sizeof(proshade_double) * ( band - bandIter ) );
333  idataptr += band - bandIter;
334  }
335 
336  //================================================ DONE
337  return ;
338 
339 }

◆ initialiseAllMemory()

void ProSHADE_internal_sphericalHarmonics::initialiseAllMemory ( proshade_unsign  band,
proshade_double *&  inputReal,
proshade_double *&  inputImag,
proshade_double *&  outputReal,
proshade_double *&  outputImag,
double *&  shWeights,
double **&  tableSpace,
double *&  tableSpaceHelper,
fftw_complex *&  workspace,
proshade_double *&  rres,
proshade_double *&  ires,
proshade_double *&  fltres,
proshade_double *&  scratchpad,
fftw_plan &  fftPlan,
fftw_plan &  dctPlan 
)

This function initialises all the memory required for spherical harmonics computation.

This function takes on all the memory allocation and filling in all data required later for the spherical harmonics computation by the SOFT2.0 library.

Parameters
[in]bandThe bandwidth to which the computation will be done.
[in]inputRealThe real part of the input will be copied here.
[in]inputImagThe immaginary part of the input will be copied here.
[in]outputRealThe real part of the output will be saved here.
[in]outputImagThe immaginary part of the output will be saved here.
[in]shWeightsThe weights for spherical harmonics computation will be stored here.
[in]tableSpaceThis space is required by SOFT for pre-computing values into this table.
[in]tableSpaceHelperThis space is required by SOFT for proper computation of the table space.
[in]workspaceThe space where multiple minor results are saved by SOFT2.0.
[in]rresPointer to where the real part of the results will be temporarily saved.
[in]iresPointer to where the imaginary part of the results will be temporarily saved.
[in]fltresPointer to where the temporary bandwise results will be saved.
[in]scratchpadPointer to extra space which is internally used (but not allocated) by SOFT2.0.
[in]fftPlanpointer to the variable where the Fourier transform should be set.
[in]dctPlanpointer to the variable where the 1D r2r Fourier transform should be set.

Definition at line 217 of file ProSHADE_sphericalHarmonics.cpp.

218 {
219  //================================================ Initialise local variables
220  proshade_unsign oneDim = band * 2;
221 
222  //================================================ Allocate memory for local pointers
223  allocateComputationMemory ( band, inputReal, inputImag, outputReal, outputImag, shWeights, tableSpaceHelper, workspace );
224 
225  //================================================ Within workspace pointers
226  placeWithinWorkspacePointers ( workspace, oneDim, rres, ires, fltres, scratchpad );
227 
228  //================================================ Generate Seminaive and naive tables for Legendre Polynomials
229  tableSpace = SemiNaive_Naive_Pml_Table ( static_cast< int > ( band ), static_cast< int > ( band ), tableSpaceHelper, reinterpret_cast<double*> ( workspace ) );
230 
231  //================================================ Make weights for spherical transform
232  makeweights ( static_cast< int > ( band ), shWeights );
233 
234  //================================================ Initialize FFTW Plans
235  initialiseFFTWPlans ( band, fftPlan, dctPlan, inputReal, inputImag, rres, ires, scratchpad );
236 
237  //================================================ Done
238  return ;
239 
240 }

◆ initialiseFFTWPlans()

void ProSHADE_internal_sphericalHarmonics::initialiseFFTWPlans ( proshade_unsign  band,
fftw_plan &  fftPlan,
fftw_plan &  dctPlan,
proshade_double *&  inputReal,
proshade_double *&  inputImag,
proshade_double *&  rres,
proshade_double *&  ires,
proshade_double *&  scratchpad 
)

This function initialises the FFTW plans.

This function initialises the FFTW plans for the spherical harmonics computations as required by the SOFT2.0 library.

Parameters
[in]bandThe bandwidth to which the computation will be done.
[in]fftPlanpointer to the variable where the Fourier transform should be set.
[in]dctPlanpointer to the variable where the 1D r2r Fourier transform should be set.
[in]inputRealpointer to the array containing (or which will contain) the input real values.
[in]inputImagpointer to the array containing (or which will contain) the input imaginary values.
[in]rrespointer to the array where the real values of result should be saved.
[in]irespointer to the array where the imaginary values of result should be saved.
[in]scratchpadpointer to the array where temporary results will be saved.

Definition at line 113 of file ProSHADE_sphericalHarmonics.cpp.

114 {
115  //================================================ Initialize fft plan along phi angles
116  fftw_iodim dims[1];
117  fftw_iodim howmany_dims[1];
118 
119  int rank = 1;
120  int howmany_rank = 1;
121 
122  dims[0].n = static_cast<int> ( band * 2 );
123  dims[0].is = 1;
124  dims[0].os = static_cast<int> ( band * 2 );
125 
126  howmany_dims[0].n = static_cast<int> ( band * 2 );
127  howmany_dims[0].is = static_cast<int> ( band * 2 );
128  howmany_dims[0].os = 1;
129 
130  //================================================ Plan fft transform
131  fftPlan = fftw_plan_guru_split_dft ( rank,
132  dims,
133  howmany_rank,
134  howmany_dims,
135  inputReal,
136  inputImag,
137  rres,
138  ires,
139  FFTW_ESTIMATE );
140 
141  //================================================ Initialize dct plan for SHT
142  dctPlan = fftw_plan_r2r_1d ( static_cast<int> ( band * 2 ),
143  scratchpad,
144  scratchpad + static_cast<int> ( band * 2 ),
145  FFTW_REDFT10,
146  FFTW_ESTIMATE ) ;
147 
148  //================================================ Done
149  return ;
150 
151 }

◆ initialSplitDiscreteTransform()

void ProSHADE_internal_sphericalHarmonics::initialSplitDiscreteTransform ( proshade_unsign  oneDim,
proshade_double *&  inputReal,
proshade_double *&  inputImag,
proshade_double *&  rres,
proshade_double *&  ires,
proshade_double *  mappedData,
fftw_plan &  fftPlan,
proshade_double  normCoeff 
)

This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer.

This function takes the already initialised and prepared pointers and values and proceeds to load the data into the proper places and compute the split discrete Fourier transform, thus preparing for the spherical transform to be done.

Parameters
[in]oneDimThis is the size of any dimension of the transform (2 * bandwidth).
[in]inputRealPointer to array which should be subjected to the transform (real part).
[in]inputRealPointer to array which should be subjected to the transform (imaginary part).
[in]rresPointer to array where the transform results will be saved (real part).
[in]iresPointer to array where the transform results will be saved (imaginary part).
[in]mappedDataPointer to the data which should be decomposed.
[in]fftPlanThe prepared plan which states how the transform will be done as set by the initialiseFFTWPlans() function.
[in]normCoeffThe transform normalisation factor.

Definition at line 257 of file ProSHADE_sphericalHarmonics.cpp.

258 {
259  //================================================ Load mapped data to decomposition array
260  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( oneDim * oneDim ); iter++ )
261  {
262  inputReal[iter] = mappedData[iter];
263  inputImag[iter] = 0.0;
264  }
265 
266  //================================================ Execute fft plan along phi
267  fftw_execute_split_dft ( fftPlan, inputReal, inputImag, rres, ires ) ;
268 
269  //================================================ Normalize
270  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( oneDim * oneDim ); iter++ )
271  {
272  rres[iter] *= normCoeff;
273  ires[iter] *= normCoeff;
274  }
275 
276  //================================================ Done
277  return ;
278 
279 }

◆ placeWithinWorkspacePointers()

void ProSHADE_internal_sphericalHarmonics::placeWithinWorkspacePointers ( fftw_complex *&  workspace,
proshade_unsign  oDim,
proshade_double *&  rres,
proshade_double *&  ires,
proshade_double *&  fltres,
proshade_double *&  scratchpad 
)

This function takes the workspace pointer and correctly places the other internal pointers.

This is a simple helper function, which places the internal workspace pointers to the correct addresses as required by the SOFT2.0 dependency.

Parameters
[in]workspacePointer to the allocated workspace, within which the other pointers will be placed.
[in]oDimThe size of the single transform dimension (twice the transform bandwidth).
[in]rresPointer to where the real part of the results will be temporarily saved.
[in]iresPointer to where the imaginary part of the results will be temporarily saved.
[in]fltresPointer to where the temporary bandwise results will be saved.
[in]scratchpadPointer to extra space which is internally used (but not allocated) by SOFT2.0.

Definition at line 86 of file ProSHADE_sphericalHarmonics.cpp.

87 {
88  //================================================ Place pointers as required by SOFT2.0
89  rres = reinterpret_cast<proshade_double*> ( workspace );
90  ires = rres + ( oDim * oDim );
91  fltres = ires + ( oDim * oDim );
92  scratchpad = fltres + ( oDim / 2 );
93 
94  //================================================ Done
95  return ;
96 
97 }

◆ releaseSphericalMemory()

void ProSHADE_internal_sphericalHarmonics::releaseSphericalMemory ( proshade_double *&  inputReal,
proshade_double *&  inputImag,
proshade_double *&  outputReal,
proshade_double *&  outputImag,
double *&  tableSpaceHelper,
double **&  tableSpace,
double *&  shWeights,
fftw_complex *&  workspace,
fftw_plan &  fftPlan,
fftw_plan &  dctPlan 
)

This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer.

This function takes all the pointers required for the spherical harmonics computation and releases all the memory.

Parameters
[in]inputRealpointer to the array that contained the input real values to be freed.
[in]inputImagpointer to the array that contained the input imaginary values to be freed.
[in]outputRealpointer to the array that contained the output real values to be freed.
[in]outputImagpointer to the array that contained the output imaginary values to be freed.
[in]tableSpaceHelperpointer to the helper array for Legendre polynomials values to be freed.
[in]tableSpacepointer to the array of Legendre polynomials values to be freed.
[in]shWeightspointer to the array spherical harmonics weighhts to be freed.
[in]workspacepointer to the array for miscellaneous temporary results to be freed.
[in]fftPlanpointer to the variable where the Fourier transform was done to be freed.
[in]dctPlanpointer to the variable where the 1D r2r Fourier transform was done to be freed.

Definition at line 170 of file ProSHADE_sphericalHarmonics.cpp.

171 {
172  //================================================ Release all memory related to SH
173  delete[] inputReal;
174  delete[] inputImag;
175  delete[] outputReal;
176  delete[] outputImag;
177  delete[] tableSpaceHelper;
178  delete[] shWeights;
179  delete[] workspace;
180 
181  //================================================ Set pointers to NULL
182  tableSpaceHelper = nullptr;
183  tableSpace = nullptr;
184  shWeights = nullptr;
185  workspace = nullptr;
186 
187  //================================================ Delete fftw plans
188  fftw_destroy_plan ( dctPlan );
189  fftw_destroy_plan ( fftPlan );
190 
191  //================================================ Done
192  return ;
193 
194 }
ProSHADE_internal_sphericalHarmonics::initialiseAllMemory
void initialiseAllMemory(proshade_unsign band, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&outputReal, proshade_double *&outputImag, double *&shWeights, double **&tableSpace, double *&tableSpaceHelper, fftw_complex *&workspace, proshade_double *&rres, proshade_double *&ires, proshade_double *&fltres, proshade_double *&scratchpad, fftw_plan &fftPlan, fftw_plan &dctPlan)
This function initialises all the memory required for spherical harmonics computation.
Definition: ProSHADE_sphericalHarmonics.cpp:217
ProSHADE_internal_sphericalHarmonics::allocateComputationMemory
void allocateComputationMemory(proshade_unsign band, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&outputReal, proshade_double *&outputImag, double *&shWeights, double *&tableSpaceHelper, fftw_complex *&workspace)
This function determines the integration order for the between spheres integration.
Definition: ProSHADE_sphericalHarmonics.cpp:39
ProSHADE_internal_sphericalHarmonics::releaseSphericalMemory
void releaseSphericalMemory(proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&outputReal, proshade_double *&outputImag, double *&tableSpaceHelper, double **&tableSpace, double *&shWeights, fftw_complex *&workspace, fftw_plan &fftPlan, fftw_plan &dctPlan)
This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer.
Definition: ProSHADE_sphericalHarmonics.cpp:170
ProSHADE_internal_sphericalHarmonics::initialSplitDiscreteTransform
void initialSplitDiscreteTransform(proshade_unsign oneDim, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&rres, proshade_double *&ires, proshade_double *mappedData, fftw_plan &fftPlan, proshade_double normCoeff)
This function computes the spherical harmonics of a aingle shell, saving them in supplied pointer.
Definition: ProSHADE_sphericalHarmonics.cpp:257
ProSHADE_internal_sphericalHarmonics::computeSphericalTransformCoeffs
void computeSphericalTransformCoeffs(proshade_unsign band, proshade_double *&rdataptr, proshade_double *&idataptr, proshade_double *&outputReal, proshade_double *&outputImag, proshade_double *&rres, proshade_double *&ires, proshade_double *&fltres, proshade_double *&scratchpad, double **&tablePml, double *&shWeights, fftw_plan &dctPlan)
This function takes the split discrete transform and proceeds to complete the spherical harmonics dec...
Definition: ProSHADE_sphericalHarmonics.cpp:300
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:67
ProSHADE_internal_sphericalHarmonics::placeWithinWorkspacePointers
void placeWithinWorkspacePointers(fftw_complex *&workspace, proshade_unsign oDim, proshade_double *&rres, proshade_double *&ires, proshade_double *&fltres, proshade_double *&scratchpad)
This function takes the workspace pointer and correctly places the other internal pointers.
Definition: ProSHADE_sphericalHarmonics.cpp:86
ProSHADE_internal_sphericalHarmonics::applyCondonShortleyPhase
void applyCondonShortleyPhase(proshade_unsign band, proshade_double *outputReal, proshade_double *outputImag, proshade_complex *&shArray)
This is the final step in computing the full spherical harmonics decomposition of the input data.
Definition: ProSHADE_sphericalHarmonics.cpp:352
ProSHADE_internal_sphericalHarmonics::initialiseFFTWPlans
void initialiseFFTWPlans(proshade_unsign band, fftw_plan &fftPlan, fftw_plan &dctPlan, proshade_double *&inputReal, proshade_double *&inputImag, proshade_double *&rres, proshade_double *&ires, proshade_double *&scratchpad)
This function initialises the FFTW plans.
Definition: ProSHADE_sphericalHarmonics.cpp:113