ProSHADE  0.7.6.2 (DEC 2021)
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", static_cast<void (ProSHADE_internal_data::ProSHADE_data::*)(gemmi::Structure, proshade_unsign, ProSHADE_settings*)> (&ProSHADE_internal_data::ProSHADE_data::readInStructure), "This function initialises the basic ProSHADE_data variables and reads in a single structure from Gemmi co-ordinate object.", pybind11::arg ( "gemmiStruct" ), pybind11::arg ( "inputO" ), pybind11::arg ( "settings" ) )
91 
92  .def ( "readInStructure",
93  [] ( 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 )
94  {
95  //== Sanity check
96  pybind11::buffer_info maskArr_buf = maskArr.request();
97  pybind11::buffer_info weightsArr_buf = weightsArr.request();
98 
99  //== Check for mask
100  if ( ( maskArr_buf.shape.at(0) != 0 ) && ( maskArr_buf.shape.at(1) != 0 ) && ( maskArr_buf.shape.at(2) != 0 ) )
101  {
102  //== Is number of dimensions correct?
103  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 ); }
104 
105  //== Mask was given! Create the array
106  proshade_double* mskArr = new proshade_double[maskArr_buf.shape.at(0)*maskArr_buf.shape.at(1)*maskArr_buf.shape.at(2)];
107  ProSHADE_internal_misc::checkMemoryAllocation ( mskArr, __FILE__, __LINE__, __func__ );
108 
109  //== Copy values
110  proshade_double* arrStart = static_cast< proshade_double* > ( maskArr_buf.ptr );
111  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]; }
112 
113  //== Check for weights
114  if ( ( weightsArr_buf.shape.at(0) != 0 ) && ( weightsArr_buf.shape.at(1) != 0 ) && ( weightsArr_buf.shape.at(2) != 0 ) )
115  {
116  //== Is number of dimensions correct?
117  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 ); }
118 
119  //== Weights were given! Create array
120  proshade_double* wghArr = new proshade_double[weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2)];
121  ProSHADE_internal_misc::checkMemoryAllocation ( wghArr, __FILE__, __LINE__, __func__ );
122 
123  //== Copy values
124  proshade_double* arr2Start = static_cast< proshade_double* > ( weightsArr_buf.ptr );
125  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]; }
126 
127  //== Call C++ function with both mask and weights
128  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) ) );
129 
130  //== Release memory
131  delete[] wghArr;
132  }
133  else
134  {
135  //== No weights given. Call C++ function with only mask
136  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 );
137  }
138 
139  //== Release memory
140  delete[] mskArr;
141  }
142  else
143  {
144  //== Mask was empty! Check for weights
145  if ( ( weightsArr_buf.shape.at(0) != 0 ) && ( weightsArr_buf.shape.at(1) != 0 ) && ( weightsArr_buf.shape.at(2) != 0 ) )
146  {
147  //== Is number of dimensions correct?
148  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 ); }
149 
150  //== Weights were given! Create array
151  proshade_double* wghArr = new proshade_double[weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2)];
152  ProSHADE_internal_misc::checkMemoryAllocation ( wghArr, __FILE__, __LINE__, __func__ );
153 
154  //== Copy values
155  proshade_double* arr2Start = static_cast< proshade_double* > ( weightsArr_buf.ptr );
156  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]; }
157 
158  //== Call C++ function with only weights
159  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) ) );
160 
161  //== Release memory
162  delete[] wghArr;
163  }
164  else
165  {
166  //== No weights either. Call C++ function with no mask and no weights
167  self.readInStructure ( fName, inputO, settings, nullptr, 0, 0, 0, nullptr, 0, 0, 0 );
168  }
169  }
170 
171  //== Done
172  return ;
173  }, "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 > () )
174  .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 )
175  .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 )
176  .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 )
177  .def ( "getMap",
178  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < proshade_double >
179  {
180  //== Copy the values
181  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { self.xDimIndices, self.yDimIndices, self.zDimIndices }, // Shape
182  { self.yDimIndices * self.zDimIndices * sizeof(proshade_double),
183  self.zDimIndices * sizeof(proshade_double),
184  sizeof(proshade_double) }, // C-stype strides
185  self.internalMap ); // Data
186 
187  //== Done
188  return ( retArr );
189  } )
190 
191  //============================================ Data processing functions
192  .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" ) )
193  .def ( "invertMirrorMap", &ProSHADE_internal_data::ProSHADE_data::invertMirrorMap, "Function for inverting the map to its mirror image.", pybind11::arg ( "settings" ) )
194  .def ( "normaliseMap", &ProSHADE_internal_data::ProSHADE_data::normaliseMap, "Function for normalising the map values to mean 0 and sd 1.", pybind11::arg ( "settings" ) )
195  .def ( "maskMap", &ProSHADE_internal_data::ProSHADE_data::maskMap, "Function for computing the map mask using blurring and X IQRs from median.", pybind11::arg ( "settings" ) )
196  .def ( "reSampleMap", &ProSHADE_internal_data::ProSHADE_data::reSampleMap, "This function changes the internal map sampling to conform to particular resolution value.", pybind11::arg ( "settings" ) )
197  .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" ) )
198  .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" ) )
199  .def ( "removePhaseInormation", &ProSHADE_internal_data::ProSHADE_data::removePhaseInormation, "This function removes phase from the map, effectively converting it to Patterson map.", pybind11::arg ( "settings" ) )
200  .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" ) )
201  .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" ) )
202  .def ( "getReBoxBoundaries",
203  [] ( ProSHADE_internal_data::ProSHADE_data &self,ProSHADE_settings* settings ) -> pybind11::array_t < proshade_signed >
204  {
205  //== Allocate output memory
206  proshade_signed* retVals = new proshade_signed[6];
207  ProSHADE_internal_misc::checkMemoryAllocation ( retVals, __FILE__, __LINE__, __func__ );
208 
209  //== If same bounds as first one are required, test if possible and return these instead
210  if ( settings->useSameBounds && ( self.inputOrder != 0 ) )
211  {
212  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->forceBounds[iter]; }
213  }
214  //== In this case, bounds need to be found de novo
215  else
216  {
217  //== Find the non-zero bounds
219  static_cast< proshade_signed > ( self.xDimIndices ),
220  static_cast< proshade_signed > ( self.yDimIndices ),
221  static_cast< proshade_signed > ( self.zDimIndices ),
222  retVals );
223 
224  //== Add the extra space
225  ProSHADE_internal_mapManip::addExtraBoundSpace ( self.xDimIndices, self.yDimIndices, self.zDimIndices,
226  self.xDimSize, self.yDimSize, self.zDimSize, retVals, settings->boundsExtraSpace );
227 
228  //== Beautify boundaries
229  ProSHADE_internal_mapManip::beautifyBoundaries ( retVals, self.xDimIndices, self.yDimIndices, self.zDimIndices, settings->boundsSimilarityThreshold );
230 
231  //== If need be, save boundaries to be used for all other structure
232  if ( settings->useSameBounds && ( self.inputOrder == 0 ) ) { for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->forceBounds[iter] = retVals[iter]; } }
233  }
234 
235  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
236  pybind11::capsule pyCapsuleReBox2 ( retVals, []( void *f ) { proshade_signed* foo = reinterpret_cast< proshade_signed* > ( f ); delete foo; } );
237 
238  //== Copy the value
239  pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 }, // Shape
240  { sizeof(proshade_signed) }, // C-stype strides
241  retVals, // Data
242  pyCapsuleReBox2 ); // Capsule
243 
244  //== Done
245  return ( retArr );
246  }, "This function finds the boundaries enclosing positive map values and adds some extra space." )
247  .def ( "createNewMapFromBounds",
248  [] ( ProSHADE_internal_data::ProSHADE_data &self, pybind11::array_t < proshade_signed > bounds, ProSHADE_internal_data::ProSHADE_data* newStr, ProSHADE_settings* settings )
249  {
250  //== Allocate memory for bounds copy
251  proshade_signed* newBounds = new proshade_signed[6];
252  ProSHADE_internal_misc::checkMemoryAllocation ( newBounds, __FILE__, __LINE__, __func__ );
253 
254  //== Copy to the pointer
255  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
256 
257  //== Run C++ function
258  self.createNewMapFromBounds ( settings, newStr, newBounds );
259 
260  //== Done
261  return ;
262  }, "This function creates a new structure from the calling structure and new bounds values." )
263 
264  //============================================ Data sphere mapping functions
265  .def ( "mapToSpheres", &ProSHADE_internal_data::ProSHADE_data::mapToSpheres, "This function converts the internal map onto a set of concentric spheres.", pybind11::arg ( "settings" ) )
266  .def ( "getSpherePositions", &ProSHADE_internal_data::ProSHADE_data::getSpherePositions, "This function determines the sphere positions (radii) for sphere mapping.", pybind11::arg ( "settings" ) )
267  .def ( "computeSphericalHarmonics", &ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics, "This function computes the spherical harmonics decomposition for the whole structure.", pybind11::arg ( "settings" ) )
268 
269  //============================================ Accessor functions
270  .def ( "getXDimSize", &ProSHADE_internal_data::ProSHADE_data::getXDimSize, "This function allows access to the map size in angstroms along the X axis." )
271  .def ( "getYDimSize", &ProSHADE_internal_data::ProSHADE_data::getYDimSize, "This function allows access to the map size in angstroms along the Y axis." )
272  .def ( "getZDimSize", &ProSHADE_internal_data::ProSHADE_data::getZDimSize, "This function allows access to the map size in angstroms along the Z axis." )
273  .def ( "getXDim", &ProSHADE_internal_data::ProSHADE_data::getXDim, "This function allows access to the map size in indices along the X axis." )
274  .def ( "getYDim", &ProSHADE_internal_data::ProSHADE_data::getYDim, "This function allows access to the map size in indices along the Y axis." )
275  .def ( "getZDim", &ProSHADE_internal_data::ProSHADE_data::getZDim, "This function allows access to the map size in indices along the Z axis." )
276 
277  //============================================ Symmetry related functions
278  .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" ) )
279  .def ( "detectSymmetryInStructure",
281  {
282  //== Call the appropriate C++ function
283  self.detectSymmetryFromAngleAxisSpace ( settings, &settings->detectedSymmetry, &settings->allDetectedCAxes );
284  }, "This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg ( "settings" ) )
285  .def ( "getRecommendedSymmetryType", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryType, "This function simply returns the detected recommended symmetry type.", pybind11::arg ( "settings" ) )
286  .def ( "getRecommendedSymmetryFold", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryFold, "This function simply returns the detected recommended symmetry fold.", pybind11::arg ( "settings" ) )
287  .def ( "getRecommendedSymmetryAxes",
288  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
289  {
290  //== Allocate memory for the numpy values
291  float* npVals = new float[static_cast<unsigned int> ( settings->detectedSymmetry.size() * 7 )];
292  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
293 
294  //== Copy values
295  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] ); } }
296 
297  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
298  pybind11::capsule pyCapsuleStrRecSym ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
299 
300  //== Copy the value
301  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( settings->detectedSymmetry.size() ), static_cast<int> ( 7 ) }, // Shape
302  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
303  npVals, // Data
304  pyCapsuleStrRecSym ); // Capsule
305 
306  //== Done
307  return ( retArr );
308  }, "This function returns the recommended symmetry axes as a 2D numpy array." )
309  .def ( "getAllCSyms",
310  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
311  {
312  //== Allocate memory for the numpy values
313  float* npVals = new float[static_cast<unsigned int> ( settings->allDetectedCAxes.size() * 7 )];
314  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
315 
316  //== Copy values
317  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) ); } }
318 
319  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
320  pybind11::capsule pyCapsuleStrSymList ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
321 
322  //== Copy the value
323  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( settings->allDetectedCAxes.size() ), 7 }, // Shape
324  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
325  npVals, // Data
326  pyCapsuleStrSymList ); // Capsule
327 
328  //== Done
329  return ( retArr );
330  }, "This function returns all symmetry axes as a 2D numpy array." )
331  .def ( "getNonCSymmetryAxesIndices",
332  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::dict
333  {
334  //== Initialise local variables
335  pybind11::dict retDict;
336  pybind11::list dList, tList, oList, iList;
337 
338  //== Fill in the D list
339  for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( settings->allDetectedDAxes.size() ); dIt++ )
340  {
341  pybind11::list memList;
342  for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( settings->allDetectedDAxes.at(dIt).size() ); memIt++ )
343  {
344  memList.append ( settings->allDetectedDAxes.at(dIt).at(memIt) );
345  }
346  dList.append ( memList );
347  }
348 
349  //== Fill in the T list
350  for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( settings->allDetectedTAxes.size() ); tIt++ )
351  {
352  tList.append ( settings->allDetectedTAxes.at ( tIt ) );
353  }
354 
355  //== Fill in the O list
356  for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( settings->allDetectedOAxes.size() ); oIt++ )
357  {
358  oList.append ( settings->allDetectedOAxes.at ( oIt ) );
359  }
360 
361  //== Fill in the T list
362  for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( settings->allDetectedIAxes.size() ); iIt++ )
363  {
364  iList.append ( settings->allDetectedIAxes.at ( iIt ) );
365  }
366 
367  //== Save the lists to the dict
368  retDict[ pybind11::handle ( pybind11::str ( "D" ).ptr ( ) ) ] = dList;
369  retDict[ pybind11::handle ( pybind11::str ( "T" ).ptr ( ) ) ] = tList;
370  retDict[ pybind11::handle ( pybind11::str ( "O" ).ptr ( ) ) ] = oList;
371  retDict[ pybind11::handle ( pybind11::str ( "I" ).ptr ( ) ) ] = iList;
372 
373  //== Done
374  return ( retDict );
375  }, "This function returns array of non-C axes indices." )
376  .def ( "getAllGroupElements",
377  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings, pybind11::array_t < int > axList, std::string groupType, proshade_double matrixTolerance ) -> pybind11::list
378  {
379  //== Sanity check
380  pybind11::buffer_info buf = axList.request();
381  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 ); }
382 
383  //== Convert to vector of unsigns
384  std::vector< proshade_unsign > axesList;
385  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) ) ); }
386 
387  //== Get the results
388  std::vector < std::vector< proshade_double > > vals = self.getAllGroupElements ( settings, axesList, groupType, matrixTolerance );
389 
390  //== Initialise return variable
391  pybind11::list retList;
392 
393  //== Copy values
394  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
395  {
396  //== Allocate memory for the numpy values
397  float* npVals = new float[static_cast<unsigned int> ( 9 )];
398  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
399 
400  //== Copy values to memory
401  for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] = static_cast< float > ( vals.at(iter).at(it) ); }
402 
403  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
404  pybind11::capsule pyCapsuleGrpEl ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
405 
406  //== Copy the value
407  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 }, // Shape
408  { 3 * sizeof(float), sizeof(float) }, // C-stype strides
409  npVals, // Data
410  pyCapsuleGrpEl ); // Capsule
411 
412  //== Save matrix
413  retList.append ( retArr );
414  }
415 
416  //== Done
417  return ( retList );
418  }, "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 )
419 
420  .def ( "getMapCOMProcessChange",
421  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < float >
422  {
423  //== Get the values
424  std::vector< proshade_double > vals = self.getMapCOMProcessChange ();
425 
426  //== Allocate memory for the numpy values
427  float* npVals = new float[static_cast<unsigned int> ( 3 )];
428  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
429 
430  //== Copy values
431  for ( proshade_unsign iter = 0; iter < 3; iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
432 
433  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
434  pybind11::capsule pyCapsuleSymShiftDat ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
435 
436  //== Copy the value
437  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( vals.size() ) }, // Shape
438  { sizeof(float) }, // C-stype strides
439  npVals, // Data
440  pyCapsuleSymShiftDat ); // Capsule
441 
442  //== Done
443  return ( retArr );
444  }, "This function returns the shift in Angstrom applied to the internal map representation in order to align its COM with the centre of box." )
445 
446  //============================================ Overlay related functions
447  .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" ) )
448  .def ( "getBestRotationMapPeaksEulerAngles",
449  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
450  {
451  //== Get values
452  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
453 
454  //== Allocate memory for the numpy values
455  float* npVals = new float[static_cast<unsigned int> (vals.size())];
456  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
457 
458  //== Copy values
459  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
460 
461  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
462  pybind11::capsule pyCapsuleEuPeak ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
463 
464  //== Copy the value
465  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<unsigned int> (vals.size()) }, // Shape
466  { sizeof(float) }, // C-stype strides
467  npVals, // Data
468  pyCapsuleEuPeak ); // Capsule (C++ destructor, basically)
469 
470  //== Done
471  return ( retArr );
472  }, "This function returns a vector of three floats, the three Euler angles of the best peak in the rotation map.", pybind11::arg ( "settings" ) )
473  .def ( "getBestRotationMapPeaksRotationMatrix",
474  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < proshade_double >
475  {
476  //== Get values
477  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
478 
479  //== Convert Euler ZXZ to matrix
480  proshade_double* retMat = new proshade_double[9];
481  ProSHADE_internal_misc::checkMemoryAllocation ( retMat, __FILE__, __LINE__, __func__ );
482  ProSHADE_internal_maths::getRotationMatrixFromEulerZXZAngles ( vals.at(0), vals.at(1), vals.at(2), retMat );
483 
484  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
485  pybind11::capsule pyCapsuleRMPeak ( retMat, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
486 
487  //== Copy the value
488  pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 }, // Shape
489  { 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
490  retMat, // Data
491  pyCapsuleRMPeak ); // Capsule
492 
493  //== Done
494  return ( retArr );
495  }, "This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg ( "settings" ) )
496  .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" ) )
497  .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" ) )
498  .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" ) )
499  .def ( "computeTranslationMap", &ProSHADE_internal_data::ProSHADE_data::computeTranslationMap, "This function does the computation of the translation map and saves results internally.", pybind11::arg ( "staticStructure" ) )
500  .def ( "getOverlayTranslations",
501  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_internal_data::ProSHADE_data* staticStructure ) -> pybind11::dict
502  {
503  //== Get values
504  std::vector< proshade_double > vals = self.getBestTranslationMapPeaksAngstrom ( staticStructure );
505 
506  //== Initialise variables
507  pybind11::dict retDict;
508  pybind11::list rotCen, toOverlay;
509 
510  //== Copy values
511  rotCen.append ( self.originalPdbRotCenX );
512  rotCen.append ( self.originalPdbRotCenY );
513  rotCen.append ( self.originalPdbRotCenZ );
514 
515  toOverlay.append ( self.originalPdbTransX );
516  toOverlay.append ( self.originalPdbTransY );
517  toOverlay.append ( self.originalPdbTransZ );
518 
519  //== Save results to return dict
520  retDict[ pybind11::handle ( pybind11::str ( "centreOfRotation" ).ptr ( ) ) ] = rotCen;
521  retDict[ pybind11::handle ( pybind11::str ( "rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
522 
523  //== Done
524  return ( retDict );
525  }, "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" ) )
526  .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" ) )
527 
528  //============================================ Internal arrays access functions
529  .def ( "findSHIndex",
530  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed shell, proshade_signed band, proshade_signed order ) -> proshade_signed
531  {
532  //== Get value
533  proshade_signed index = seanindex ( static_cast< int > ( order ),
534  static_cast< int > ( band ),
535  static_cast< int > ( self.spheres[shell]->getLocalBandwidth() ) );
536 
537  //== Done
538  return ( index );
539  }, "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." )
540  .def ( "getSphericalHarmonics",
541  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex<proshade_double> >
542  {
543  //== Allocate memory for the numpy values
544  std::complex<proshade_double>* npVals = new std::complex<proshade_double>[static_cast<proshade_unsign> ( self.noSpheres * pow( self.maxShellBand, 2.0 ) )];
545  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
546 
547  //== Fill with zeroes as the matrix will be sparse
548  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 ); }
549 
550  //== Initialise variables
551  proshade_signed pyPosSH;
552  proshade_signed pyPos;
553 
554  //== Copy data to new memory
555  for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> ( self.noSpheres ); shIt++ )
556  {
557  for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> ( self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
558  {
559  for ( proshade_signed order = -bnd; order <= bnd; order++ )
560  {
561  pyPosSH = ( static_cast< proshade_signed > ( shIt ) * static_cast< proshade_signed > ( std::pow ( self.maxShellBand, 2 ) ) );
562  pyPos = seanindex ( static_cast< int > ( order ),
563  static_cast< int > ( bnd ),
564  static_cast< int > ( self.spheres[shIt]->getLocalBandwidth() ) );
565  npVals[pyPosSH+pyPos].real ( self.sphericalHarmonics[shIt][pyPos][0] );
566  npVals[pyPosSH+pyPos].imag ( self.sphericalHarmonics[shIt][pyPos][1] );
567  }
568  }
569  }
570 
571  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
572  pybind11::capsule pyCapsuleSHs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
573 
574  //== Copy the value
575  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
576  { static_cast<int> ( self.noSpheres ), static_cast<int> ( pow ( self.maxShellBand, 2.0 ) ) },
577  { sizeof ( std::complex < proshade_double > ) * static_cast< proshade_unsign > ( std::pow ( static_cast< proshade_unsign > ( self.maxShellBand ), 2 ) ), sizeof ( std::complex < proshade_double > ) },
578  npVals,
579  pyCapsuleSHs );
580 
581  //== Done
582  return ( retArr );
583  }, "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." )
584  .def ( "getEMatrix",
585  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
586  {
587  //== Allocate memory for the numpy values
588  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) )];
589  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
590 
591  //== Allocate local variables
592  proshade_signed index = 0;
593 
594  //== Fill with zeroes as the matrix will be sparse
595  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 ); }
596 
597  //== Copy data to new memory
598  proshade_double emReal, emImag;
599  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxShellBand ); bandIter++ )
600  {
601  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
602  {
603  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
604  {
605  index = order2 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * bandIter );
606  self.getEMatrixValue( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( order1 ), static_cast< proshade_unsign > ( order2 ), &emReal, &emImag );
607  npVals[index].real ( emReal );
608  npVals[index].imag ( emImag );
609  }
610  }
611  }
612 
613  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
614  pybind11::capsule pyCapsuleEMs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
615 
616  //== Create the output object
617  proshade_unsign orderSize = ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 );
618  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
619  { self.maxShellBand, orderSize, orderSize },
620  { sizeof ( std::complex < proshade_double > ) * orderSize * orderSize,
621  sizeof ( std::complex < proshade_double > ) * orderSize,
622  sizeof ( std::complex < proshade_double > ) },
623  npVals,
624  pyCapsuleEMs );
625 
626  //== Done
627  return ( retArr );
628  }, "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." )
629  .def ( "getSO3Coefficients",
630  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
631  {
632  //== Allocate memory for the numpy values
633  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) )];
634  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
635 
636  //== Allocate local variables
637  proshade_signed index = 0;
638 
639  //== Fill with zeroes as the matrix will be sparse
640  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 ); }
641 
642  //== Copy data to new memory
643  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxShellBand ); bandIter++ )
644  {
645  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
646  {
647  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
648  {
649  index = order2 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * bandIter );
650  npVals[index].real ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][0] );
651  npVals[index].imag ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][1] );
652  }
653  }
654  }
655 
656  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
657  pybind11::capsule pyCapsuleSOCoeffs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
658 
659  //== Create the output object
660  proshade_unsign orderSize = ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 );
661  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
662  { self.maxShellBand, orderSize, orderSize },
663  { sizeof ( std::complex < proshade_double > ) * orderSize * orderSize,
664  sizeof ( std::complex < proshade_double > ) * orderSize,
665  sizeof ( std::complex < proshade_double > ) },
666  npVals,
667  pyCapsuleSOCoeffs );
668 
669  //== Done
670  return ( retArr );
671  }, "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." )
672  .def ( "getRotationFunctionMap",
673  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
674  {
675  //== Allocate memory for the numpy values
676  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) )];
677  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
678 
679  //== Copy the values
680  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] ); }
681 
682  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
683  pybind11::capsule pyCapsuleRotMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
684 
685  //== Create python readable object with the correct memory access
686  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
687  { ( self.maxShellBand * 2 ), ( self.maxShellBand * 2 ), ( self.maxShellBand * 2 ) }, // Shape
688  { ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) * sizeof(std::complex < proshade_double >),
689  ( self.maxShellBand * 2 ) * sizeof(std::complex < proshade_double >),
690  sizeof(std::complex < proshade_double >) }, // C-stype strides
691  npVals, // Data
692  pyCapsuleRotMap ); // Capsule (destructor)
693 
694  //== Done
695  return ( retArr );
696  }, "This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
697  .def ( "getRotationMatrixFromSOFTCoordinates",
698  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed xPos, proshade_signed yPos, proshade_signed zPos ) -> pybind11::array_t < proshade_double >
699  {
700  //== Allocate memory for the numpy values
701  proshade_double* npVals = new proshade_double[9];
702  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
703 
704  //== Initialise local variables
705  proshade_double eulA, eulB, eulG;
706 
707  //== Compute the Euler angles from SOFT position
708  ProSHADE_internal_maths::getEulerZXZFromSOFTPosition ( static_cast< proshade_signed > ( self.maxShellBand ), xPos, yPos, zPos, &eulA, &eulB, &eulG );
709 
710  //== Compute the rotation matrix
712 
713  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
714  pybind11::capsule pyCapsuleRMSoft ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
715 
716  //== Create python readable object with the correct memory access
717  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
718  { 3 * sizeof(proshade_double), sizeof(proshade_double) },
719  npVals,
720  pyCapsuleRMSoft );
721 
722  //== Done
723  return ( retArr );
724  }, "This function converts a given rotation function map position onto the corresponding rotation matrix.", pybind11::arg ( "xPos" ), pybind11::arg ( "yPos" ), pybind11::arg ( "zPos" ) )
725  .def ( "getTranslationFunctionMap",
726  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
727  {
728  //== Allocate memory for the numpy values
729  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.getXDim() * self.getYDim() * self.getZDim() )];
730  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
731 
732  //== Copy the values
733  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] ); }
734 
735  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
736  pybind11::capsule pyCapsuleTrsMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
737 
738  //== Create python readable object with the correct memory access
739  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
740  { self.getXDim(), self.getYDim(), self.getZDim() }, // Shape
741  { self.getYDim() * self.getZDim() * sizeof(std::complex < proshade_double >),
742  self.getZDim() * sizeof(std::complex < proshade_double >),
743  sizeof(std::complex < proshade_double >) }, // C-stype strides
744  npVals, // Data
745  pyCapsuleTrsMap ); // Capsule (destructor)
746 
747  //== Done
748  return ( retArr );
749  }, "This function returns the translation function as a three-dimensional map of complex numbers." )
750 
751  //============================================ Member variables
752  .def_readwrite ( "fileName", &ProSHADE_internal_data::ProSHADE_data::fileName )
753  .def_readwrite ( "xDimSize", &ProSHADE_internal_data::ProSHADE_data::xDimSize )
754  .def_readwrite ( "yDimSize", &ProSHADE_internal_data::ProSHADE_data::yDimSize )
755  .def_readwrite ( "zDimSize", &ProSHADE_internal_data::ProSHADE_data::zDimSize )
756  .def_readwrite ( "aAngle", &ProSHADE_internal_data::ProSHADE_data::aAngle )
757  .def_readwrite ( "bAngle", &ProSHADE_internal_data::ProSHADE_data::bAngle )
758  .def_readwrite ( "cAngle", &ProSHADE_internal_data::ProSHADE_data::cAngle )
759  .def_readwrite ( "xDimIndices", &ProSHADE_internal_data::ProSHADE_data::xDimIndices )
760  .def_readwrite ( "yDimIndices", &ProSHADE_internal_data::ProSHADE_data::yDimIndices )
761  .def_readwrite ( "zDimIndices", &ProSHADE_internal_data::ProSHADE_data::zDimIndices )
762  .def_readwrite ( "xGridIndices", &ProSHADE_internal_data::ProSHADE_data::xGridIndices )
763  .def_readwrite ( "yGridIndices", &ProSHADE_internal_data::ProSHADE_data::yGridIndices )
764  .def_readwrite ( "zGridIndices", &ProSHADE_internal_data::ProSHADE_data::zGridIndices )
765  .def_readwrite ( "xAxisOrder", &ProSHADE_internal_data::ProSHADE_data::xAxisOrder )
766  .def_readwrite ( "yAxisOrder", &ProSHADE_internal_data::ProSHADE_data::yAxisOrder )
767  .def_readwrite ( "zAxisOrder", &ProSHADE_internal_data::ProSHADE_data::zAxisOrder )
768  .def_readwrite ( "xAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::xAxisOrigin )
769  .def_readwrite ( "yAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::yAxisOrigin )
770  .def_readwrite ( "zAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::zAxisOrigin )
771  .def_readwrite ( "xCom", &ProSHADE_internal_data::ProSHADE_data::xCom )
772  .def_readwrite ( "yCom", &ProSHADE_internal_data::ProSHADE_data::yCom )
773  .def_readwrite ( "zCom", &ProSHADE_internal_data::ProSHADE_data::zCom )
774  .def_readwrite ( "xFrom", &ProSHADE_internal_data::ProSHADE_data::xFrom )
775  .def_readwrite ( "yFrom", &ProSHADE_internal_data::ProSHADE_data::yFrom )
776  .def_readwrite ( "zFrom", &ProSHADE_internal_data::ProSHADE_data::zFrom )
777  .def_readwrite ( "xTo", &ProSHADE_internal_data::ProSHADE_data::xTo )
778  .def_readwrite ( "yTo", &ProSHADE_internal_data::ProSHADE_data::yTo )
779  .def_readwrite ( "zTo", &ProSHADE_internal_data::ProSHADE_data::zTo )
780  .def_readwrite ( "spherePos", &ProSHADE_internal_data::ProSHADE_data::spherePos )
781  .def_readwrite ( "noSpheres", &ProSHADE_internal_data::ProSHADE_data::noSpheres )
782  .def_readwrite ( "maxShellBand", &ProSHADE_internal_data::ProSHADE_data::maxShellBand )
783  .def_readwrite ( "isEmpty", &ProSHADE_internal_data::ProSHADE_data::isEmpty )
784  .def_readwrite ( "inputOrder", &ProSHADE_internal_data::ProSHADE_data::inputOrder )
785 
786  //============================================ Description
787  .def ( "__repr__", [] ( ) { return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );
788 
789  //================================================ Export extra symmetry elements functions
790  pyProSHADE.def ( "computeGroupElementsForGroup",
791  [] ( proshade_double x, proshade_double y, proshade_double z, proshade_unsign fold ) -> pybind11::array_t < proshade_double >
792  {
793  //== Get the results
794  std::vector<std::vector< proshade_double > > retVec = ProSHADE_internal_data::computeGroupElementsForGroup ( x, y, z, static_cast< proshade_signed > ( fold ) );
795 
796  //== Allocate memory for the numpy values
797  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
798  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
799 
800  //== Copy the values
801  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); } }
802 
803  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
804  pybind11::capsule pyCapsuleGrEls ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
805 
806  //== Create python readable object with the correct memory access
807  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
808  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
809  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
810  npVals, // Data
811  pyCapsuleGrEls ); // Capsule (destructor)
812 
813  //== Done
814  return ( retArr );
815  }, "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" ) );
816  pyProSHADE.def ( "joinElementsFromDifferentGroups",
817  [] ( pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > first,
818  pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > second,
819  proshade_double matrixTolerance,
820  bool combine ) -> pybind11::array_t < proshade_double >
821  {
822  //== Get array buffers
823  pybind11::buffer_info buf1 = first.request();
824  pybind11::buffer_info buf2 = second.request();
825 
826  //== Sanity check
827  if ( buf1.ndim != 3 || buf2.ndim != 3 )
828  {
829  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;
830  return ( pybind11::array_t < proshade_double > () );
831  }
832 
833  //== Convert input to C++ vectors
834  std::vector< std::vector< proshade_double > > fVec, sVec;
835 
836  proshade_double* dataPtr1 = reinterpret_cast < proshade_double* > ( buf1.ptr );
837  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf1.shape.at(0) ); elIt++ )
838  {
839  std::vector< proshade_double > rotMat;
840  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf1.shape.at(1) ); rowIt++ )
841  {
842  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf1.shape.at(2) ); colIt++ )
843  {
844  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 )] );
845  }
846  }
848  }
849 
850  proshade_double* dataPtr2 = reinterpret_cast < proshade_double* > ( buf2.ptr );
851  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf2.shape.at(0) ); elIt++ )
852  {
853  std::vector< proshade_double > rotMat;
854  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf2.shape.at(1) ); rowIt++ )
855  {
856  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf2.shape.at(2) ); colIt++ )
857  {
858  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 )] );
859  }
860  }
862  }
863 
864 
865  //== Get the results
866  std::vector< std::vector< proshade_double > > retVec = ProSHADE_internal_data::joinElementsFromDifferentGroups ( &fVec, &sVec, matrixTolerance, combine );
867 
868  //== Allocate memory for the numpy values
869  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
870  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
871 
872  //== Copy the values
873  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); } }
874 
875  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
876  pybind11::capsule pyCapsuleGrElsCombo ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
877 
878  //== Create python readable object with the correct memory access
879  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
880  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
881  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
882  npVals, // Data
883  pyCapsuleGrElsCombo ); // Capsule (destructor)
884 
885  //== Done
886  return ( retArr );
887  }, "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" ) );
888  pyProSHADE.def ( "getRotationFunctionSpheres",
889  [] ( ProSHADE_internal_data::ProSHADE_data *dataObj, int fold ) -> pybind11::array_t < proshade_double >
890  {
891 
892  //== Check if rotation function was computed (sanity check)
893  if ( dataObj->getInvSO3Coeffs() == nullptr )
894  {
895  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;
896  return ( pybind11::array_t < proshade_double > () );
897  }
898 
899  //== Check that fold makes sense
900  if ( fold <= 1 )
901  {
902  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;
903  return ( pybind11::array_t < proshade_double > () );
904  }
905 
906  //== Compute the sphere mapping
907  std::vector<ProSHADE_internal_spheres::ProSHADE_rotFun_sphere*> sphereMappedRotFun;
908 
909  //== Convert rotation function to only the required angle-axis space spheres and find all peaks
910  proshade_double soughtAngle = 0.0;
911  for ( proshade_double angIt = 1.0; angIt < static_cast < proshade_double > ( fold ); angIt += 1.0 )
912  {
913  //== Figure the angles to form the symmetry
914  soughtAngle = angIt * ( 2.0 * M_PI / static_cast<proshade_double> ( fold ) );
915 
916  //== Create the angle-axis sphere with correct radius (angle)
917  sphereMappedRotFun.emplace_back ( new ProSHADE_internal_spheres::ProSHADE_rotFun_sphere ( soughtAngle,
918  M_PI / static_cast< proshade_double > ( dataObj->getMaxBand() ),
919  dataObj->getMaxBand() * 2,
920  soughtAngle,
921  static_cast<proshade_unsign> ( angIt - 1.0 ) ) );
922 
923  //== Interpolate rotation function onto the sphere
924  sphereMappedRotFun.at( static_cast < size_t > ( angIt - 1.0 ))->interpolateSphereValues ( dataObj->getInvSO3Coeffs ( ) );
925  }
926 
927  //== Save values to pointer
928  proshade_double* npVals = new proshade_double[(dataObj->getMaxBand() * 2) * (dataObj->getMaxBand() * 2) * (fold - 1)];
929  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
930 
931  //== Copy value to the pointer
932  std::vector< std::vector < proshade_double > > vls;
933  size_t valIter = 0;
934  for ( size_t sphIt = 0; sphIt < sphereMappedRotFun.size(); sphIt++ )
935  {
936  vls.clear();
937  vls = sphereMappedRotFun.at(sphIt)->getCopyOfValues ();
938 
939  for ( size_t lonIt = 0; lonIt < (dataObj->getMaxBand() * 2); lonIt++ )
940  {
941  for ( size_t latIt = 0; latIt < (dataObj->getMaxBand() * 2); latIt++ )
942  {
943  valIter = latIt + static_cast< size_t > ( dataObj->getMaxBand() * 2 ) * ( lonIt + static_cast< size_t > ( dataObj->getMaxBand() * 2 ) * sphIt );
944  npVals[valIter] = vls.at(lonIt).at(latIt);
945  }
946  }
947  }
948 
949  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
950  pybind11::capsule pyCapsuleRFMapSph ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
951 
952  //== Create python readable object with the correct memory access
953  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
954  { static_cast< int > ( sphereMappedRotFun.size() ), static_cast< int > (dataObj->getMaxBand() * 2), static_cast< int > (dataObj->getMaxBand() * 2) }, // Shape
955  { (dataObj->getMaxBand() * 2) * (dataObj->getMaxBand() * 2) * sizeof(proshade_double), (dataObj->getMaxBand() * 2) * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
956  npVals, // Data
957  pyCapsuleRFMapSph ); // Capsule (destructor)
958 
959  //== Done
960  return ( retArr );
961  }, "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)." );
962 }
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:1181
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:1847
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:3919
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:3092
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:4395
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:93
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_settings::allDetectedDAxes
std::vector< std::vector< proshade_unsign > > allDetectedDAxes
The vector of all detected dihedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:156
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:3899
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:3909
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:1746
ProSHADE_internal_data::ProSHADE_data::invertMirrorMap
void invertMirrorMap(ProSHADE_settings *settings)
Function for inverting the map to its mirror image.
Definition: ProSHADE_data.cpp:1175
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:1020
ProSHADE_settings::allDetectedTAxes
std::vector< proshade_unsign > allDetectedTAxes
The vector of all detected tetrahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:157
ProSHADE_settings::allDetectedOAxes
std::vector< proshade_unsign > allDetectedOAxes
The vector of all detected octahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:158
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:1520
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::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:1140
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:3929
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::readInStructure
void readInStructure(std::string fName, proshade_unsign inputO, ProSHADE_settings *settings, proshade_double *maskArr=nullptr, proshade_unsign maskXDim=0, proshade_unsign maskYDim=0, proshade_unsign maskZDim=0, proshade_double *weightsArr=nullptr, proshade_unsign weigXDim=0, proshade_unsign weigYDim=0, proshade_unsign weigZDim=0)
This function initialises the basic ProSHADE_data variables and reads in a single structure.
Definition: ProSHADE_data.cpp:509
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:1803
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:959
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:892
ProSHADE_internal_data::ProSHADE_data::getMaxBand
proshade_unsign getMaxBand(void)
This function returns the maximum band value for the object.
Definition: ProSHADE_data.cpp:3611
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:1075
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:868
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:996
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:3939
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:1276
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:1124
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:1695
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:3949
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:94
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:1966
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:1440
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:68
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_internal_data::ProSHADE_data::getInvSO3Coeffs
proshade_complex * getInvSO3Coeffs(void)
This function allows access to the inverse SO(3) coefficients array.
Definition: ProSHADE_data.cpp:3847
ProSHADE_settings::allDetectedIAxes
std::vector< proshade_unsign > allDetectedIAxes
The vector of all detected icosahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:159
ProSHADE_settings::forceBounds
proshade_signed * forceBounds
These will be the boundaries to be forced upon the map.
Definition: ProSHADE_settings.hpp:96
ProSHADE_settings::detectedSymmetry
std::vector< proshade_double * > detectedSymmetry
The vector of detected symmetry axes.
Definition: ProSHADE_settings.hpp:154
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:95
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:2939
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:4384
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:1229
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:155
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:3657
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:1585