ProSHADE  0.7.6.0 (JUL 2021)
Protein Shape Detection
pyProSHADE_data.cpp
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", &ProSHADE_internal_data::ProSHADE_data::readInStructure, "This function initialises the basic ProSHADE_data variables and reads in a single structure.", pybind11::arg ( "fname" ), pybind11::arg ( "inputO" ), pybind11::arg ( "settings" ) )
91  .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 )
92  .def ( "writePdb", &ProSHADE_internal_data::ProSHADE_data::writePdb, "This function writes out the PDB formatted file coresponding to the structure.", 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 ( "firstModel" ) = true )
93  .def ( "getMap",
94  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < proshade_double >
95  {
96  //== Copy the values
97  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { self.xDimIndices, self.yDimIndices, self.zDimIndices }, // Shape
98  { self.yDimIndices * self.zDimIndices * sizeof(proshade_double),
99  self.zDimIndices * sizeof(proshade_double),
100  sizeof(proshade_double) }, // C-stype strides
101  self.internalMap ); // Data
102 
103  //== Done
104  return ( retArr );
105  } )
106 
107  //============================================ Data processing functions
108  .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" ) )
109  .def ( "invertMirrorMap", &ProSHADE_internal_data::ProSHADE_data::invertMirrorMap, "Function for inverting the map to its mirror image.", pybind11::arg ( "settings" ) )
110  .def ( "normaliseMap", &ProSHADE_internal_data::ProSHADE_data::normaliseMap, "Function for normalising the map values to mean 0 and sd 1.", pybind11::arg ( "settings" ) )
111  .def ( "maskMap", &ProSHADE_internal_data::ProSHADE_data::maskMap, "Function for computing the map mask using blurring and X IQRs from median.", pybind11::arg ( "settings" ) )
112  .def ( "reSampleMap", &ProSHADE_internal_data::ProSHADE_data::reSampleMap, "This function changes the internal map sampling to conform to particular resolution value.", pybind11::arg ( "settings" ) )
113  .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" ) )
114  .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" ) )
115  .def ( "removePhaseInormation", &ProSHADE_internal_data::ProSHADE_data::removePhaseInormation, "This function removes phase from the map, effectively converting it to Patterson map.", pybind11::arg ( "settings" ) )
116  .def ( "getReBoxBoundaries",
117  [] ( ProSHADE_internal_data::ProSHADE_data &self,ProSHADE_settings* settings ) -> pybind11::array_t < proshade_signed >
118  {
119  //== Allocate output memory
120  proshade_signed* retVals = new proshade_signed[6];
121  ProSHADE_internal_misc::checkMemoryAllocation ( retVals, __FILE__, __LINE__, __func__ );
122 
123  //== If same bounds as first one are required, test if possible and return these instead
124  if ( settings->useSameBounds && ( self.inputOrder != 0 ) )
125  {
126  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->forceBounds[iter]; }
127  }
128  //== In this case, bounds need to be found de novo
129  else
130  {
131  //== Find the non-zero bounds
133  static_cast< proshade_signed > ( self.xDimIndices ),
134  static_cast< proshade_signed > ( self.yDimIndices ),
135  static_cast< proshade_signed > ( self.zDimIndices ),
136  retVals );
137 
138  //== Add the extra space
139  ProSHADE_internal_mapManip::addExtraBoundSpace ( self.xDimIndices, self.yDimIndices, self.zDimIndices,
140  self.xDimSize, self.yDimSize, self.zDimSize, retVals, settings->boundsExtraSpace );
141 
142  //== Beautify boundaries
143  ProSHADE_internal_mapManip::beautifyBoundaries ( retVals, self.xDimIndices, self.yDimIndices, self.zDimIndices, settings->boundsSimilarityThreshold );
144 
145  //== If need be, save boundaries to be used for all other structure
146  if ( settings->useSameBounds && ( self.inputOrder == 0 ) ) { for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->forceBounds[iter] = retVals[iter]; } }
147  }
148 
149  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
150  pybind11::capsule pyCapsuleReBox2 ( retVals, []( void *f ) { proshade_signed* foo = reinterpret_cast< proshade_signed* > ( f ); delete foo; } );
151 
152  //== Copy the value
153  pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 }, // Shape
154  { sizeof(proshade_signed) }, // C-stype strides
155  retVals, // Data
156  pyCapsuleReBox2 ); // Capsule
157 
158  //== Done
159  return ( retArr );
160  }, "This function finds the boundaries enclosing positive map values and adds some extra space." )
161  .def ( "createNewMapFromBounds",
162  [] ( ProSHADE_internal_data::ProSHADE_data &self, pybind11::array_t < proshade_signed > bounds, ProSHADE_internal_data::ProSHADE_data* newStr, ProSHADE_settings* settings )
163  {
164  //== Allocate memory for bounds copy
165  proshade_signed* newBounds = new proshade_signed[6];
166  ProSHADE_internal_misc::checkMemoryAllocation ( newBounds, __FILE__, __LINE__, __func__ );
167 
168  //== Copy to the pointer
169  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
170 
171  //== Run C++ function
172  self.createNewMapFromBounds ( settings, newStr, newBounds );
173 
174  //== Done
175  return ;
176  }, "This function creates a new structure from the calling structure and new bounds values." )
177 
178  //============================================ Data sphere mapping functions
179  .def ( "mapToSpheres", &ProSHADE_internal_data::ProSHADE_data::mapToSpheres, "This function converts the internal map onto a set of concentric spheres.", pybind11::arg ( "settings" ) )
180  .def ( "getSpherePositions", &ProSHADE_internal_data::ProSHADE_data::getSpherePositions, "This function determines the sphere positions (radii) for sphere mapping.", pybind11::arg ( "settings" ) )
181  .def ( "computeSphericalHarmonics", &ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics, "This function computes the spherical harmonics decomposition for the whole structure.", pybind11::arg ( "settings" ) )
182 
183  //============================================ Accessor functions
184  .def ( "getXDimSize", &ProSHADE_internal_data::ProSHADE_data::getXDimSize, "This function allows access to the map size in angstroms along the X axis." )
185  .def ( "getYDimSize", &ProSHADE_internal_data::ProSHADE_data::getYDimSize, "This function allows access to the map size in angstroms along the Y axis." )
186  .def ( "getZDimSize", &ProSHADE_internal_data::ProSHADE_data::getZDimSize, "This function allows access to the map size in angstroms along the Z axis." )
187  .def ( "getXDim", &ProSHADE_internal_data::ProSHADE_data::getXDim, "This function allows access to the map size in indices along the X axis." )
188  .def ( "getYDim", &ProSHADE_internal_data::ProSHADE_data::getYDim, "This function allows access to the map size in indices along the Y axis." )
189  .def ( "getZDim", &ProSHADE_internal_data::ProSHADE_data::getZDim, "This function allows access to the map size in indices along the Z axis." )
190 
191  //============================================ Symmetry related functions
192  .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" ) )
193  .def ( "detectSymmetryInStructure",
195  {
196  //== Call the appropriate C++ function
197  if ( settings->usePeakSearchInRotationFunctionSpace )
198  {
199  //== Detect point groups in the angle-axis space
200  self.detectSymmetryFromAngleAxisSpace ( settings, &settings->detectedSymmetry, &settings->allDetectedCAxes );
201  }
202  else
203  {
204  //== Detect symmetry using the peak detection in rotation function space
205  self.detectSymmetryInStructure ( settings, &settings->detectedSymmetry, &settings->allDetectedCAxes );
206  }
207  }, "This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg ( "settings" ) )
208  .def ( "getRecommendedSymmetryType", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryType, "This function simply returns the detected recommended symmetry type.", pybind11::arg ( "settings" ) )
209  .def ( "getRecommendedSymmetryFold", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryFold, "This function simply returns the detected recommended symmetry fold.", pybind11::arg ( "settings" ) )
210  .def ( "getRecommendedSymmetryAxes",
211  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
212  {
213  //== Allocate memory for the numpy values
214  float* npVals = new float[static_cast<unsigned int> ( settings->detectedSymmetry.size() * 7 )];
215  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
216 
217  //== Copy values
218  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->detectedSymmetry.size() ); iter++ ) { for ( proshade_unsign it = 0; it < 7; it++ ) { npVals[(iter*7)+it] = static_cast< float > ( settings->detectedSymmetry.at(iter)[it] ); } }
219 
220  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
221  pybind11::capsule pyCapsuleStrRecSym ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
222 
223  //== Copy the value
224  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( settings->detectedSymmetry.size() ), static_cast<int> ( 7 ) }, // Shape
225  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
226  npVals, // Data
227  pyCapsuleStrRecSym ); // Capsule
228 
229  //== Done
230  return ( retArr );
231  }, "This function returns the recommended symmetry axes as a 2D numpy array." )
232  .def ( "getAllCSyms",
233  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
234  {
235  //== Allocate memory for the numpy values
236  float* npVals = new float[static_cast<unsigned int> ( settings->allDetectedCAxes.size() * 7 )];
237  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
238 
239  //== Copy values
240  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->allDetectedCAxes.size() ); iter++ ) { for ( proshade_unsign it = 0; it < 7; it++ ) { npVals[(iter*7)+it] = static_cast< float > ( settings->allDetectedCAxes.at(iter).at(it) ); } }
241 
242  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
243  pybind11::capsule pyCapsuleStrSymList ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
244 
245  //== Copy the value
246  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( settings->allDetectedCAxes.size() ), 7 }, // Shape
247  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
248  npVals, // Data
249  pyCapsuleStrSymList ); // Capsule
250 
251  //== Done
252  return ( retArr );
253  }, "This function returns all symmetry axes as a 2D numpy array." )
254  .def ( "getNonCSymmetryAxesIndices",
255  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::dict
256  {
257  //== Initialise local variables
258  pybind11::dict retDict;
259  pybind11::list dList, tList, oList, iList;
260 
261  //== Fill in the D list
262  for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( settings->allDetectedDAxes.size() ); dIt++ )
263  {
264  pybind11::list memList;
265  for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( settings->allDetectedDAxes.at(dIt).size() ); memIt++ )
266  {
267  memList.append ( settings->allDetectedDAxes.at(dIt).at(memIt) );
268  }
269  dList.append ( memList );
270  }
271 
272  //== Fill in the T list
273  for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( settings->allDetectedTAxes.size() ); tIt++ )
274  {
275  tList.append ( settings->allDetectedTAxes.at ( tIt ) );
276  }
277 
278  //== Fill in the O list
279  for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( settings->allDetectedOAxes.size() ); oIt++ )
280  {
281  oList.append ( settings->allDetectedOAxes.at ( oIt ) );
282  }
283 
284  //== Fill in the T list
285  for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( settings->allDetectedIAxes.size() ); iIt++ )
286  {
287  iList.append ( settings->allDetectedIAxes.at ( iIt ) );
288  }
289 
290  //== Save the lists to the dict
291  retDict[ pybind11::handle ( pybind11::str ( "D" ).ptr ( ) ) ] = dList;
292  retDict[ pybind11::handle ( pybind11::str ( "T" ).ptr ( ) ) ] = tList;
293  retDict[ pybind11::handle ( pybind11::str ( "O" ).ptr ( ) ) ] = oList;
294  retDict[ pybind11::handle ( pybind11::str ( "I" ).ptr ( ) ) ] = iList;
295 
296  //== Done
297  return ( retDict );
298  }, "This function returns array of non-C axes indices." )
299  .def ( "getAllGroupElements",
300  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings, pybind11::array_t < int > axList, std::string groupType, proshade_double matrixTolerance ) -> pybind11::list
301  {
302  //== Sanity check
303  pybind11::buffer_info buf = axList.request();
304  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 ); }
305 
306  //== Convert to vector of unsigns
307  std::vector< proshade_unsign > axesList;
308  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) ) ); }
309 
310  //== Get the results
311  std::vector < std::vector< proshade_double > > vals = self.getAllGroupElements ( settings, axesList, groupType, matrixTolerance );
312 
313  //== Initialise return variable
314  pybind11::list retList;
315 
316  //== Copy values
317  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
318  {
319  //== Allocate memory for the numpy values
320  float* npVals = new float[static_cast<unsigned int> ( 9 )];
321  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
322 
323  //== Copy values to memory
324  for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] = static_cast< float > ( vals.at(iter).at(it) ); }
325 
326  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
327  pybind11::capsule pyCapsuleGrpEl ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
328 
329  //== Copy the value
330  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 }, // Shape
331  { 3 * sizeof(float), sizeof(float) }, // C-stype strides
332  npVals, // Data
333  pyCapsuleGrpEl ); // Capsule
334 
335  //== Save matrix
336  retList.append ( retArr );
337  }
338 
339  //== Done
340  return ( retList );
341  }, "This function returns the group elements as rotation matrices of any point group described by the detected axes.", pybind11::arg ( "settings" ), pybind11::arg ( "axList" ), pybind11::arg ( "groupType" ) = "", pybind11::arg( "matrixTolerance" ) = 0.05 )
342 
343  .def ( "getMapCOMProcessChange",
344  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < float >
345  {
346  //== Get the values
347  std::vector< proshade_double > vals = self.getMapCOMProcessChange ();
348 
349  //== Allocate memory for the numpy values
350  float* npVals = new float[static_cast<unsigned int> ( 3 )];
351  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
352 
353  //== Copy values
354  for ( proshade_unsign iter = 0; iter < 3; iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
355 
356  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
357  pybind11::capsule pyCapsuleSymShiftDat ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
358 
359  //== Copy the value
360  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( vals.size() ) }, // Shape
361  { sizeof(float) }, // C-stype strides
362  npVals, // Data
363  pyCapsuleSymShiftDat ); // Capsule
364 
365  //== Done
366  return ( retArr );
367  }, "This function returns the shift in Angstrom applied to the internal map representation in order to align its COM with the centre of box." )
368 
369  //============================================ Overlay related functions
370  .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" ) )
371  .def ( "getBestRotationMapPeaksEulerAngles",
372  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
373  {
374  //== Get values
375  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
376 
377  //== Allocate memory for the numpy values
378  float* npVals = new float[static_cast<unsigned int> (vals.size())];
379  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
380 
381  //== Copy values
382  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
383 
384  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
385  pybind11::capsule pyCapsuleEuPeak ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
386 
387  //== Copy the value
388  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<unsigned int> (vals.size()) }, // Shape
389  { sizeof(float) }, // C-stype strides
390  npVals, // Data
391  pyCapsuleEuPeak ); // Capsule (C++ destructor, basically)
392 
393  //== Done
394  return ( retArr );
395  }, "This function returns a vector of three floats, the three Euler angles of the best peak in the rotation map.", pybind11::arg ( "settings" ) )
396  .def ( "getBestRotationMapPeaksRotationMatrix",
397  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < proshade_double >
398  {
399  //== Get values
400  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
401 
402  //== Convert Euler ZXZ to matrix
403  proshade_double* retMat = new proshade_double[9];
404  ProSHADE_internal_misc::checkMemoryAllocation ( retMat, __FILE__, __LINE__, __func__ );
405  ProSHADE_internal_maths::getRotationMatrixFromEulerZXZAngles ( vals.at(0), vals.at(1), vals.at(2), retMat );
406 
407  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
408  pybind11::capsule pyCapsuleRMPeak ( retMat, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
409 
410  //== Copy the value
411  pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 }, // Shape
412  { 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
413  retMat, // Data
414  pyCapsuleRMPeak ); // Capsule
415 
416  //== Done
417  return ( retArr );
418  }, "This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg ( "settings" ) )
419  .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" ) )
420  .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" ) )
421  .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" ) )
422  .def ( "computeTranslationMap", &ProSHADE_internal_data::ProSHADE_data::computeTranslationMap, "This function does the computation of the translation map and saves results internally.", pybind11::arg ( "staticStructure" ) )
423  .def ( "getOverlayTranslations",
424  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_internal_data::ProSHADE_data* staticStructure , proshade_double eulA, proshade_double eulB, proshade_double eulG) -> pybind11::dict
425  {
426  //== Get values
427  std::vector< proshade_double > vals = self.getBestTranslationMapPeaksAngstrom ( staticStructure, eulA, eulB, eulG );
428 
429  //== Initialise variables
430  pybind11::dict retDict;
431  pybind11::list rotCen, toOverlay;
432 
433  //== Copy values
434  rotCen.append ( self.originalPdbRotCenX );
435  rotCen.append ( self.originalPdbRotCenY );
436  rotCen.append ( self.originalPdbRotCenZ );
437 
438  toOverlay.append ( self.originalPdbTransX );
439  toOverlay.append ( self.originalPdbTransY );
440  toOverlay.append ( self.originalPdbTransZ );
441 
442  //== Save results to return dict
443  retDict[ pybind11::handle ( pybind11::str ( "centreOfRotation" ).ptr ( ) ) ] = rotCen;
444  retDict[ pybind11::handle ( pybind11::str ( "rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
445 
446  //== Done
447  return ( retDict );
448  }, "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" ), pybind11::arg ( "eulA" ), pybind11::arg ( "eulB" ), pybind11::arg ( "eulG" ) )
449  .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" ) )
450 
451  //============================================ Internal arrays access functions
452  .def ( "findSHIndex",
453  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed shell, proshade_signed band, proshade_signed order ) -> proshade_signed
454  {
455  //== Get value
456  proshade_signed index = seanindex ( static_cast< int > ( order ),
457  static_cast< int > ( band ),
458  static_cast< int > ( self.spheres[shell]->getLocalBandwidth() ) );
459 
460  //== Done
461  return ( index );
462  }, "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." )
463  .def ( "getSphericalHarmonics",
464  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex<proshade_double> >
465  {
466  //== Allocate memory for the numpy values
467  std::complex<proshade_double>* npVals = new std::complex<proshade_double>[static_cast<proshade_unsign> ( self.noSpheres * pow( self.maxShellBand, 2.0 ) )];
468  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
469 
470  //== Fill with zeroes as the matrix will be sparse
471  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 ); }
472 
473  //== Initialise variables
474  proshade_signed pyPosSH;
475  proshade_signed pyPos;
476 
477  //== Copy data to new memory
478  for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> ( self.noSpheres ); shIt++ )
479  {
480  for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> ( self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
481  {
482  for ( proshade_signed order = -bnd; order <= bnd; order++ )
483  {
484  pyPosSH = ( static_cast< proshade_signed > ( shIt ) * static_cast< proshade_signed > ( std::pow ( self.maxShellBand, 2 ) ) );
485  pyPos = seanindex ( static_cast< int > ( order ),
486  static_cast< int > ( bnd ),
487  static_cast< int > ( self.spheres[shIt]->getLocalBandwidth() ) );
488  npVals[pyPosSH+pyPos].real ( self.sphericalHarmonics[shIt][pyPos][0] );
489  npVals[pyPosSH+pyPos].imag ( self.sphericalHarmonics[shIt][pyPos][1] );
490  }
491  }
492  }
493 
494  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
495  pybind11::capsule pyCapsuleSHs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
496 
497  //== Copy the value
498  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
499  { static_cast<int> ( self.noSpheres ), static_cast<int> ( pow ( self.maxShellBand, 2.0 ) ) },
500  { sizeof ( std::complex < proshade_double > ) * static_cast< proshade_unsign > ( std::pow ( static_cast< proshade_unsign > ( self.maxShellBand ), 2 ) ), sizeof ( std::complex < proshade_double > ) },
501  npVals,
502  pyCapsuleSHs );
503 
504  //== Done
505  return ( retArr );
506  }, "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 maxShellBand**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." )
507  .def ( "getEMatrix",
508  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
509  {
510  //== Allocate memory for the numpy values
511  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) )];
512  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
513 
514  //== Allocate local variables
515  proshade_signed index = 0;
516 
517  //== Fill with zeroes as the matrix will be sparse
518  for ( proshade_unsign iter = 0; iter < static_cast < proshade_unsign > ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
519 
520  //== Copy data to new memory
521  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxShellBand ); bandIter++ )
522  {
523  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
524  {
525  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
526  {
527  index = order2 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * bandIter );
528  npVals[index].real ( self.eMatrices[bandIter][order1][order2][0] );
529  npVals[index].imag ( self.eMatrices[bandIter][order1][order2][1] );
530  }
531  }
532  }
533 
534  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
535  pybind11::capsule pyCapsuleEMs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
536 
537  //== Create the output object
538  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
539  { self.maxShellBand, ( ( self.maxShellBand * 2 ) + 1 ), ( ( self.maxShellBand * 2 ) + 1 ) },
540  { sizeof ( std::complex < proshade_double > ) * static_cast< proshade_unsign > ( ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 ) * ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 ) ),
541  sizeof ( std::complex < proshade_double > ) * ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1,
542  sizeof ( std::complex < proshade_double > ) },
543  npVals,
544  pyCapsuleEMs );
545 
546  //== Done
547  return ( retArr );
548  }, "This function returns the 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 simply work, the order indexing starts from 0 - i.e. array[1][0][0] means band 1 ; order1 = -1 and order2 = -1." )
549  .def ( "getSO3Coefficients",
550  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
551  {
552  //== Allocate memory for the numpy values
553  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) )];
554  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
555 
556  //== Allocate local variables
557  proshade_signed index = 0;
558 
559  //== Fill with zeroes as the matrix will be sparse
560  for ( proshade_unsign iter = 0; iter < static_cast < proshade_unsign > ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) ); iter++ ) { npVals[iter].real ( 0.0 ); npVals[iter].imag ( 0.0 ); }
561 
562  //== Copy data to new memory
563  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxShellBand ); bandIter++ )
564  {
565  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
566  {
567  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
568  {
569  index = order2 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * bandIter );
570  npVals[index].real ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][0] );
571  npVals[index].imag ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][1] );
572  }
573  }
574  }
575 
576  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
577  pybind11::capsule pyCapsuleSOCoeffs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
578 
579  //== Create the output object
580  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
581  { self.maxShellBand, ( ( self.maxShellBand * 2 ) + 1 ), ( ( self.maxShellBand * 2 ) + 1 ) },
582  { sizeof ( std::complex < proshade_double > ) * ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 ) * ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 ),
583  sizeof ( std::complex < proshade_double > ) * ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1,
584  sizeof ( std::complex < proshade_double > ) },
585  npVals,
586  pyCapsuleSOCoeffs );
587 
588  //== Done
589  return ( retArr );
590  }, "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." )
591  .def ( "getRotationFunctionMap",
592  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
593  {
594  //== Allocate memory for the numpy values
595  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) )];
596  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
597 
598  //== Copy the values
599  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) ); iter++ ) { npVals[iter].real ( self.so3CoeffsInverse[iter][0] ); npVals[iter].imag ( self.so3CoeffsInverse[iter][1] ); }
600 
601  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
602  pybind11::capsule pyCapsuleRotMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
603 
604  //== Create python readable object with the correct memory access
605  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
606  { ( self.maxShellBand * 2 ), ( self.maxShellBand * 2 ), ( self.maxShellBand * 2 ) }, // Shape
607  { ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) * sizeof(std::complex < proshade_double >),
608  ( self.maxShellBand * 2 ) * sizeof(std::complex < proshade_double >),
609  sizeof(std::complex < proshade_double >) }, // C-stype strides
610  npVals, // Data
611  pyCapsuleRotMap ); // Capsule (destructor)
612 
613  //== Done
614  return ( retArr );
615  }, "This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
616  .def ( "getRotationMatrixFromSOFTCoordinates",
617  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed xPos, proshade_signed yPos, proshade_signed zPos ) -> pybind11::array_t < proshade_double >
618  {
619  //== Allocate memory for the numpy values
620  proshade_double* npVals = new proshade_double[9];
621  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
622 
623  //== Initialise local variables
624  proshade_double eulA, eulB, eulG;
625 
626  //== Compute the Euler angles from SOFT position
627  ProSHADE_internal_maths::getEulerZXZFromSOFTPosition ( static_cast< proshade_signed > ( self.maxShellBand ), xPos, yPos, zPos, &eulA, &eulB, &eulG );
628 
629  //== Compute the rotation matrix
631 
632  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
633  pybind11::capsule pyCapsuleRMSoft ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
634 
635  //== Create python readable object with the correct memory access
636  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
637  { 3 * sizeof(proshade_double), sizeof(proshade_double) },
638  npVals,
639  pyCapsuleRMSoft );
640 
641  //== Done
642  return ( retArr );
643  }, "This function converts a given rotation function map position onto the corresponding rotation matrix.", pybind11::arg ( "xPos" ), pybind11::arg ( "yPos" ), pybind11::arg ( "zPos" ) )
644  .def ( "getTranslationFunctionMap",
645  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
646  {
647  //== Allocate memory for the numpy values
648  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.getXDim() * self.getYDim() * self.getZDim() )];
649  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
650 
651  //== Copy the values
652  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] ); }
653 
654  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
655  pybind11::capsule pyCapsuleTrsMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
656 
657  //== Create python readable object with the correct memory access
658  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
659  { self.getXDim(), self.getYDim(), self.getZDim() }, // Shape
660  { self.getYDim() * self.getZDim() * sizeof(std::complex < proshade_double >),
661  self.getZDim() * sizeof(std::complex < proshade_double >),
662  sizeof(std::complex < proshade_double >) }, // C-stype strides
663  npVals, // Data
664  pyCapsuleTrsMap ); // Capsule (destructor)
665 
666  //== Done
667  return ( retArr );
668  }, "This function returns the translation function as a three-dimensional map of complex numbers." )
669 
670  //============================================ Member variables
671  .def_readwrite ( "fileName", &ProSHADE_internal_data::ProSHADE_data::fileName )
672  .def_readwrite ( "xDimSize", &ProSHADE_internal_data::ProSHADE_data::xDimSize )
673  .def_readwrite ( "yDimSize", &ProSHADE_internal_data::ProSHADE_data::yDimSize )
674  .def_readwrite ( "zDimSize", &ProSHADE_internal_data::ProSHADE_data::zDimSize )
675  .def_readwrite ( "aAngle", &ProSHADE_internal_data::ProSHADE_data::aAngle )
676  .def_readwrite ( "bAngle", &ProSHADE_internal_data::ProSHADE_data::bAngle )
677  .def_readwrite ( "cAngle", &ProSHADE_internal_data::ProSHADE_data::cAngle )
678  .def_readwrite ( "xDimIndices", &ProSHADE_internal_data::ProSHADE_data::xDimIndices )
679  .def_readwrite ( "yDimIndices", &ProSHADE_internal_data::ProSHADE_data::yDimIndices )
680  .def_readwrite ( "zDimIndices", &ProSHADE_internal_data::ProSHADE_data::zDimIndices )
681  .def_readwrite ( "xGridIndices", &ProSHADE_internal_data::ProSHADE_data::xGridIndices )
682  .def_readwrite ( "yGridIndices", &ProSHADE_internal_data::ProSHADE_data::yGridIndices )
683  .def_readwrite ( "zGridIndices", &ProSHADE_internal_data::ProSHADE_data::zGridIndices )
684  .def_readwrite ( "xAxisOrder", &ProSHADE_internal_data::ProSHADE_data::xAxisOrder )
685  .def_readwrite ( "yAxisOrder", &ProSHADE_internal_data::ProSHADE_data::yAxisOrder )
686  .def_readwrite ( "zAxisOrder", &ProSHADE_internal_data::ProSHADE_data::zAxisOrder )
687  .def_readwrite ( "xAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::xAxisOrigin )
688  .def_readwrite ( "yAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::yAxisOrigin )
689  .def_readwrite ( "zAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::zAxisOrigin )
690  .def_readwrite ( "xCom", &ProSHADE_internal_data::ProSHADE_data::xCom )
691  .def_readwrite ( "yCom", &ProSHADE_internal_data::ProSHADE_data::yCom )
692  .def_readwrite ( "zCom", &ProSHADE_internal_data::ProSHADE_data::zCom )
693  .def_readwrite ( "xFrom", &ProSHADE_internal_data::ProSHADE_data::xFrom )
694  .def_readwrite ( "yFrom", &ProSHADE_internal_data::ProSHADE_data::yFrom )
695  .def_readwrite ( "zFrom", &ProSHADE_internal_data::ProSHADE_data::zFrom )
696  .def_readwrite ( "xTo", &ProSHADE_internal_data::ProSHADE_data::xTo )
697  .def_readwrite ( "yTo", &ProSHADE_internal_data::ProSHADE_data::yTo )
698  .def_readwrite ( "zTo", &ProSHADE_internal_data::ProSHADE_data::zTo )
699  .def_readwrite ( "spherePos", &ProSHADE_internal_data::ProSHADE_data::spherePos )
700  .def_readwrite ( "noSpheres", &ProSHADE_internal_data::ProSHADE_data::noSpheres )
701  .def_readwrite ( "maxShellBand", &ProSHADE_internal_data::ProSHADE_data::maxShellBand )
702  .def_readwrite ( "isEmpty", &ProSHADE_internal_data::ProSHADE_data::isEmpty )
703  .def_readwrite ( "inputOrder", &ProSHADE_internal_data::ProSHADE_data::inputOrder )
704 
705  //============================================ Description
706  .def ( "__repr__", [] ( ) { return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );
707 
708  //================================================ Export extra symmetry elements functions
709  pyProSHADE.def ( "computeGroupElementsForGroup",
710  [] ( proshade_double x, proshade_double y, proshade_double z, proshade_unsign fold ) -> pybind11::array_t < proshade_double >
711  {
712  //== Get the results
713  std::vector<std::vector< proshade_double > > retVec = ProSHADE_internal_data::computeGroupElementsForGroup ( x, y, z, static_cast< proshade_signed > ( fold ) );
714 
715  //== Allocate memory for the numpy values
716  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
717  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
718 
719  //== Copy the values
720  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); } }
721 
722  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
723  pybind11::capsule pyCapsuleGrEls ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
724 
725  //== Create python readable object with the correct memory access
726  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
727  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
728  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
729  npVals, // Data
730  pyCapsuleGrEls ); // Capsule (destructor)
731 
732  //== Done
733  return ( retArr );
734  }, "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" ) );
735  pyProSHADE.def ( "joinElementsFromDifferentGroups",
736  [] ( pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > first,
737  pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > second,
738  proshade_double matrixTolerance,
739  bool combine ) -> pybind11::array_t < proshade_double >
740  {
741  //== Get array buffers
742  pybind11::buffer_info buf1 = first.request();
743  pybind11::buffer_info buf2 = second.request();
744 
745  //== Sanity check
746  if ( buf1.ndim != 3 || buf2.ndim != 3 )
747  {
748  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;
749  return ( pybind11::array_t < proshade_double > () );
750  }
751 
752  //== Convert input to C++ vectors
753  std::vector< std::vector< proshade_double > > fVec, sVec;
754 
755  proshade_double* dataPtr1 = reinterpret_cast < proshade_double* > ( buf1.ptr );
756  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf1.shape.at(0) ); elIt++ )
757  {
758  std::vector< proshade_double > rotMat;
759  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf1.shape.at(1) ); rowIt++ )
760  {
761  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf1.shape.at(2) ); colIt++ )
762  {
763  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 )] );
764  }
765  }
767  }
768 
769  proshade_double* dataPtr2 = reinterpret_cast < proshade_double* > ( buf2.ptr );
770  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf2.shape.at(0) ); elIt++ )
771  {
772  std::vector< proshade_double > rotMat;
773  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf2.shape.at(1) ); rowIt++ )
774  {
775  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf2.shape.at(2) ); colIt++ )
776  {
777  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 )] );
778  }
779  }
781  }
782 
783 
784  //== Get the results
785  std::vector< std::vector< proshade_double > > retVec = ProSHADE_internal_data::joinElementsFromDifferentGroups ( &fVec, &sVec, matrixTolerance, combine );
786 
787  //== Allocate memory for the numpy values
788  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
789  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
790 
791  //== Copy the values
792  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); } }
793 
794  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
795  pybind11::capsule pyCapsuleGrElsCombo ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
796 
797  //== Create python readable object with the correct memory access
798  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
799  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
800  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
801  npVals, // Data
802  pyCapsuleGrElsCombo ); // Capsule (destructor)
803 
804  //== Done
805  return ( retArr );
806  }, "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" ) );
807 }
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:1139
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:1756
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:3929
ProSHADE_internal_data::ProSHADE_data::rotateMapRealSpaceInPlace
void 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:982
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:3133
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::getRecommendedSymmetryFold
proshade_unsign getRecommendedSymmetryFold(ProSHADE_settings *settings)
This function simply returns the detected recommended symmetry fold.
Definition: ProSHADE_data.cpp:4405
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:90
ProSHADE_settings::allDetectedDAxes
std::vector< std::vector< proshade_unsign > > allDetectedDAxes
The vector of all detected dihedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:146
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:3909
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:3919
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_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:1655
ProSHADE_internal_data::ProSHADE_data::invertMirrorMap
void invertMirrorMap(ProSHADE_settings *settings)
Function for inverting the map to its mirror image.
Definition: ProSHADE_data.cpp:1101
ProSHADE_internal_maths::getRotationMatrixFromEulerZXZAngles
void getRotationMatrixFromEulerZXZAngles(proshade_double eulerAlpha, proshade_double eulerBeta, proshade_double eulerGamma, proshade_double *matrix)
Function to find the rotation matrix from Euler angles (ZXZ convention).
Definition: ProSHADE_maths.cpp:1011
ProSHADE_settings::allDetectedTAxes
std::vector< proshade_unsign > allDetectedTAxes
The vector of all detected tetrahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:147
ProSHADE_settings::allDetectedOAxes
std::vector< proshade_unsign > allDetectedOAxes
The vector of all detected octahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:148
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:1438
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::getXDim
proshade_unsign getXDim(void)
This function allows access to the map size in indices along the X axis.
Definition: ProSHADE_data.cpp:3939
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:756
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:1712
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:398
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::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:948
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:1040
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:3949
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:1202
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:1082
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:1612
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::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:3959
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:91
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, bool firstModel=true)
This function writes out the PDB formatted file coresponding to the structure so that its COM is at s...
Definition: ProSHADE_data.cpp:1013
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:1916
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:1366
ProSHADE_internal_misc::checkMemoryAllocation
void checkMemoryAllocation(chVar checkVar, std::string fileP, unsigned int lineP, std::string funcP, std::string infoP="This error may occurs when ProSHADE requests memory to be\n : allocated to it and this operation fails. This could\n : happen when not enough memory is available, either due to\n : other processes using a lot of memory, or when the machine\n : does not have sufficient memory available. Re-run to see\n : if this problem persists.")
Checks if memory was allocated properly.
Definition: ProSHADE_misc.hpp:67
ProSHADE_internal_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_maths::getEulerZXZFromSOFTPosition
void getEulerZXZFromSOFTPosition(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 (ZXZ convention) from index position in the inverse SOFT map.
Definition: ProSHADE_maths.cpp:963
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_settings::allDetectedIAxes
std::vector< proshade_unsign > allDetectedIAxes
The vector of all detected icosahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:149
ProSHADE_settings::forceBounds
proshade_signed * forceBounds
These will be the boundaries to be forced upon the map.
Definition: ProSHADE_settings.hpp:93
ProSHADE_settings::detectedSymmetry
std::vector< proshade_double * > detectedSymmetry
The vector of detected symmetry axes.
Definition: ProSHADE_settings.hpp:144
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:92
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:2980
ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryType
std::string getRecommendedSymmetryType(ProSHADE_settings *settings)
This function simply returns the detected recommended symmetry type.
Definition: ProSHADE_data.cpp:4394
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_settings::usePeakSearchInRotationFunctionSpace
bool usePeakSearchInRotationFunctionSpace
This variable switch decides whether symmetry detection will be done using peak search in rotation fu...
Definition: ProSHADE_settings.hpp:129
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:1155
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_settings::allDetectedCAxes
std::vector< std::vector< proshade_double > > allDetectedCAxes
The vector of all detected cyclic symmetry axes.
Definition: ProSHADE_settings.hpp:145
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:3676
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:618
ProSHADE_internal_data::ProSHADE_data::readInStructure
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
Definition: ProSHADE_data.cpp:501
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:1503