ProSHADE  0.7.5.1 (JAN 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  tableSpaceHelper = new proshade_double [static_cast<proshade_unsign> ( Reduced_Naive_TableSize ( band, band ) +
55  Reduced_SpharmonicTableSize ( band, band ) )];
56  workspace = new fftw_complex [( 8 * band * band ) + ( 10 * band )];
57 
58  //================================================ Check memory allocation success
59  ProSHADE_internal_misc::checkMemoryAllocation ( inputReal, __FILE__, __LINE__, __func__ );
60  ProSHADE_internal_misc::checkMemoryAllocation ( inputImag, __FILE__, __LINE__, __func__ );
61  ProSHADE_internal_misc::checkMemoryAllocation ( outputReal, __FILE__, __LINE__, __func__ );
62  ProSHADE_internal_misc::checkMemoryAllocation ( outputImag, __FILE__, __LINE__, __func__ );
63  ProSHADE_internal_misc::checkMemoryAllocation ( shWeights, __FILE__, __LINE__, __func__ );
64  ProSHADE_internal_misc::checkMemoryAllocation ( tableSpaceHelper, __FILE__, __LINE__, __func__ );
65  ProSHADE_internal_misc::checkMemoryAllocation ( workspace, __FILE__, __LINE__, __func__ );
66 
67  //================================================ Done
68  return ;
69 
70 }

◆ 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 351 of file ProSHADE_sphericalHarmonics.cpp.

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

◆ 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 393 of file ProSHADE_sphericalHarmonics.cpp.

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

◆ 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 299 of file ProSHADE_sphericalHarmonics.cpp.

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

◆ 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 216 of file ProSHADE_sphericalHarmonics.cpp.

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

◆ 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 111 of file ProSHADE_sphericalHarmonics.cpp.

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

◆ 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 256 of file ProSHADE_sphericalHarmonics.cpp.

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

◆ 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 84 of file ProSHADE_sphericalHarmonics.cpp.

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

◆ 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 168 of file ProSHADE_sphericalHarmonics.cpp.

169 {
170  //================================================ Release all memory related to SH
171  delete[] inputReal;
172  delete[] inputImag;
173  delete[] outputReal;
174  delete[] outputImag;
175 
176  delete[] tableSpaceHelper;
177  delete[] tableSpace;
178  delete[] shWeights;
179 
180  fftw_free ( workspace );
181 
182  tableSpaceHelper = NULL;
183  tableSpace = NULL;
184  shWeights = NULL;
185  workspace = NULL;
186 
187  fftw_destroy_plan ( dctPlan );
188  fftw_destroy_plan ( fftPlan );
189 
190  //================================================ Done
191  return ;
192 
193 }
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:216
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:168
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:256
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:299
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:65
ProSHADE_internal_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:84
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:351
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:111