ProSHADE  0.7.6.6 (JUL 2022)
Protein Shape Detection
pyProSHADE_data.cpp
Go to the documentation of this file.
1 
22 //==================================================== Include PyBind11 header
23 #include <pybind11/pybind11.h>
24 #include <pybind11/stl.h>
25 #include <pybind11/numpy.h>
26 
27 //==================================================== Add the ProSHADE_settings and ProSHADE_run classes to the PyBind11 module
28 void add_dataClass ( pybind11::module& pyProSHADE )
29 {
30  //================================================ Export the ProSHADE_settings class
31  pybind11::class_ < ProSHADE_internal_data::ProSHADE_data > ( pyProSHADE, "ProSHADE_data" )
32 
33  //============================================ Constructors (destructors do not need wrappers???)
34  .def ( pybind11::init ( ) )
35  .def ( pybind11::init ( [] ( std::string strName, pybind11::array_t < float, pybind11::array::c_style | pybind11::array::forcecast > mapData, proshade_single xDmSz, proshade_single yDmSz, proshade_single zDmSz, proshade_unsign xDmInd, proshade_unsign yDmInd, proshade_unsign zDmInd, proshade_signed xFr, proshade_signed yFr, proshade_signed zFr, proshade_signed xT, proshade_signed yT, proshade_signed zT, proshade_unsign inputO )
36  {
37  //== Find the array size
38  pybind11::buffer_info buf = mapData.request();
39  proshade_unsign len = static_cast< proshade_unsign > ( buf.size );
40 
41  //== Allocate memory
42  double* npVals = new double[static_cast<proshade_unsign> ( len )];
43  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
44 
45  //== Copy from numpy to ProSHADE (I do not want to be directly changing the python memory, that screams trouble)
46  if ( buf.ndim == 1 )
47  {
48  for ( proshade_unsign iter = 0; iter < len; iter++ ) { npVals[iter] = static_cast < double > ( mapData.at(iter) ); }
49  }
50  else if ( buf.ndim == 3 )
51  {
52  float* dataPtr = reinterpret_cast < float* > ( buf.ptr );
53  for ( proshade_unsign xIt = 0; xIt < static_cast<proshade_unsign> ( buf.shape.at(0) ); xIt++ )
54  {
55  for ( proshade_unsign yIt = 0; yIt < static_cast<proshade_unsign> ( buf.shape.at(1) ); yIt++ )
56  {
57  for ( proshade_unsign zIt = 0; zIt < static_cast<proshade_unsign> ( buf.shape.at(2) ); zIt++ )
58  {
59  npVals[zIt + static_cast< proshade_unsign > ( buf.shape.at(2) ) * ( yIt + static_cast< proshade_unsign > ( buf.shape.at(1) ) * xIt )] = static_cast < double > ( dataPtr[zIt + static_cast< proshade_unsign > ( buf.shape.at(2) ) * ( yIt + static_cast< proshade_unsign > ( buf.shape.at(1) ) * xIt )] );
60  }
61  }
62  }
63  }
64  else
65  {
66  std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The ProSHADE_data class constructor ( ProSHADE_settings, str, numpy.ndarray, float, float, float, ... ) only supports the third argument input array in the 1D or 3D numpy.ndarray format. The supplied array has " << buf.ndim << " dims. Terminating..." << std::endl;
67  exit ( EXIT_FAILURE );
68  }
69 
70  //== Call the ProSHADE_data constructor
71  return new ProSHADE_internal_data::ProSHADE_data ( strName,
72  npVals,
73  static_cast<int> ( len ),
74  xDmSz,
75  yDmSz,
76  zDmSz,
77  xDmInd,
78  yDmInd,
79  zDmInd,
80  xFr,
81  yFr,
82  zFr,
83  xT,
84  yT,
85  zT,
86  inputO );
87  } ) )
88 
89  //============================================ Data I/O functions
90  .def ( "readInStructure",
91  [] ( ProSHADE_internal_data::ProSHADE_data &self, std::string fName, proshade_unsign inputO, ProSHADE_settings* settings, pybind11::array_t < proshade_double > maskArr, pybind11::array_t < proshade_double > weightsArr )
92  {
93  //== Sanity check
94  pybind11::buffer_info maskArr_buf = maskArr.request();
95  pybind11::buffer_info weightsArr_buf = weightsArr.request();
96 
97  //== Check for mask
98  if ( ( maskArr_buf.shape.at(0) != 0 ) && ( maskArr_buf.shape.at(1) != 0 ) && ( maskArr_buf.shape.at(2) != 0 ) )
99  {
100  //== Is number of dimensions correct?
101  if ( maskArr_buf.ndim != 3 ) { std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The fourth argument to readInStructure() must be a 3D numpy array or empty and the dimensions must match!" << std::endl; exit ( EXIT_FAILURE ); }
102 
103  //== Mask was given! Create the array
104  proshade_double* mskArr = new proshade_double[maskArr_buf.shape.at(0)*maskArr_buf.shape.at(1)*maskArr_buf.shape.at(2)];
105  ProSHADE_internal_misc::checkMemoryAllocation ( mskArr, __FILE__, __LINE__, __func__ );
106 
107  //== Copy values
108  proshade_double* arrStart = static_cast< proshade_double* > ( maskArr_buf.ptr );
109  for ( size_t iter = 0; iter < static_cast< size_t > ( maskArr_buf.shape.at(0)*maskArr_buf.shape.at(1)*maskArr_buf.shape.at(2) ); iter++ ) { mskArr[iter] = arrStart[iter]; }
110 
111  //== Check for weights
112  if ( ( weightsArr_buf.shape.at(0) != 0 ) && ( weightsArr_buf.shape.at(1) != 0 ) && ( weightsArr_buf.shape.at(2) != 0 ) )
113  {
114  //== Is number of dimensions correct?
115  if ( weightsArr_buf.ndim != 3 ) { std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The fifth argument to readInStructure() must be a 3D numpy array or empty and the dimensions must match!" << std::endl; exit ( EXIT_FAILURE ); }
116 
117  //== Weights were given! Create array
118  proshade_double* wghArr = new proshade_double[weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2)];
119  ProSHADE_internal_misc::checkMemoryAllocation ( wghArr, __FILE__, __LINE__, __func__ );
120 
121  //== Copy values
122  proshade_double* arr2Start = static_cast< proshade_double* > ( weightsArr_buf.ptr );
123  for ( size_t iter = 0; iter < static_cast< size_t > ( weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2) ); iter++ ) { wghArr[iter] = arr2Start[iter]; }
124 
125  //== Call C++ function with both mask and weights
126  self.readInStructure ( fName, inputO, settings, mskArr, static_cast< proshade_unsign > ( maskArr_buf.shape.at(0) ), static_cast< proshade_unsign > ( maskArr_buf.shape.at(1) ), static_cast< proshade_unsign > ( maskArr_buf.shape.at(2) ), wghArr, static_cast< proshade_unsign > ( weightsArr_buf.shape.at(0) ), static_cast< proshade_unsign > ( weightsArr_buf.shape.at(1) ), static_cast< proshade_unsign > ( weightsArr_buf.shape.at(2) ) );
127 
128  //== Release memory
129  delete[] wghArr;
130  }
131  else
132  {
133  //== No weights given. Call C++ function with only mask
134  self.readInStructure ( fName, inputO, settings, mskArr, static_cast< proshade_unsign > ( maskArr_buf.shape.at(0) ), static_cast< proshade_unsign > ( maskArr_buf.shape.at(1) ), static_cast< proshade_unsign > ( maskArr_buf.shape.at(2) ), nullptr, 0, 0, 0 );
135  }
136 
137  //== Release memory
138  delete[] mskArr;
139  }
140  else
141  {
142  //== Mask was empty! Check for weights
143  if ( ( weightsArr_buf.shape.at(0) != 0 ) && ( weightsArr_buf.shape.at(1) != 0 ) && ( weightsArr_buf.shape.at(2) != 0 ) )
144  {
145  //== Is number of dimensions correct?
146  if ( weightsArr_buf.ndim != 3 ) { std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The fifth argument to readInStructure() must be a 3D numpy array or empty and the dimensions must match!" << std::endl; exit ( EXIT_FAILURE ); }
147 
148  //== Weights were given! Create array
149  proshade_double* wghArr = new proshade_double[weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2)];
150  ProSHADE_internal_misc::checkMemoryAllocation ( wghArr, __FILE__, __LINE__, __func__ );
151 
152  //== Copy values
153  proshade_double* arr2Start = static_cast< proshade_double* > ( weightsArr_buf.ptr );
154  for ( size_t iter = 0; iter < static_cast< size_t > ( weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2) ); iter++ ) { wghArr[iter] = arr2Start[iter]; }
155 
156  //== Call C++ function with only weights
157  self.readInStructure ( fName, inputO, settings, nullptr, 0, 0, 0, wghArr, static_cast< proshade_unsign > ( weightsArr_buf.shape.at(0) ), static_cast< proshade_unsign > ( weightsArr_buf.shape.at(1) ), static_cast< proshade_unsign > ( weightsArr_buf.shape.at(2) ) );
158 
159  //== Release memory
160  delete[] wghArr;
161  }
162  else
163  {
164  //== No weights either. Call C++ function with no mask and no weights
165  self.readInStructure ( fName, inputO, settings, nullptr, 0, 0, 0, nullptr, 0, 0, 0 );
166  }
167  }
168 
169  //== Done
170  return ;
171  }, "This function returns the group elements as rotation matrices of any point group described by the detected axes.", pybind11::arg ( "fName" ), pybind11::arg ( "inputO" ), pybind11::arg ( "settings" ), pybind11::arg( "maskArr" ) = pybind11::array_t < proshade_double > (), pybind11::arg( "weightsArr" ) = pybind11::array_t < proshade_double > () )
172  .def ( "writeMap", &ProSHADE_internal_data::ProSHADE_data::writeMap, "Function for writing out the internal structure representation in MRC MAP format.", pybind11::arg ( "fname" ), pybind11::arg ( "title" ) = "Created by ProSHADE and written by GEMMI", pybind11::arg ( "mode" ) = 2 )
173  .def ( "writePdb", &ProSHADE_internal_data::ProSHADE_data::writePdb, "This function writes out the co-ordinates file with ProSHADE type rotation and translation applied.", pybind11::arg ( "fname" ), pybind11::arg ( "euA" ) = 0.0, pybind11::arg ( "euB" ) = 0.0, pybind11::arg ( "euG" ) = 0.0, pybind11::arg ( "trsX" ) = 0.0, pybind11::arg ( "trsY" ) = 0.0, pybind11::arg ( "trsZ" ) = 0.0, pybind11::arg ( "rotX" ) = 0.0, pybind11::arg ( "rotY" ) = 0.0, pybind11::arg ( "rotZ" ) = 0.0, pybind11::arg ( "firstModel" ) = true )
174  .def ( "writeGemmi", &ProSHADE_internal_data::ProSHADE_data::writeGemmi, "This function writes out the gemmi::Structure object with ProSHADE type rotation and translation applied.", pybind11::arg ( "fname" ), pybind11::arg ( "gemmiStruct" ), pybind11::arg ( "euA" ) = 0.0, pybind11::arg ( "euB" ) = 0.0, pybind11::arg ( "euG" ) = 0.0, pybind11::arg ( "trsX" ) = 0.0, pybind11::arg ( "trsY" ) = 0.0, pybind11::arg ( "trsZ" ) = 0.0, pybind11::arg ( "rotX" ) = 0.0, pybind11::arg ( "rotY" ) = 0.0, pybind11::arg ( "rotZ" ) = 0.0, pybind11::arg ( "firstModel" ) = true )
175  .def ( "getMap",
176  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < proshade_double >
177  {
178  //== Copy the values
179  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { self.xDimIndices, self.yDimIndices, self.zDimIndices }, // Shape
180  { self.yDimIndices * self.zDimIndices * sizeof(proshade_double),
181  self.zDimIndices * sizeof(proshade_double),
182  sizeof(proshade_double) }, // C-stype strides
183  self.internalMap ); // Data
184 
185  //== Done
186  return ( retArr );
187  } )
188 
189  //============================================ Data processing functions
190  .def ( "processInternalMap", &ProSHADE_internal_data::ProSHADE_data::processInternalMap, "This function simply clusters several map manipulating functions which should be called together. These include centering, phase removal, normalisation, adding extra space, etc.", pybind11::arg ( "settings" ) )
191  .def ( "invertMirrorMap", &ProSHADE_internal_data::ProSHADE_data::invertMirrorMap, "Function for inverting the map to its mirror image.", pybind11::arg ( "settings" ) )
192  .def ( "normaliseMap", &ProSHADE_internal_data::ProSHADE_data::normaliseMap, "Function for normalising the map values to mean 0 and sd 1.", pybind11::arg ( "settings" ) )
193  .def ( "maskMap", &ProSHADE_internal_data::ProSHADE_data::maskMap, "Function for computing the map mask using blurring and X IQRs from median.", pybind11::arg ( "settings" ) )
194  .def ( "reSampleMap", &ProSHADE_internal_data::ProSHADE_data::reSampleMap, "This function changes the internal map sampling to conform to particular resolution value.", pybind11::arg ( "settings" ) )
195  .def ( "centreMapOnCOM", &ProSHADE_internal_data::ProSHADE_data::centreMapOnCOM, "This function shits the map so that its COM is in the centre of the map.", pybind11::arg ( "settings" ) )
196  .def ( "addExtraSpace", &ProSHADE_internal_data::ProSHADE_data::addExtraSpace, "This function increases the size of the map so that it can add empty space around it.", pybind11::arg ( "settings" ) )
197  .def ( "removePhaseInormation", &ProSHADE_internal_data::ProSHADE_data::removePhaseInormation, "This function removes phase from the map, effectively converting it to Patterson map.", pybind11::arg ( "settings" ) )
198  .def ( "shiftToBoxCentre" , &ProSHADE_internal_data::ProSHADE_data::shiftToBoxCentre, "This function shifts the internal map so that its centre of the box is at required position.", pybind11::arg ( "settings" ) )
199  .def ( "shiftToRotationCentre", &ProSHADE_internal_data::ProSHADE_data::shiftToRotationCentre, "This function shifts the internal map so that its rotation centre is at the centre of the box.", pybind11::arg ( "settings" ) )
200  .def ( "getReBoxBoundaries",
201  [] ( ProSHADE_internal_data::ProSHADE_data &self,ProSHADE_settings* settings ) -> pybind11::array_t < proshade_signed >
202  {
203  //== Allocate output memory
204  proshade_signed* retVals = new proshade_signed[6];
205  ProSHADE_internal_misc::checkMemoryAllocation ( retVals, __FILE__, __LINE__, __func__ );
206 
207  //== If same bounds as first one are required, test if possible and return these instead
208  if ( settings->useSameBounds && ( self.inputOrder != 0 ) )
209  {
210  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->forceBounds[iter]; }
211  }
212  //== In this case, bounds need to be found de novo
213  else
214  {
215  //== Find the non-zero bounds
217  static_cast< proshade_signed > ( self.xDimIndices ),
218  static_cast< proshade_signed > ( self.yDimIndices ),
219  static_cast< proshade_signed > ( self.zDimIndices ),
220  retVals );
221 
222  //== Add the extra space
223  ProSHADE_internal_mapManip::addExtraBoundSpace ( self.xDimIndices, self.yDimIndices, self.zDimIndices,
224  self.xDimSize, self.yDimSize, self.zDimSize, retVals, settings->boundsExtraSpace );
225 
226  //== Beautify boundaries
227  ProSHADE_internal_mapManip::beautifyBoundaries ( retVals, self.xDimIndices, self.yDimIndices, self.zDimIndices, settings->boundsSimilarityThreshold );
228 
229  //== If need be, save boundaries to be used for all other structure
230  if ( settings->useSameBounds && ( self.inputOrder == 0 ) ) { for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->forceBounds[iter] = retVals[iter]; } }
231  }
232 
233  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
234  pybind11::capsule pyCapsuleReBox2 ( retVals, []( void *f ) { proshade_signed* foo = reinterpret_cast< proshade_signed* > ( f ); delete foo; } );
235 
236  //== Copy the value
237  pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 }, // Shape
238  { sizeof(proshade_signed) }, // C-stype strides
239  retVals, // Data
240  pyCapsuleReBox2 ); // Capsule
241 
242  //== Done
243  return ( retArr );
244  }, "This function finds the boundaries enclosing positive map values and adds some extra space." )
245  .def ( "createNewMapFromBounds",
246  [] ( ProSHADE_internal_data::ProSHADE_data &self, pybind11::array_t < proshade_signed > bounds, ProSHADE_internal_data::ProSHADE_data* newStr, ProSHADE_settings* settings )
247  {
248  //== Allocate memory for bounds copy
249  proshade_signed* newBounds = new proshade_signed[6];
250  ProSHADE_internal_misc::checkMemoryAllocation ( newBounds, __FILE__, __LINE__, __func__ );
251 
252  //== Copy to the pointer
253  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
254 
255  //== Run C++ function
256  self.createNewMapFromBounds ( settings, newStr, newBounds );
257 
258  //== Done
259  return ;
260  }, "This function creates a new structure from the calling structure and new bounds values." )
261 
262  //============================================ Data sphere mapping functions
263  .def ( "mapToSpheres", &ProSHADE_internal_data::ProSHADE_data::mapToSpheres, "This function converts the internal map onto a set of concentric spheres.", pybind11::arg ( "settings" ) )
264  .def ( "getSpherePositions", &ProSHADE_internal_data::ProSHADE_data::getSpherePositions, "This function determines the sphere positions (radii) for sphere mapping.", pybind11::arg ( "settings" ) )
265  .def ( "computeSphericalHarmonics", &ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics, "This function computes the spherical harmonics decomposition for the whole structure.", pybind11::arg ( "settings" ) )
266 
267  //============================================ Accessor functions
268  .def ( "getXDimSize", &ProSHADE_internal_data::ProSHADE_data::getXDimSize, "This function allows access to the map size in angstroms along the X axis." )
269  .def ( "getYDimSize", &ProSHADE_internal_data::ProSHADE_data::getYDimSize, "This function allows access to the map size in angstroms along the Y axis." )
270  .def ( "getZDimSize", &ProSHADE_internal_data::ProSHADE_data::getZDimSize, "This function allows access to the map size in angstroms along the Z axis." )
271  .def ( "getXDim", &ProSHADE_internal_data::ProSHADE_data::getXDim, "This function allows access to the map size in indices along the X axis." )
272  .def ( "getYDim", &ProSHADE_internal_data::ProSHADE_data::getYDim, "This function allows access to the map size in indices along the Y axis." )
273  .def ( "getZDim", &ProSHADE_internal_data::ProSHADE_data::getZDim, "This function allows access to the map size in indices along the Z axis." )
274 
275  //============================================ Symmetry related functions
276  .def ( "computeRotationFunction", &ProSHADE_internal_data::ProSHADE_data::computeRotationFunction, "This function computes the self-rotation function for this structure and stores it internally in the ProSHADE_data object.", pybind11::arg ( "settings" ) )
277  .def ( "detectSymmetryInStructure",
279  {
280  //== Call the appropriate C++ function
281  self.detectSymmetryFromAngleAxisSpace ( settings );
282  }, "This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg ( "settings" ) )
283  .def ( "reRunSymmetryDetectionThreshold",
284  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings, proshade_double threshold )
285  {
286  //== Prepare FSC computation memory and variables
287  fftw_complex* fCoeffsCut;
288  proshade_double **bindata, *fscByBin;
289  proshade_signed *cutIndices, *binCounts, noBins, cutXDim, cutYDim, cutZDim;
290  self.prepareFSCFourierMemory ( cutIndices, fCoeffsCut, &noBins, bindata, binCounts, fscByBin, settings->requestedResolution, &cutXDim, &cutYDim, &cutZDim );
291 
292  //== Run detection with threshold
293  self.determineRecommendedSymmetry ( settings, threshold, cutIndices, fCoeffsCut, noBins, bindata, binCounts, fscByBin, cutXDim, cutYDim, cutZDim );
294 
295  //== Release FSC computation memory
296  for (size_t binIt = 0; binIt < static_cast< size_t > ( noBins ); binIt++ ) { delete[] bindata[binIt]; }
297  delete[] bindata;
298  delete[] binCounts;
299  delete[] cutIndices;
300  fftw_free ( fCoeffsCut );
301  }, "This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg ( "settings" ), pybind11::arg ( "threshold" ) )
302  .def ( "getRecommendedSymmetryType", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryType, "This function simply returns the detected recommended symmetry type." )
303  .def ( "getRecommendedSymmetryFold", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryFold, "This function simply returns the detected recommended symmetry fold." )
304  .def ( "getRecommendedSymmetryAxes",
305  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < float >
306  {
307  //== Allocate memory for the numpy values
308  float* npVals = new float[static_cast<unsigned int> ( self.recommendedSymmetryValues.size() * 7 )];
309  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
310 
311  //== Copy values
312  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( self.recommendedSymmetryValues.size() ); iter++ ) { for ( proshade_unsign it = 0; it < 7; it++ ) { npVals[(iter*7)+it] = static_cast< float > ( self.recommendedSymmetryValues.at(iter)[it] ); } }
313 
314  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
315  pybind11::capsule pyCapsuleStrRecSym ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
316 
317  //== Copy the value
318  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( self.recommendedSymmetryValues.size() ), static_cast<int> ( 7 ) }, // Shape
319  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
320  npVals, // Data
321  pyCapsuleStrRecSym ); // Capsule
322 
323  //== Done
324  return ( retArr );
325  }, "This function returns the recommended symmetry axes as a 2D numpy array." )
326  .def ( "getAllCSyms",
327  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < float >
328  {
329  //== Allocate memory for the numpy values
330  float* npVals = new float[static_cast<unsigned int> ( self.cyclicSymmetries.size() * 7 )];
331  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
332 
333  //== Copy values
334  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( self.cyclicSymmetries.size() ); iter++ ) { for ( proshade_unsign it = 0; it < 7; it++ ) { npVals[(iter*7)+it] = static_cast< float > ( self.cyclicSymmetries.at(iter)[it] ); } }
335 
336  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
337  pybind11::capsule pyCapsuleStrSymList ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
338 
339  //== Copy the value
340  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( self.cyclicSymmetries.size() ), 7 }, // Shape
341  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
342  npVals, // Data
343  pyCapsuleStrSymList ); // Capsule
344 
345  //== Done
346  return ( retArr );
347  }, "This function returns all symmetry axes as a 2D numpy array." )
348  .def ( "getNonCSymmetryAxesIndices",
349  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::dict
350  {
351  //== Initialise local variables
352  pybind11::dict retDict;
353  pybind11::list dList, tList, oList, iList;
354 
355  //== Fill in the D list
356  for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( self.dihedralSymmetries.size() ); dIt++ )
357  {
358  pybind11::list memList;
359  for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( self.dihedralSymmetries.at(dIt).size() ); memIt++ )
360  {
361  pybind11::list axList;
362  for ( size_t vIt = 0; vIt < 7; vIt++ ) { axList.append ( self.dihedralSymmetries.at(dIt).at(memIt)[vIt] ); }
363  memList.append ( axList );
364  }
365  dList.append ( memList );
366  }
367 
368  //== Fill in the T list
369  for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( self.tetrahedralSymmetries.size() ); tIt++ )
370  {
371  pybind11::list axList;
372  for ( size_t vIt = 0; vIt < 7; vIt++ ) { axList.append ( self.tetrahedralSymmetries.at(tIt)[vIt] ); }
373  tList.append ( axList );
374  }
375 
376  //== Fill in the O list
377  for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( self.octahedralSymmetries.size() ); oIt++ )
378  {
379  pybind11::list axList;
380  for ( size_t vIt = 0; vIt < 7; vIt++ ) { axList.append ( self.octahedralSymmetries.at(oIt)[vIt] ); }
381  oList.append ( axList );
382  }
383 
384  //== Fill in the T list
385  for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( self.icosahedralSymmetries.size() ); iIt++ )
386  {
387  pybind11::list axList;
388  for ( size_t vIt = 0; vIt < 7; vIt++ ) { axList.append ( self.icosahedralSymmetries.at(iIt)[vIt] ); }
389  iList.append ( axList );
390  }
391 
392  //== Save the lists to the dict
393  retDict[ pybind11::handle ( pybind11::str ( "D" ).ptr ( ) ) ] = dList;
394  retDict[ pybind11::handle ( pybind11::str ( "T" ).ptr ( ) ) ] = tList;
395  retDict[ pybind11::handle ( pybind11::str ( "O" ).ptr ( ) ) ] = oList;
396  retDict[ pybind11::handle ( pybind11::str ( "I" ).ptr ( ) ) ] = iList;
397 
398  //== Done
399  return ( retDict );
400  }, "This function returns array of non-C axes indices." )
401  .def ( "getAllGroupElements",
402  [] ( ProSHADE_internal_data::ProSHADE_data &self, pybind11::array_t < int > axList, std::string groupType, proshade_double matrixTolerance ) -> pybind11::list
403  {
404  //== Sanity check
405  pybind11::buffer_info buf = axList.request();
406  if ( buf.ndim != 1 ) { std::cerr << "!!! ProSHADE PYTHON MODULE ERROR !!! The second argument to getAllGroupElements() must be a 1D numpy array stating the indices of the axes forming the group!" << std::endl; exit ( EXIT_FAILURE ); }
407 
408  //== Convert to vector of unsigns
409  std::vector< proshade_unsign > axesList;
410  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( axList.size() ); iter++ ) { ProSHADE_internal_misc::addToUnsignVector ( &axesList, static_cast< proshade_unsign > ( axList.at(iter) ) ); }
411 
412  //== Get the results
413  std::vector < std::vector< proshade_double > > vals = self.getAllGroupElements ( axesList, groupType, matrixTolerance );
414 
415  //== Initialise return variable
416  pybind11::list retList;
417 
418  //== Copy values
419  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
420  {
421  //== Allocate memory for the numpy values
422  float* npVals = new float[static_cast<unsigned int> ( 9 )];
423  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
424 
425  //== Copy values to memory
426  for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] = static_cast< float > ( vals.at(iter).at(it) ); }
427 
428  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
429  pybind11::capsule pyCapsuleGrpEl ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
430 
431  //== Copy the value
432  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 }, // Shape
433  { 3 * sizeof(float), sizeof(float) }, // C-stype strides
434  npVals, // Data
435  pyCapsuleGrpEl ); // Capsule
436 
437  //== Save matrix
438  retList.append ( retArr );
439  }
440 
441  //== Done
442  return ( retList );
443  }, "This function returns the group elements as rotation matrices of any point group described by the detected axes.", pybind11::arg ( "axList" ), pybind11::arg ( "groupType" ) = "", pybind11::arg( "matrixTolerance" ) = 0.05 )
444 
445  .def ( "getMapCOMProcessChange",
446  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < float >
447  {
448  //== Get the values
449  std::vector< proshade_double > vals = self.getMapCOMProcessChange ();
450 
451  //== Allocate memory for the numpy values
452  float* npVals = new float[static_cast<unsigned int> ( 3 )];
453  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
454 
455  //== Copy values
456  for ( proshade_unsign iter = 0; iter < 3; iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
457 
458  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
459  pybind11::capsule pyCapsuleSymShiftDat ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
460 
461  //== Copy the value
462  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( vals.size() ) }, // Shape
463  { sizeof(float) }, // C-stype strides
464  npVals, // Data
465  pyCapsuleSymShiftDat ); // Capsule
466 
467  //== Done
468  return ( retArr );
469  }, "This function returns the shift in Angstrom applied to the internal map representation in order to align its COM with the centre of box." )
470 
471  //============================================ Overlay related functions
472  .def ( "getOverlayRotationFunction", &ProSHADE_internal_data::ProSHADE_data::getOverlayRotationFunction, "This function computes the overlay rotation function (i.e. the correlation function in SO(3) space).", pybind11::arg ( "settings" ), pybind11::arg ( "obj2" ) )
473  .def ( "getBestRotationMapPeaksEulerAngles",
474  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
475  {
476  //== Get values
477  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
478 
479  //== Allocate memory for the numpy values
480  float* npVals = new float[static_cast<unsigned int> (vals.size())];
481  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
482 
483  //== Copy values
484  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
485 
486  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
487  pybind11::capsule pyCapsuleEuPeak ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
488 
489  //== Copy the value
490  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<unsigned int> (vals.size()) }, // Shape
491  { sizeof(float) }, // C-stype strides
492  npVals, // Data
493  pyCapsuleEuPeak ); // Capsule (C++ destructor, basically)
494 
495  //== Done
496  return ( retArr );
497  }, "This function returns a vector of three floats, the three Euler angles of the best peak in the rotation map.", pybind11::arg ( "settings" ) )
498  .def ( "getBestRotationMapPeaksRotationMatrix",
499  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < proshade_double >
500  {
501  //== Get values
502  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
503 
504  //== Convert Euler ZYZ to matrix
505  proshade_double* retMat = new proshade_double[9];
506  ProSHADE_internal_misc::checkMemoryAllocation ( retMat, __FILE__, __LINE__, __func__ );
507  ProSHADE_internal_maths::getRotationMatrixFromEulerZYZAngles ( vals.at(0), vals.at(1), vals.at(2), retMat );
508 
509  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
510  pybind11::capsule pyCapsuleRMPeak ( retMat, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
511 
512  //== Copy the value
513  pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 }, // Shape
514  { 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
515  retMat, // Data
516  pyCapsuleRMPeak ); // Capsule
517 
518  //== Done
519  return ( retArr );
520  }, "This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg ( "settings" ) )
521  .def ( "rotateMapReciprocalSpace", &ProSHADE_internal_data::ProSHADE_data::rotateMapReciprocalSpace, "This function rotates a map based on the given Euler angles.", pybind11::arg ( "settings" ), pybind11::arg ( "eulerAlpha" ), pybind11::arg ( "eulerBeta" ), pybind11::arg ( "eulerGamma" ) )
522  .def ( "rotateMapRealSpaceInPlace", &ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpaceInPlace, "This function rotates a map based on the given Euler angles in real space using interpolation.", pybind11::arg ( "eulerAlpha" ), pybind11::arg ( "eulerBeta" ), pybind11::arg ( "eulerGamma" ) )
523  .def ( "zeroPaddToDims", &ProSHADE_internal_data::ProSHADE_data::zeroPaddToDims, "This function changes the size of a structure to fit the supplied new limits.", pybind11::arg ( "xDimMax" ), pybind11::arg ( "yDimMax" ), pybind11::arg ( "zDimMax" ) )
524  .def ( "computeTranslationMap", &ProSHADE_internal_data::ProSHADE_data::computeTranslationMap, "This function does the computation of the translation map and saves results internally.", pybind11::arg ( "staticStructure" ) )
525  .def ( "getOverlayTranslations",
526  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_internal_data::ProSHADE_data* staticStructure ) -> pybind11::dict
527  {
528  //== Get values
529  std::vector< proshade_double > vals = self.getBestTranslationMapPeaksAngstrom ( staticStructure );
530 
531  //== Initialise variables
532  pybind11::dict retDict;
533  pybind11::list rotCen, toOverlay;
534 
535  //== Copy values
536  rotCen.append ( self.originalPdbRotCenX );
537  rotCen.append ( self.originalPdbRotCenY );
538  rotCen.append ( self.originalPdbRotCenZ );
539 
540  toOverlay.append ( self.originalPdbTransX );
541  toOverlay.append ( self.originalPdbTransY );
542  toOverlay.append ( self.originalPdbTransZ );
543 
544  //== Save results to return dict
545  retDict[ pybind11::handle ( pybind11::str ( "centreOfRotation" ).ptr ( ) ) ] = rotCen;
546  retDict[ pybind11::handle ( pybind11::str ( "rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
547 
548  //== Done
549  return ( retDict );
550  }, "This function returns the vector from optimal rotation centre to origin and the optimal overlay translation vector. These two vectors allow overlaying the inputs (see documentation for details on how the two vectors should be used).", pybind11::arg ( "staticStructure" ) )
551  .def ( "translateMap", &ProSHADE_internal_data::ProSHADE_data::translateMap, "This function translates the map by a given number of Angstroms along the three axes. Please note the translation happens firstly to the whole map box and only the translation remainder that cannot be achieved by moving the box will be corrected for using reciprocal space translation within the box.", pybind11::arg ( "trsX" ), pybind11::arg ( "trsY" ), pybind11::arg ( "trsZ" ) )
552 
553  //============================================ Internal arrays access functions
554  .def ( "findSHIndex",
555  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed shell, proshade_signed band, proshade_signed order ) -> proshade_signed
556  {
557  //== Get value
558  proshade_signed index = seanindex ( static_cast< int > ( order ),
559  static_cast< int > ( band ),
560  static_cast< int > ( self.spheres[shell]->getLocalBandwidth() ) );
561 
562  //== Done
563  return ( index );
564  }, "This function finds the correct index for given shell, band and order in the spherical harmonics array. Please note that the order is expected in range -band <= 0 <- band and NOT from 0 to ( 2 * band ) + 1." )
565  .def ( "getSphericalHarmonics",
566  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex<proshade_double> >
567  {
568  //== Allocate memory for the numpy values
569  std::complex<proshade_double>* npVals = new std::complex<proshade_double>[static_cast<proshade_unsign> ( self.noSpheres * pow( self.maxShellBand, 2.0 ) )];
570  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
571 
572  //== Fill with zeroes as the matrix will be sparse
573  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( self.noSpheres * pow( self.maxShellBand, 2.0 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
574 
575  //== Initialise variables
576  proshade_signed pyPosSH;
577  proshade_signed pyPos;
578 
579  //== Copy data to new memory
580  for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> ( self.noSpheres ); shIt++ )
581  {
582  for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> ( self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
583  {
584  for ( proshade_signed order = -bnd; order <= bnd; order++ )
585  {
586  pyPosSH = ( static_cast< proshade_signed > ( shIt ) * static_cast< proshade_signed > ( std::pow ( self.maxShellBand, 2 ) ) );
587  pyPos = seanindex ( static_cast< int > ( order ),
588  static_cast< int > ( bnd ),
589  static_cast< int > ( self.spheres[shIt]->getLocalBandwidth() ) );
590  npVals[pyPosSH+pyPos].real ( self.sphericalHarmonics[shIt][pyPos][0] );
591  npVals[pyPosSH+pyPos].imag ( self.sphericalHarmonics[shIt][pyPos][1] );
592  }
593  }
594  }
595 
596  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
597  pybind11::capsule pyCapsuleSHs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
598 
599  //== Copy the value
600  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
601  { static_cast<int> ( self.noSpheres ), static_cast<int> ( pow ( self.maxShellBand, 2.0 ) ) },
602  { sizeof ( std::complex < proshade_double > ) * static_cast< proshade_unsign > ( std::pow ( static_cast< proshade_unsign > ( self.maxShellBand ), 2 ) ), sizeof ( std::complex < proshade_double > ) },
603  npVals,
604  pyCapsuleSHs );
605 
606  //== Done
607  return ( retArr );
608  }, "This function returns a 2D numpy array of complex numbers representing the spherical harmonics computed for the structure. The first dimension of the array is the spheres (i.e. each sphere has its own array) and the second dimension is the band and order combination as given by the findSHIndex() function. Please note that while each sphere can have different number of spherical harmonics coefficients (depending on the settings.progressiveSphereMapping value), this array uses maxEMatDim**2 values to make sure the length are equal for all spheres. To avoid issues, please use the findSHIndex() to correctly find the index of a particular shell, band, order combination." )
609  .def ( "getEMatrix",
610  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
611  {
612  //== Allocate memory for the numpy values
613  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxEMatDim * ( ( self.maxEMatDim * 2 ) + 1 ) * ( ( self.maxEMatDim * 2 ) + 1 ) )];
614  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
615 
616  //== Allocate local variables
617  proshade_signed index = 0;
618 
619  //== Fill with zeroes as the matrix will be sparse
620  for ( proshade_unsign iter = 0; iter < static_cast < proshade_unsign > ( self.maxEMatDim * ( ( self.maxEMatDim * 2 ) + 1 ) * ( ( self.maxEMatDim * 2 ) + 1 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
621 
622  //== Copy data to new memory
623  proshade_double emReal, emImag;
624  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxEMatDim ); bandIter++ )
625  {
626  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
627  {
628  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
629  {
630  index = order2 + ( ( static_cast< proshade_signed > ( self.maxEMatDim ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxEMatDim ) * 2 ) + 1 ) * bandIter );
631  self.getEMatrixValue( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( order1 ), static_cast< proshade_unsign > ( order2 ), &emReal, &emImag );
632  npVals[index].real ( emReal );
633  npVals[index].imag ( emImag );
634  }
635  }
636  }
637 
638  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
639  pybind11::capsule pyCapsuleEMs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
640 
641  //== Create the output object
642  proshade_unsign orderSize = ( ( static_cast< proshade_unsign > ( self.maxEMatDim ) * 2 ) + 1 );
643  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
644  { self.maxEMatDim, orderSize, orderSize },
645  { sizeof ( std::complex < proshade_double > ) * orderSize * orderSize,
646  sizeof ( std::complex < proshade_double > ) * orderSize,
647  sizeof ( std::complex < proshade_double > ) },
648  npVals,
649  pyCapsuleEMs );
650 
651  //== Done
652  return ( retArr );
653  }, "This function returns the weighted E matrix values (these are the integral over all spheres of ( c1^(l,m) * c2^(l,m') ) values) obtained when rotation function or self-rotation function are computed. The returned array has three dimensions, first being the band, second being the order1 and third being the order2. Please note that as negative indexing does not work, the order indexing starts from 0 - i.e. array[1][0][0] means band 1 ; order1 = -1 and order2 = -1." )
654  .def ( "getSO3Coefficients",
655  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
656  {
657  //== Allocate memory for the numpy values
658  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxEMatDim * ( ( self.maxEMatDim * 2 ) + 1 ) * ( ( self.maxEMatDim * 2 ) + 1 ) )];
659  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
660 
661  //== Allocate local variables
662  proshade_signed index = 0;
663 
664  //== Fill with zeroes as the matrix will be sparse
665  for ( proshade_unsign iter = 0; iter < static_cast < proshade_unsign > ( self.maxEMatDim * ( ( self.maxEMatDim * 2 ) + 1 ) * ( ( self.maxEMatDim * 2 ) + 1 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
666 
667  //== Copy data to new memory
668  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxEMatDim ); bandIter++ )
669  {
670  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
671  {
672  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
673  {
674  index = order2 + ( ( static_cast< proshade_signed > ( self.maxEMatDim ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxEMatDim ) * 2 ) + 1 ) * bandIter );
675  npVals[index].real ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][0] );
676  npVals[index].imag ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][1] );
677  }
678  }
679  }
680 
681  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
682  pybind11::capsule pyCapsuleSOCoeffs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
683 
684  //== Create the output object
685  proshade_unsign orderSize = ( ( static_cast< proshade_unsign > ( self.maxEMatDim ) * 2 ) + 1 );
686  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
687  { self.maxEMatDim, orderSize, orderSize },
688  { sizeof ( std::complex < proshade_double > ) * orderSize * orderSize,
689  sizeof ( std::complex < proshade_double > ) * orderSize,
690  sizeof ( std::complex < proshade_double > ) },
691  npVals,
692  pyCapsuleSOCoeffs );
693 
694  //== Done
695  return ( retArr );
696  }, "This function returns the SO(3) coefficient values (these are normalised E matrix values with sign modification) obtained when rotation function or self-rotation function are computed. The returned array has three dimensions, first being the band, second being the order1 and third being the order2. Please note that as negative indexing does not simply work, the order indexing starts from 0 - i.e. array[1][0][0] means band 1 ; order1 = -1 and order2 = -1." )
697  .def ( "getRotationFunctionMap",
698  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
699  {
700  //== Allocate memory for the numpy values
701  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( ( self.maxEMatDim * 2 ) * ( self.maxEMatDim * 2 ) * ( self.maxEMatDim * 2 ) )];
702  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
703 
704  //== Copy the values
705  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( ( self.maxEMatDim * 2 ) * ( self.maxEMatDim * 2 ) * ( self.maxEMatDim * 2 ) ); iter++ ) { npVals[iter].real ( self.so3CoeffsInverse[iter][0] ); npVals[iter].imag ( self.so3CoeffsInverse[iter][1] ); }
706 
707  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
708  pybind11::capsule pyCapsuleRotMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
709 
710  //== Create python readable object with the correct memory access
711  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
712  { ( self.maxEMatDim * 2 ), ( self.maxEMatDim * 2 ), ( self.maxEMatDim * 2 ) }, // Shape
713  { ( self.maxEMatDim * 2 ) * ( self.maxEMatDim * 2 ) * sizeof(std::complex < proshade_double >),
714  ( self.maxEMatDim * 2 ) * sizeof(std::complex < proshade_double >),
715  sizeof(std::complex < proshade_double >) }, // C-stype strides
716  npVals, // Data
717  pyCapsuleRotMap ); // Capsule (destructor)
718 
719  //== Done
720  return ( retArr );
721  }, "This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
722  .def ( "getRotationMatrixFromSOFTCoordinates",
723  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed xPos, proshade_signed yPos, proshade_signed zPos ) -> pybind11::array_t < proshade_double >
724  {
725  //== Allocate memory for the numpy values
726  proshade_double* npVals = new proshade_double[9];
727  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
728 
729  //== Initialise local variables
730  proshade_double eulA, eulB, eulG;
731 
732  //== Compute the Euler angles from SOFT position
733  ProSHADE_internal_maths::getEulerZYZFromSOFTPosition ( static_cast< proshade_signed > ( self.maxEMatDim ), xPos, yPos, zPos, &eulA, &eulB, &eulG );
734 
735  //== Compute the rotation matrix
737 
738  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
739  pybind11::capsule pyCapsuleRMSoft ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
740 
741  //== Create python readable object with the correct memory access
742  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
743  { 3 * sizeof(proshade_double), sizeof(proshade_double) },
744  npVals,
745  pyCapsuleRMSoft );
746 
747  //== Done
748  return ( retArr );
749  }, "This function converts a given rotation function map position onto the corresponding rotation matrix.", pybind11::arg ( "xPos" ), pybind11::arg ( "yPos" ), pybind11::arg ( "zPos" ) )
750  .def ( "getTranslationFunctionMap",
751  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
752  {
753  //== Allocate memory for the numpy values
754  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.getXDim() * self.getYDim() * self.getZDim() )];
755  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
756 
757  //== Copy the values
758  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( self.getXDim() * self.getYDim() * self.getZDim() ); iter++ ) { npVals[iter].real ( self.translationMap[iter][0] ); npVals[iter].imag ( self.translationMap[iter][1] ); }
759 
760  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
761  pybind11::capsule pyCapsuleTrsMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
762 
763  //== Create python readable object with the correct memory access
764  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
765  { self.getXDim(), self.getYDim(), self.getZDim() }, // Shape
766  { self.getYDim() * self.getZDim() * sizeof(std::complex < proshade_double >),
767  self.getZDim() * sizeof(std::complex < proshade_double >),
768  sizeof(std::complex < proshade_double >) }, // C-stype strides
769  npVals, // Data
770  pyCapsuleTrsMap ); // Capsule (destructor)
771 
772  //== Done
773  return ( retArr );
774  }, "This function returns the translation function as a three-dimensional map of complex numbers." )
775 
776  //============================================ Member variables
777  .def_readwrite ( "fileName", &ProSHADE_internal_data::ProSHADE_data::fileName )
778  .def_readwrite ( "xDimSize", &ProSHADE_internal_data::ProSHADE_data::xDimSize )
779  .def_readwrite ( "yDimSize", &ProSHADE_internal_data::ProSHADE_data::yDimSize )
780  .def_readwrite ( "zDimSize", &ProSHADE_internal_data::ProSHADE_data::zDimSize )
781  .def_readwrite ( "aAngle", &ProSHADE_internal_data::ProSHADE_data::aAngle )
782  .def_readwrite ( "bAngle", &ProSHADE_internal_data::ProSHADE_data::bAngle )
783  .def_readwrite ( "cAngle", &ProSHADE_internal_data::ProSHADE_data::cAngle )
784  .def_readwrite ( "xDimIndices", &ProSHADE_internal_data::ProSHADE_data::xDimIndices )
785  .def_readwrite ( "yDimIndices", &ProSHADE_internal_data::ProSHADE_data::yDimIndices )
786  .def_readwrite ( "zDimIndices", &ProSHADE_internal_data::ProSHADE_data::zDimIndices )
787  .def_readwrite ( "xGridIndices", &ProSHADE_internal_data::ProSHADE_data::xGridIndices )
788  .def_readwrite ( "yGridIndices", &ProSHADE_internal_data::ProSHADE_data::yGridIndices )
789  .def_readwrite ( "zGridIndices", &ProSHADE_internal_data::ProSHADE_data::zGridIndices )
790  .def_readwrite ( "xAxisOrder", &ProSHADE_internal_data::ProSHADE_data::xAxisOrder )
791  .def_readwrite ( "yAxisOrder", &ProSHADE_internal_data::ProSHADE_data::yAxisOrder )
792  .def_readwrite ( "zAxisOrder", &ProSHADE_internal_data::ProSHADE_data::zAxisOrder )
793  .def_readwrite ( "xAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::xAxisOrigin )
794  .def_readwrite ( "yAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::yAxisOrigin )
795  .def_readwrite ( "zAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::zAxisOrigin )
796  .def_readwrite ( "xCom", &ProSHADE_internal_data::ProSHADE_data::xCom )
797  .def_readwrite ( "yCom", &ProSHADE_internal_data::ProSHADE_data::yCom )
798  .def_readwrite ( "zCom", &ProSHADE_internal_data::ProSHADE_data::zCom )
799  .def_readwrite ( "xFrom", &ProSHADE_internal_data::ProSHADE_data::xFrom )
800  .def_readwrite ( "yFrom", &ProSHADE_internal_data::ProSHADE_data::yFrom )
801  .def_readwrite ( "zFrom", &ProSHADE_internal_data::ProSHADE_data::zFrom )
802  .def_readwrite ( "xTo", &ProSHADE_internal_data::ProSHADE_data::xTo )
803  .def_readwrite ( "yTo", &ProSHADE_internal_data::ProSHADE_data::yTo )
804  .def_readwrite ( "zTo", &ProSHADE_internal_data::ProSHADE_data::zTo )
805  .def_readwrite ( "spherePos", &ProSHADE_internal_data::ProSHADE_data::spherePos )
806  .def_readwrite ( "noSpheres", &ProSHADE_internal_data::ProSHADE_data::noSpheres )
807  .def_readwrite ( "maxShellBand", &ProSHADE_internal_data::ProSHADE_data::maxShellBand )
808  .def_readwrite ( "maxEMatDim", &ProSHADE_internal_data::ProSHADE_data::maxEMatDim )
809  .def_readwrite ( "isEmpty", &ProSHADE_internal_data::ProSHADE_data::isEmpty )
810  .def_readwrite ( "inputOrder", &ProSHADE_internal_data::ProSHADE_data::inputOrder )
811 
812  //============================================ Description
813  .def ( "__repr__", [] ( ) { return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );
814 
815  //================================================ Export extra symmetry elements functions
816  pyProSHADE.def ( "computeGroupElementsForGroup",
817  [] ( proshade_double x, proshade_double y, proshade_double z, proshade_unsign fold ) -> pybind11::array_t < proshade_double >
818  {
819  //== Get the results
820  std::vector<std::vector< proshade_double > > retVec = ProSHADE_internal_data::computeGroupElementsForGroup ( x, y, z, static_cast< proshade_signed > ( fold ) );
821 
822  //== Allocate memory for the numpy values
823  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
824  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
825 
826  //== Copy the values
827  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( retVec.size() ); iter++ ) { for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[(iter*9)+it] = retVec.at(iter).at(it); } }
828 
829  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
830  pybind11::capsule pyCapsuleGrEls ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
831 
832  //== Create python readable object with the correct memory access
833  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
834  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
835  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
836  npVals, // Data
837  pyCapsuleGrEls ); // Capsule (destructor)
838 
839  //== Done
840  return ( retArr );
841  }, "This function takes the axis and fold and computes all resulting group elements (as rotation matrices), returning them in numpy.ndarray.", pybind11::arg ( "x" ), pybind11::arg ( "y" ), pybind11::arg ( "z" ), pybind11::arg ( "fold" ) );
842  pyProSHADE.def ( "joinElementsFromDifferentGroups",
843  [] ( pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > first,
844  pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > second,
845  proshade_double matrixTolerance,
846  bool combine ) -> pybind11::array_t < proshade_double >
847  {
848  //== Get array buffers
849  pybind11::buffer_info buf1 = first.request();
850  pybind11::buffer_info buf2 = second.request();
851 
852  //== Sanity check
853  if ( buf1.ndim != 3 || buf2.ndim != 3 )
854  {
855  std::cerr << "Function joinElementsFromDifferentGroups() arguments first and second should be numpy.ndarrays with 3 dimensions indexed as follos: first[elementNumber][elementRotationMatrixRow][elementRotationMatrixColumn] - the same format as returned by the computeGroupElementsForGroup() function." << std::endl;
856  return ( pybind11::array_t < proshade_double > () );
857  }
858 
859  //== Convert input to C++ vectors
860  std::vector< std::vector< proshade_double > > fVec, sVec;
861 
862  proshade_double* dataPtr1 = reinterpret_cast < proshade_double* > ( buf1.ptr );
863  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf1.shape.at(0) ); elIt++ )
864  {
865  std::vector< proshade_double > rotMat;
866  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf1.shape.at(1) ); rowIt++ )
867  {
868  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf1.shape.at(2) ); colIt++ )
869  {
870  ProSHADE_internal_misc::addToDoubleVector ( &rotMat, dataPtr1[colIt + static_cast< proshade_unsign > ( buf1.shape.at(2) ) * ( rowIt + static_cast< proshade_unsign > ( buf1.shape.at(1) ) * elIt )] );
871  }
872  }
874  }
875 
876  proshade_double* dataPtr2 = reinterpret_cast < proshade_double* > ( buf2.ptr );
877  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf2.shape.at(0) ); elIt++ )
878  {
879  std::vector< proshade_double > rotMat;
880  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf2.shape.at(1) ); rowIt++ )
881  {
882  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf2.shape.at(2) ); colIt++ )
883  {
884  ProSHADE_internal_misc::addToDoubleVector ( &rotMat, dataPtr2[colIt + static_cast< proshade_unsign > ( buf2.shape.at(2) ) * ( rowIt + static_cast< proshade_unsign > ( buf2.shape.at(1) ) * elIt )] );
885  }
886  }
888  }
889 
890 
891  //== Get the results
892  std::vector< std::vector< proshade_double > > retVec = ProSHADE_internal_data::joinElementsFromDifferentGroups ( &fVec, &sVec, matrixTolerance, combine );
893 
894  //== Allocate memory for the numpy values
895  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
896  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
897 
898  //== Copy the values
899  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( retVec.size() ); iter++ ) { for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[(iter*9)+it] = retVec.at(iter).at(it); } }
900 
901  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
902  pybind11::capsule pyCapsuleGrElsCombo ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
903 
904  //== Create python readable object with the correct memory access
905  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
906  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
907  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
908  npVals, // Data
909  pyCapsuleGrElsCombo ); // Capsule (destructor)
910 
911  //== Done
912  return ( retArr );
913  }, "This function takes the axis and fold and computes all resulting group elements (as rotation matrices), returning them in numpy.ndarray.", pybind11::arg ( "first" ), pybind11::arg ( "second" ), pybind11::arg ( "matrixTolerance" ), pybind11::arg ( "combine" ) );
914  pyProSHADE.def ( "getRotationFunctionSpheres",
915  [] ( ProSHADE_internal_data::ProSHADE_data *dataObj, int fold ) -> pybind11::array_t < proshade_double >
916  {
917 
918  //== Check if rotation function was computed (sanity check)
919  if ( dataObj->getInvSO3Coeffs() == nullptr )
920  {
921  std::cerr << "Function getRotationFunctionSpheres() was called before rotation function was computed. Please compute the rotation function by calling the computeRotationFunction() function first. Returning empty array." << std::endl;
922  return ( pybind11::array_t < proshade_double > () );
923  }
924 
925  //== Check that fold makes sense
926  if ( fold <= 1 )
927  {
928  std::cerr << "The fold for which the rotation function should be mapped is not in the allowed range (i.e. 2+). Returning empty array." << std::endl;
929  return ( pybind11::array_t < proshade_double > () );
930  }
931 
932  //== Compute the sphere mapping
933  std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*> sphereMappedRotFun;
934 
935  //== Convert rotation function to only the required angle-axis space spheres and find all peaks
936  proshade_double soughtAngle = 0.0;
937  for ( proshade_double angIt = 1.0; angIt < static_cast < proshade_double > ( fold ); angIt += 1.0 )
938  {
939  //== Figure the angles to form the symmetry
940  soughtAngle = angIt * ( 2.0 * M_PI / static_cast<proshade_double> ( fold ) );
941 
942  //== Create the angle-axis sphere with correct radius (angle)
943  sphereMappedRotFun.emplace_back ( new ProSHADE_internal_spheres::ProSHADE_rotFun_sphere ( soughtAngle,
944  M_PI / static_cast< proshade_double > ( dataObj->getMaxBand() ),
945  dataObj->getMaxBand ( ) * 2,
946  dataObj->getEMatDim ( ) * 2,
947  soughtAngle,
948  static_cast<proshade_unsign> ( angIt - 1.0 ) ) );
949 
950  //== Interpolate rotation function onto the sphere
951  sphereMappedRotFun.at( static_cast < size_t > ( angIt - 1.0 ))->interpolateSphereValues ( dataObj->getInvSO3Coeffs ( ) );
952  }
953 
954  //== Save values to pointer
955  proshade_double* npVals = new proshade_double[(dataObj->getMaxBand() * 2) * (dataObj->getMaxBand() * 2) * (static_cast< proshade_unsign > ( fold ) - 1)];
956  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
957 
958  //== Copy value to the pointer
959  std::vector< std::vector < proshade_double > > vls;
960  size_t valIter = 0;
961  for ( size_t sphIt = 0; sphIt < sphereMappedRotFun.size(); sphIt++ )
962  {
963  vls.clear();
964  vls = sphereMappedRotFun.at(sphIt)->getCopyOfValues ();
965 
966  for ( size_t lonIt = 0; lonIt < (dataObj->getMaxBand() * 2); lonIt++ )
967  {
968  for ( size_t latIt = 0; latIt < (dataObj->getMaxBand() * 2); latIt++ )
969  {
970  valIter = latIt + static_cast< size_t > ( dataObj->getMaxBand() * 2 ) * ( lonIt + static_cast< size_t > ( dataObj->getMaxBand() * 2 ) * sphIt );
971  npVals[valIter] = vls.at(lonIt).at(latIt);
972  }
973  }
974  }
975 
976  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
977  pybind11::capsule pyCapsuleRFMapSph ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
978 
979  //== Create python readable object with the correct memory access
980  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
981  { static_cast< int > ( sphereMappedRotFun.size() ), static_cast< int > (dataObj->getMaxBand() * 2), static_cast< int > (dataObj->getMaxBand() * 2) }, // Shape
982  { (dataObj->getMaxBand() * 2) * (dataObj->getMaxBand() * 2) * sizeof(proshade_double), (dataObj->getMaxBand() * 2) * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
983  npVals, // Data
984  pyCapsuleRFMapSph ); // Capsule (destructor)
985 
986  //== Done
987  return ( retArr );
988  }, "This function allows access to the (self-)rotation function mapping to spheres with the angle being the radius and the longitude/latitude angles being the axis of the rotation for a particular fold (i.e. for all angles required by the fold except for the angle 0)." );
989 }
ProSHADE_internal_mapManip::addExtraBoundSpace
void addExtraBoundSpace(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_single xAngs, proshade_single yAngs, proshade_single zAngs, proshade_signed *&bounds, proshade_single extraSpace)
This function takes a set of bounds and adds a given number of Angstroms to them.
Definition: ProSHADE_mapManip.cpp:1193
ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics
void computeSphericalHarmonics(ProSHADE_settings *settings)
This function computes the spherical harmonics decomposition for the whole structure.
Definition: ProSHADE_data.cpp:1895
ProSHADE_internal_data::ProSHADE_data::getOverlayRotationFunction
void getOverlayRotationFunction(ProSHADE_settings *settings, ProSHADE_internal_data::ProSHADE_data *obj2)
This function computes the overlay rotation function (i.e. the correlation function in SO(3) space).
Definition: ProSHADE_overlay.cpp:35
ProSHADE_internal_data::ProSHADE_data::fileName
std::string fileName
This is the original file from which the data were obtained.
Definition: ProSHADE_data.hpp:52
ProSHADE_internal_data::ProSHADE_data::zFrom
proshade_signed zFrom
This is the starting index along the z axis.
Definition: ProSHADE_data.hpp:112
ProSHADE_internal_data::ProSHADE_data::xDimIndices
proshade_unsign xDimIndices
This is the size of the map cell x dimension in indices.
Definition: ProSHADE_data.hpp:65
ProSHADE_internal_data::ProSHADE_data::getZDimSize
proshade_single getZDimSize(void)
This function allows access to the map size in angstroms along the Z axis.
Definition: ProSHADE_data.cpp:4064
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_data::ProSHADE_data::xGridIndices
proshade_unsign xGridIndices
As far as I know, this is identical to the xDimIndices.
Definition: ProSHADE_data.hpp:68
ProSHADE_internal_data::ProSHADE_data::computeRotationFunction
void computeRotationFunction(ProSHADE_settings *settings)
This function computes the self-rotation function for this structure.
Definition: ProSHADE_symmetry.cpp:41
ProSHADE_internal_data::joinElementsFromDifferentGroups
std::vector< std::vector< proshade_double > > joinElementsFromDifferentGroups(std::vector< std::vector< proshade_double > > *first, std::vector< std::vector< proshade_double > > *second, proshade_double matrixTolerance, bool combine)
This function joins two group element lists using only unique elements.
Definition: ProSHADE_data.cpp:3064
ProSHADE_internal_data::ProSHADE_data::xAxisOrder
proshade_unsign xAxisOrder
This is the order of the x axis.
Definition: ProSHADE_data.hpp:71
ProSHADE_internal_data::ProSHADE_data::cAngle
proshade_single cAngle
This is the angle c of the map cell in degrees.
Definition: ProSHADE_data.hpp:64
ProSHADE_internal_data::ProSHADE_data::yCom
proshade_double yCom
The COM of the map after processing along the Y-axis.
Definition: ProSHADE_data.hpp:78
ProSHADE_settings::boundsExtraSpace
proshade_single boundsExtraSpace
The number of extra angstroms to be added to all re-boxing bounds just for safety.
Definition: ProSHADE_settings.hpp:97
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_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_data::ProSHADE_data::getXDimSize
proshade_single getXDimSize(void)
This function allows access to the map size in angstroms along the X axis.
Definition: ProSHADE_data.cpp:4044
ProSHADE_internal_data::ProSHADE_data
This class contains all inputed and derived data for a single structure.
Definition: ProSHADE_data.hpp:49
ProSHADE_internal_data::ProSHADE_data::noSpheres
proshade_unsign noSpheres
The number of spheres with map projected onto them.
Definition: ProSHADE_data.hpp:119
ProSHADE_internal_data::ProSHADE_data::getYDimSize
proshade_single getYDimSize(void)
This function allows access to the map size in angstroms along the Y axis.
Definition: ProSHADE_data.cpp:4054
ProSHADE_internal_data::ProSHADE_data::inputOrder
proshade_unsign inputOrder
This value is the input order - it is useful to know for writing out files, so that they would not ov...
Definition: ProSHADE_data.hpp:140
ProSHADE_internal_data::ProSHADE_data::zDimSize
proshade_single zDimSize
This is the size of the map cell z dimension in Angstroms.
Definition: ProSHADE_data.hpp:61
ProSHADE_settings::requestedResolution
proshade_single requestedResolution
The resolution to which the calculations are to be done.
Definition: ProSHADE_settings.hpp:50
ProSHADE_internal_data::ProSHADE_data::maxShellBand
proshade_unsign maxShellBand
The maximum band for any shell of the object.
Definition: ProSHADE_data.hpp:123
ProSHADE_internal_data::ProSHADE_data::getSpherePositions
void getSpherePositions(ProSHADE_settings *settings)
This function determines the sphere positions (radii) for sphere mapping.
Definition: ProSHADE_data.cpp:1790
ProSHADE_internal_data::ProSHADE_data::invertMirrorMap
void invertMirrorMap(ProSHADE_settings *settings)
Function for inverting the map to its mirror image.
Definition: ProSHADE_data.cpp:1211
ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryFold
proshade_unsign getRecommendedSymmetryFold(void)
This function simply returns the detected recommended symmetry fold.
Definition: ProSHADE_data.cpp:4609
ProSHADE_internal_data::ProSHADE_data::centreMapOnCOM
void centreMapOnCOM(ProSHADE_settings *settings)
This function shits the map so that its COM is in the centre of the map.
Definition: ProSHADE_data.cpp:1556
ProSHADE_internal_data::ProSHADE_data::yDimIndices
proshade_unsign yDimIndices
This is the size of the map cell y dimension in indices.
Definition: ProSHADE_data.hpp:66
ProSHADE_internal_data::ProSHADE_data::yGridIndices
proshade_unsign yGridIndices
As far as I know, this is identical to the yDimIndices.
Definition: ProSHADE_data.hpp:69
ProSHADE_internal_data::ProSHADE_data::zAxisOrigin
proshade_signed zAxisOrigin
This is the origin position along the z axis.
Definition: ProSHADE_data.hpp:76
ProSHADE_internal_misc::addToDoubleVector
void addToDoubleVector(std::vector< proshade_double > *vecToAddTo, proshade_double elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:77
ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryType
std::string getRecommendedSymmetryType(void)
This function simply returns the detected recommended symmetry type.
Definition: ProSHADE_data.cpp:4600
ProSHADE_internal_data::ProSHADE_data::shiftToRotationCentre
void shiftToRotationCentre(ProSHADE_settings *settings)
Function for shifting map so that its rotation centre is at the centre of the box.
Definition: ProSHADE_data.cpp:1164
ProSHADE_internal_data::ProSHADE_data::getXDim
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
Definition: ProSHADE_data.cpp:4074
ProSHADE_internal_data::ProSHADE_data::rotateMapReciprocalSpace
void rotateMapReciprocalSpace(ProSHADE_settings *settings, proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma)
This function rotates a map based on the given Euler angles.
Definition: ProSHADE_overlay.cpp:571
ProSHADE_internal_data::ProSHADE_data::zTo
proshade_signed zTo
This is the final index along the z axis.
Definition: ProSHADE_data.hpp:115
ProSHADE_internal_data::ProSHADE_data::xTo
proshade_signed xTo
This is the final index along the x axis.
Definition: ProSHADE_data.hpp:113
ProSHADE_internal_data::ProSHADE_data::mapToSpheres
void mapToSpheres(ProSHADE_settings *settings)
This function converts the internal map onto a set of concentric spheres.
Definition: ProSHADE_data.cpp:1851
ProSHADE_internal_data::ProSHADE_data::xAxisOrigin
proshade_signed xAxisOrigin
This is the origin position along the x axis.
Definition: ProSHADE_data.hpp:74
ProSHADE_internal_data::ProSHADE_data::yTo
proshade_signed yTo
This is the final index along the y axis.
Definition: ProSHADE_data.hpp:114
ProSHADE_internal_data::ProSHADE_data::xDimSize
proshade_single xDimSize
This is the size of the map cell x dimension in Angstroms.
Definition: ProSHADE_data.hpp:59
ProSHADE_internal_data::ProSHADE_data::yAxisOrigin
proshade_signed yAxisOrigin
This is the origin position along the y axis.
Definition: ProSHADE_data.hpp:75
ProSHADE_internal_data::ProSHADE_data::computeTranslationMap
void computeTranslationMap(ProSHADE_internal_data::ProSHADE_data *obj1)
This function does the computation of the translation map and saves results internally.
Definition: ProSHADE_overlay.cpp:324
ProSHADE_internal_data::ProSHADE_data::yAxisOrder
proshade_unsign yAxisOrder
This is the order of the y axis.
Definition: ProSHADE_data.hpp:72
ProSHADE_internal_data::ProSHADE_data::writePdb
void writePdb(std::string fName, proshade_double euA=0.0, proshade_double euB=0.0, proshade_double euG=0.0, proshade_double trsX=0.0, proshade_double trsY=0.0, proshade_double trsZ=0.0, proshade_double rotX=0.0, proshade_double rotY=0.0, proshade_double rotZ=0.0, bool firstModel=true)
This function writes out the co-ordinates file with ProSHADE type rotation and translation applied.
Definition: ProSHADE_data.cpp:983
ProSHADE_internal_data::ProSHADE_data::xFrom
proshade_signed xFrom
This is the starting index along the x axis.
Definition: ProSHADE_data.hpp:110
ProSHADE_internal_data::ProSHADE_data::writeMap
void writeMap(std::string fName, std::string title="Created by ProSHADE and written by GEMMI", int mode=2)
Function for writing out the internal structure representation in MRC MAP format.
Definition: ProSHADE_data.cpp:916
ProSHADE_internal_data::ProSHADE_data::getMaxBand
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
Definition: ProSHADE_data.cpp:3748
ProSHADE_internal_data::ProSHADE_data::shiftToBoxCentre
void shiftToBoxCentre(ProSHADE_settings *settings)
Function for shifting map so that its centre of box is at required position.
Definition: ProSHADE_data.cpp:1099
ProSHADE_internal_data::ProSHADE_data::translateMap
void translateMap(proshade_double trsX, proshade_double trsY, proshade_double trsZ)
This function simply translates the map by a given number of Angstroms along the three axes.
Definition: ProSHADE_overlay.cpp:1064
ProSHADE_internal_data::ProSHADE_data::writeGemmi
void writeGemmi(std::string fName, gemmi::Structure gemmiStruct, proshade_double euA=0.0, proshade_double euB=0.0, proshade_double euG=0.0, proshade_double trsX=0.0, proshade_double trsY=0.0, proshade_double trsZ=0.0, proshade_double rotX=0.0, proshade_double rotY=0.0, proshade_double rotZ=0.0, bool firstModel=true)
This function writes out the gemmi::Structure object with ProSHADE type rotation and translation appl...
Definition: ProSHADE_data.cpp:1020
ProSHADE_internal_data::ProSHADE_data::getYDim
proshade_unsign getYDim(void)
This function allows access to the map size in indices along the Y axis.
Definition: ProSHADE_data.cpp:4084
ProSHADE_internal_data::ProSHADE_data::maskMap
void maskMap(ProSHADE_settings *settings)
Function for computing the map mask using blurring and X IQRs from median.
Definition: ProSHADE_data.cpp:1312
ProSHADE_internal_mapManip::getNonZeroBounds
void getNonZeroBounds(proshade_double *map, proshade_signed xDim, proshade_signed yDim, proshade_signed zDim, proshade_signed *&ret)
Function for finding the map boundaries enclosing positive only values.
Definition: ProSHADE_mapManip.cpp:1136
ProSHADE_settings
This class stores all the settings and is passed to the executive classes instead of a multitude of p...
Definition: ProSHADE_settings.hpp:37
ProSHADE_internal_data::ProSHADE_data::xCom
proshade_double xCom
The COM of the map after processing along the X-axis.
Definition: ProSHADE_data.hpp:77
ProSHADE_internal_data::ProSHADE_data::processInternalMap
void processInternalMap(ProSHADE_settings *settings)
This function simply clusters several other functions which should be called together.
Definition: ProSHADE_data.cpp:1739
ProSHADE_internal_data::ProSHADE_data::zCom
proshade_double zCom
The COM of the map after processing along the Z-axis.
Definition: ProSHADE_data.hpp:79
ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpaceInPlace
std::vector< proshade_double > rotateMapRealSpaceInPlace(proshade_double eulA, proshade_double eulB, proshade_double eulG)
This function rotates a map based on the given Euler angles in place.
Definition: ProSHADE_overlay.cpp:805
ProSHADE_internal_data::ProSHADE_data::spherePos
std::vector< proshade_single > spherePos
Vector of sphere radii from the centre of the map.
Definition: ProSHADE_data.hpp:118
ProSHADE_internal_data::ProSHADE_data::getZDim
proshade_unsign getZDim(void)
This function allows access to the map size in indices along the Z axis.
Definition: ProSHADE_data.cpp:4094
ProSHADE_settings::boundsSimilarityThreshold
proshade_signed boundsSimilarityThreshold
Number of indices which can be added just to make sure same size in indices is achieved.
Definition: ProSHADE_settings.hpp:98
ProSHADE_internal_data::ProSHADE_data::zAxisOrder
proshade_unsign zAxisOrder
This is the order of the z axis.
Definition: ProSHADE_data.hpp:73
ProSHADE_internal_mapManip::beautifyBoundaries
void beautifyBoundaries(proshade_signed *&bounds, proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim, proshade_signed boundsDiffThres)
Function for modifying boundaries to a mathematically more pleasant values.
Definition: ProSHADE_mapManip.cpp:1976
ProSHADE_internal_data::ProSHADE_data::reSampleMap
void reSampleMap(ProSHADE_settings *settings)
This function changes the internal map sampling to conform to particular resolution value.
Definition: ProSHADE_data.cpp:1476
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:73
ProSHADE_internal_data::ProSHADE_data::maxEMatDim
proshade_unsign maxEMatDim
The band (l) value for E matrix (i.e. the smallest of the two bands).
Definition: ProSHADE_data.hpp:124
ProSHADE_internal_data::ProSHADE_data::bAngle
proshade_single bAngle
This is the angle b of the map cell in degrees.
Definition: ProSHADE_data.hpp:63
ProSHADE_internal_data::ProSHADE_data::aAngle
proshade_single aAngle
This is the angle a of the map cell in degrees.
Definition: ProSHADE_data.hpp:62
ProSHADE_internal_misc::addToDoubleVectorVector
void addToDoubleVectorVector(std::vector< std::vector< proshade_double > > *vecToAddTo, std::vector< proshade_double > elementToAdd)
Adds the element to the vector of vectors.
Definition: ProSHADE_misc.cpp:233
ProSHADE_internal_data::ProSHADE_data::isEmpty
bool isEmpty
This variable stated whether the class contains any information.
Definition: ProSHADE_data.hpp:139
ProSHADE_internal_data::ProSHADE_data::getInvSO3Coeffs
proshade_complex * getInvSO3Coeffs(void)
This function allows access to the inverse SO(3) coefficients array.
Definition: ProSHADE_data.cpp:3992
ProSHADE_internal_data::ProSHADE_data::getEMatDim
proshade_unsign getEMatDim(void)
This function allows access to the maximum band for the E matrix.
Definition: ProSHADE_data.cpp:4014
ProSHADE_settings::forceBounds
proshade_signed * forceBounds
These will be the boundaries to be forced upon the map.
Definition: ProSHADE_settings.hpp:100
ProSHADE_internal_data::ProSHADE_data::yFrom
proshade_signed yFrom
This is the starting index along the y axis.
Definition: ProSHADE_data.hpp:111
ProSHADE_internal_misc::addToUnsignVector
void addToUnsignVector(std::vector< proshade_unsign > *vecToAddTo, proshade_unsign elementToAdd)
Adds the element to the vector.
Definition: ProSHADE_misc.cpp:99
ProSHADE_settings::useSameBounds
bool useSameBounds
Switch to say that the same boundaries as used for the first should be used for all input maps.
Definition: ProSHADE_settings.hpp:99
ProSHADE_internal_data::ProSHADE_data::zDimIndices
proshade_unsign zDimIndices
This is the size of the map cell z dimension in indices.
Definition: ProSHADE_data.hpp:67
ProSHADE_internal_data::computeGroupElementsForGroup
std::vector< std::vector< proshade_double > > computeGroupElementsForGroup(proshade_double xAx, proshade_double yAx, proshade_double zAx, proshade_signed fold)
This function computes the group elements as rotation matrices (except for the identity element) for ...
Definition: ProSHADE_data.cpp:2911
ProSHADE_internal_data::ProSHADE_data::zGridIndices
proshade_unsign zGridIndices
As far as I know, this is identical to the zDimIndices.
Definition: ProSHADE_data.hpp:70
ProSHADE_internal_data::ProSHADE_data::normaliseMap
void normaliseMap(ProSHADE_settings *settings)
Function for normalising the map values to mean 0 and sd 1..
Definition: ProSHADE_data.cpp:1265
ProSHADE_internal_data::ProSHADE_data::yDimSize
proshade_single yDimSize
This is the size of the map cell y dimension in Angstroms.
Definition: ProSHADE_data.hpp:60
ProSHADE_internal_data::ProSHADE_data::removePhaseInormation
void removePhaseInormation(ProSHADE_settings *settings)
This function removes phase from the map, effectively converting it to Patterson map.
Definition: ProSHADE_data.cpp:3794
ProSHADE_internal_data::ProSHADE_data::zeroPaddToDims
void zeroPaddToDims(proshade_unsign xDim, proshade_unsign yDim, proshade_unsign zDim)
This function changes the size of a structure to fit the supplied new limits.
Definition: ProSHADE_overlay.cpp:433
ProSHADE_internal_data::ProSHADE_data::addExtraSpace
void addExtraSpace(ProSHADE_settings *settings)
This function increases the size of the map so that it can add empty space around it.
Definition: ProSHADE_data.cpp:1629