ProSHADE  0.7.6.1 (AUG 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 ( "getReBoxBoundaries",
201  [] ( ProSHADE_internal_data::ProSHADE_data &self,ProSHADE_settings* settings ) -> pybind11::array_t < proshade_signed >
202  {
203  //== Allocate output memory
204  proshade_signed* retVals = new proshade_signed[6];
205  ProSHADE_internal_misc::checkMemoryAllocation ( retVals, __FILE__, __LINE__, __func__ );
206 
207  //== If same bounds as first one are required, test if possible and return these instead
208  if ( settings->useSameBounds && ( self.inputOrder != 0 ) )
209  {
210  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->forceBounds[iter]; }
211  }
212  //== In this case, bounds need to be found de novo
213  else
214  {
215  //== Find the non-zero bounds
217  static_cast< proshade_signed > ( self.xDimIndices ),
218  static_cast< proshade_signed > ( self.yDimIndices ),
219  static_cast< proshade_signed > ( self.zDimIndices ),
220  retVals );
221 
222  //== Add the extra space
223  ProSHADE_internal_mapManip::addExtraBoundSpace ( self.xDimIndices, self.yDimIndices, self.zDimIndices,
224  self.xDimSize, self.yDimSize, self.zDimSize, retVals, settings->boundsExtraSpace );
225 
226  //== Beautify boundaries
227  ProSHADE_internal_mapManip::beautifyBoundaries ( retVals, self.xDimIndices, self.yDimIndices, self.zDimIndices, settings->boundsSimilarityThreshold );
228 
229  //== If need be, save boundaries to be used for all other structure
230  if ( settings->useSameBounds && ( self.inputOrder == 0 ) ) { for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->forceBounds[iter] = retVals[iter]; } }
231  }
232 
233  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
234  pybind11::capsule pyCapsuleReBox2 ( retVals, []( void *f ) { proshade_signed* foo = reinterpret_cast< proshade_signed* > ( f ); delete foo; } );
235 
236  //== Copy the value
237  pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 }, // Shape
238  { sizeof(proshade_signed) }, // C-stype strides
239  retVals, // Data
240  pyCapsuleReBox2 ); // Capsule
241 
242  //== Done
243  return ( retArr );
244  }, "This function finds the boundaries enclosing positive map values and adds some extra space." )
245  .def ( "createNewMapFromBounds",
246  [] ( ProSHADE_internal_data::ProSHADE_data &self, pybind11::array_t < proshade_signed > bounds, ProSHADE_internal_data::ProSHADE_data* newStr, ProSHADE_settings* settings )
247  {
248  //== Allocate memory for bounds copy
249  proshade_signed* newBounds = new proshade_signed[6];
250  ProSHADE_internal_misc::checkMemoryAllocation ( newBounds, __FILE__, __LINE__, __func__ );
251 
252  //== Copy to the pointer
253  for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
254 
255  //== Run C++ function
256  self.createNewMapFromBounds ( settings, newStr, newBounds );
257 
258  //== Done
259  return ;
260  }, "This function creates a new structure from the calling structure and new bounds values." )
261 
262  //============================================ Data sphere mapping functions
263  .def ( "mapToSpheres", &ProSHADE_internal_data::ProSHADE_data::mapToSpheres, "This function converts the internal map onto a set of concentric spheres.", pybind11::arg ( "settings" ) )
264  .def ( "getSpherePositions", &ProSHADE_internal_data::ProSHADE_data::getSpherePositions, "This function determines the sphere positions (radii) for sphere mapping.", pybind11::arg ( "settings" ) )
265  .def ( "computeSphericalHarmonics", &ProSHADE_internal_data::ProSHADE_data::computeSphericalHarmonics, "This function computes the spherical harmonics decomposition for the whole structure.", pybind11::arg ( "settings" ) )
266 
267  //============================================ Accessor functions
268  .def ( "getXDimSize", &ProSHADE_internal_data::ProSHADE_data::getXDimSize, "This function allows access to the map size in angstroms along the X axis." )
269  .def ( "getYDimSize", &ProSHADE_internal_data::ProSHADE_data::getYDimSize, "This function allows access to the map size in angstroms along the Y axis." )
270  .def ( "getZDimSize", &ProSHADE_internal_data::ProSHADE_data::getZDimSize, "This function allows access to the map size in angstroms along the Z axis." )
271  .def ( "getXDim", &ProSHADE_internal_data::ProSHADE_data::getXDim, "This function allows access to the map size in indices along the X axis." )
272  .def ( "getYDim", &ProSHADE_internal_data::ProSHADE_data::getYDim, "This function allows access to the map size in indices along the Y axis." )
273  .def ( "getZDim", &ProSHADE_internal_data::ProSHADE_data::getZDim, "This function allows access to the map size in indices along the Z axis." )
274 
275  //============================================ Symmetry related functions
276  .def ( "computeRotationFunction", &ProSHADE_internal_data::ProSHADE_data::computeRotationFunction, "This function computes the self-rotation function for this structure and stores it internally in the ProSHADE_data object.", pybind11::arg ( "settings" ) )
277  .def ( "detectSymmetryInStructure",
279  {
280  //== Call the appropriate C++ function
281  self.detectSymmetryFromAngleAxisSpace ( settings, &settings->detectedSymmetry, &settings->allDetectedCAxes );
282  }, "This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg ( "settings" ) )
283  .def ( "getRecommendedSymmetryType", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryType, "This function simply returns the detected recommended symmetry type.", pybind11::arg ( "settings" ) )
284  .def ( "getRecommendedSymmetryFold", &ProSHADE_internal_data::ProSHADE_data::getRecommendedSymmetryFold, "This function simply returns the detected recommended symmetry fold.", pybind11::arg ( "settings" ) )
285  .def ( "getRecommendedSymmetryAxes",
286  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
287  {
288  //== Allocate memory for the numpy values
289  float* npVals = new float[static_cast<unsigned int> ( settings->detectedSymmetry.size() * 7 )];
290  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
291 
292  //== Copy values
293  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] ); } }
294 
295  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
296  pybind11::capsule pyCapsuleStrRecSym ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
297 
298  //== Copy the value
299  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( settings->detectedSymmetry.size() ), static_cast<int> ( 7 ) }, // Shape
300  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
301  npVals, // Data
302  pyCapsuleStrRecSym ); // Capsule
303 
304  //== Done
305  return ( retArr );
306  }, "This function returns the recommended symmetry axes as a 2D numpy array." )
307  .def ( "getAllCSyms",
308  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
309  {
310  //== Allocate memory for the numpy values
311  float* npVals = new float[static_cast<unsigned int> ( settings->allDetectedCAxes.size() * 7 )];
312  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
313 
314  //== Copy values
315  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) ); } }
316 
317  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
318  pybind11::capsule pyCapsuleStrSymList ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
319 
320  //== Copy the value
321  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( settings->allDetectedCAxes.size() ), 7 }, // Shape
322  { 7 * sizeof(float), sizeof(float) }, // C-stype strides
323  npVals, // Data
324  pyCapsuleStrSymList ); // Capsule
325 
326  //== Done
327  return ( retArr );
328  }, "This function returns all symmetry axes as a 2D numpy array." )
329  .def ( "getNonCSymmetryAxesIndices",
330  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::dict
331  {
332  //== Initialise local variables
333  pybind11::dict retDict;
334  pybind11::list dList, tList, oList, iList;
335 
336  //== Fill in the D list
337  for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( settings->allDetectedDAxes.size() ); dIt++ )
338  {
339  pybind11::list memList;
340  for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( settings->allDetectedDAxes.at(dIt).size() ); memIt++ )
341  {
342  memList.append ( settings->allDetectedDAxes.at(dIt).at(memIt) );
343  }
344  dList.append ( memList );
345  }
346 
347  //== Fill in the T list
348  for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( settings->allDetectedTAxes.size() ); tIt++ )
349  {
350  tList.append ( settings->allDetectedTAxes.at ( tIt ) );
351  }
352 
353  //== Fill in the O list
354  for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( settings->allDetectedOAxes.size() ); oIt++ )
355  {
356  oList.append ( settings->allDetectedOAxes.at ( oIt ) );
357  }
358 
359  //== Fill in the T list
360  for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( settings->allDetectedIAxes.size() ); iIt++ )
361  {
362  iList.append ( settings->allDetectedIAxes.at ( iIt ) );
363  }
364 
365  //== Save the lists to the dict
366  retDict[ pybind11::handle ( pybind11::str ( "D" ).ptr ( ) ) ] = dList;
367  retDict[ pybind11::handle ( pybind11::str ( "T" ).ptr ( ) ) ] = tList;
368  retDict[ pybind11::handle ( pybind11::str ( "O" ).ptr ( ) ) ] = oList;
369  retDict[ pybind11::handle ( pybind11::str ( "I" ).ptr ( ) ) ] = iList;
370 
371  //== Done
372  return ( retDict );
373  }, "This function returns array of non-C axes indices." )
374  .def ( "getAllGroupElements",
375  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings, pybind11::array_t < int > axList, std::string groupType, proshade_double matrixTolerance ) -> pybind11::list
376  {
377  //== Sanity check
378  pybind11::buffer_info buf = axList.request();
379  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 ); }
380 
381  //== Convert to vector of unsigns
382  std::vector< proshade_unsign > axesList;
383  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) ) ); }
384 
385  //== Get the results
386  std::vector < std::vector< proshade_double > > vals = self.getAllGroupElements ( settings, axesList, groupType, matrixTolerance );
387 
388  //== Initialise return variable
389  pybind11::list retList;
390 
391  //== Copy values
392  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
393  {
394  //== Allocate memory for the numpy values
395  float* npVals = new float[static_cast<unsigned int> ( 9 )];
396  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
397 
398  //== Copy values to memory
399  for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] = static_cast< float > ( vals.at(iter).at(it) ); }
400 
401  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
402  pybind11::capsule pyCapsuleGrpEl ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
403 
404  //== Copy the value
405  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 }, // Shape
406  { 3 * sizeof(float), sizeof(float) }, // C-stype strides
407  npVals, // Data
408  pyCapsuleGrpEl ); // Capsule
409 
410  //== Save matrix
411  retList.append ( retArr );
412  }
413 
414  //== Done
415  return ( retList );
416  }, "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 )
417 
418  .def ( "getMapCOMProcessChange",
419  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < float >
420  {
421  //== Get the values
422  std::vector< proshade_double > vals = self.getMapCOMProcessChange ();
423 
424  //== Allocate memory for the numpy values
425  float* npVals = new float[static_cast<unsigned int> ( 3 )];
426  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
427 
428  //== Copy values
429  for ( proshade_unsign iter = 0; iter < 3; iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
430 
431  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
432  pybind11::capsule pyCapsuleSymShiftDat ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
433 
434  //== Copy the value
435  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<int> ( vals.size() ) }, // Shape
436  { sizeof(float) }, // C-stype strides
437  npVals, // Data
438  pyCapsuleSymShiftDat ); // Capsule
439 
440  //== Done
441  return ( retArr );
442  }, "This function returns the shift in Angstrom applied to the internal map representation in order to align its COM with the centre of box." )
443 
444  //============================================ Overlay related functions
445  .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" ) )
446  .def ( "getBestRotationMapPeaksEulerAngles",
447  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < float >
448  {
449  //== Get values
450  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
451 
452  //== Allocate memory for the numpy values
453  float* npVals = new float[static_cast<unsigned int> (vals.size())];
454  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
455 
456  //== Copy values
457  for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] = static_cast< float > ( vals.at(iter) ); }
458 
459  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
460  pybind11::capsule pyCapsuleEuPeak ( npVals, []( void *f ) { float* foo = reinterpret_cast< float* > ( f ); delete foo; } );
461 
462  //== Copy the value
463  pybind11::array_t < float > retArr = pybind11::array_t<float> ( { static_cast<unsigned int> (vals.size()) }, // Shape
464  { sizeof(float) }, // C-stype strides
465  npVals, // Data
466  pyCapsuleEuPeak ); // Capsule (C++ destructor, basically)
467 
468  //== Done
469  return ( retArr );
470  }, "This function returns a vector of three floats, the three Euler angles of the best peak in the rotation map.", pybind11::arg ( "settings" ) )
471  .def ( "getBestRotationMapPeaksRotationMatrix",
472  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_settings* settings ) -> pybind11::array_t < proshade_double >
473  {
474  //== Get values
475  std::vector< proshade_double > vals = self.getBestRotationMapPeaksEulerAngles ( settings );
476 
477  //== Convert Euler ZXZ to matrix
478  proshade_double* retMat = new proshade_double[9];
479  ProSHADE_internal_misc::checkMemoryAllocation ( retMat, __FILE__, __LINE__, __func__ );
480  ProSHADE_internal_maths::getRotationMatrixFromEulerZXZAngles ( vals.at(0), vals.at(1), vals.at(2), retMat );
481 
482  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
483  pybind11::capsule pyCapsuleRMPeak ( retMat, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
484 
485  //== Copy the value
486  pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 }, // Shape
487  { 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
488  retMat, // Data
489  pyCapsuleRMPeak ); // Capsule
490 
491  //== Done
492  return ( retArr );
493  }, "This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg ( "settings" ) )
494  .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" ) )
495  .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" ) )
496  .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" ) )
497  .def ( "computeTranslationMap", &ProSHADE_internal_data::ProSHADE_data::computeTranslationMap, "This function does the computation of the translation map and saves results internally.", pybind11::arg ( "staticStructure" ) )
498  .def ( "getOverlayTranslations",
499  [] ( ProSHADE_internal_data::ProSHADE_data &self, ProSHADE_internal_data::ProSHADE_data* staticStructure ) -> pybind11::dict
500  {
501  //== Get values
502  std::vector< proshade_double > vals = self.getBestTranslationMapPeaksAngstrom ( staticStructure );
503 
504  //== Initialise variables
505  pybind11::dict retDict;
506  pybind11::list rotCen, toOverlay;
507 
508  //== Copy values
509  rotCen.append ( self.originalPdbRotCenX );
510  rotCen.append ( self.originalPdbRotCenY );
511  rotCen.append ( self.originalPdbRotCenZ );
512 
513  toOverlay.append ( self.originalPdbTransX );
514  toOverlay.append ( self.originalPdbTransY );
515  toOverlay.append ( self.originalPdbTransZ );
516 
517  //== Save results to return dict
518  retDict[ pybind11::handle ( pybind11::str ( "centreOfRotation" ).ptr ( ) ) ] = rotCen;
519  retDict[ pybind11::handle ( pybind11::str ( "rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
520 
521  //== Done
522  return ( retDict );
523  }, "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" ) )
524  .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" ) )
525 
526  //============================================ Internal arrays access functions
527  .def ( "findSHIndex",
528  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed shell, proshade_signed band, proshade_signed order ) -> proshade_signed
529  {
530  //== Get value
531  proshade_signed index = seanindex ( static_cast< int > ( order ),
532  static_cast< int > ( band ),
533  static_cast< int > ( self.spheres[shell]->getLocalBandwidth() ) );
534 
535  //== Done
536  return ( index );
537  }, "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." )
538  .def ( "getSphericalHarmonics",
539  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex<proshade_double> >
540  {
541  //== Allocate memory for the numpy values
542  std::complex<proshade_double>* npVals = new std::complex<proshade_double>[static_cast<proshade_unsign> ( self.noSpheres * pow( self.maxShellBand, 2.0 ) )];
543  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
544 
545  //== Fill with zeroes as the matrix will be sparse
546  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 ); }
547 
548  //== Initialise variables
549  proshade_signed pyPosSH;
550  proshade_signed pyPos;
551 
552  //== Copy data to new memory
553  for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> ( self.noSpheres ); shIt++ )
554  {
555  for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> ( self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
556  {
557  for ( proshade_signed order = -bnd; order <= bnd; order++ )
558  {
559  pyPosSH = ( static_cast< proshade_signed > ( shIt ) * static_cast< proshade_signed > ( std::pow ( self.maxShellBand, 2 ) ) );
560  pyPos = seanindex ( static_cast< int > ( order ),
561  static_cast< int > ( bnd ),
562  static_cast< int > ( self.spheres[shIt]->getLocalBandwidth() ) );
563  npVals[pyPosSH+pyPos].real ( self.sphericalHarmonics[shIt][pyPos][0] );
564  npVals[pyPosSH+pyPos].imag ( self.sphericalHarmonics[shIt][pyPos][1] );
565  }
566  }
567  }
568 
569  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
570  pybind11::capsule pyCapsuleSHs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
571 
572  //== Copy the value
573  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
574  { static_cast<int> ( self.noSpheres ), static_cast<int> ( pow ( self.maxShellBand, 2.0 ) ) },
575  { sizeof ( std::complex < proshade_double > ) * static_cast< proshade_unsign > ( std::pow ( static_cast< proshade_unsign > ( self.maxShellBand ), 2 ) ), sizeof ( std::complex < proshade_double > ) },
576  npVals,
577  pyCapsuleSHs );
578 
579  //== Done
580  return ( retArr );
581  }, "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." )
582  .def ( "getEMatrix",
583  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
584  {
585  //== Allocate memory for the numpy values
586  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) )];
587  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
588 
589  //== Allocate local variables
590  proshade_signed index = 0;
591 
592  //== Fill with zeroes as the matrix will be sparse
593  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 ); }
594 
595  //== Copy data to new memory
596  proshade_double emReal, emImag;
597  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxShellBand ); bandIter++ )
598  {
599  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
600  {
601  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
602  {
603  index = order2 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * bandIter );
604  self.getEMatrixValue( static_cast< proshade_unsign > ( bandIter ), static_cast< proshade_unsign > ( order1 ), static_cast< proshade_unsign > ( order2 ), &emReal, &emImag );
605  npVals[index].real ( emReal );
606  npVals[index].imag ( emImag );
607  }
608  }
609  }
610 
611  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
612  pybind11::capsule pyCapsuleEMs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
613 
614  //== Create the output object
615  proshade_unsign orderSize = ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 );
616  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
617  { self.maxShellBand, orderSize, orderSize },
618  { sizeof ( std::complex < proshade_double > ) * orderSize * orderSize,
619  sizeof ( std::complex < proshade_double > ) * orderSize,
620  sizeof ( std::complex < proshade_double > ) },
621  npVals,
622  pyCapsuleEMs );
623 
624  //== Done
625  return ( retArr );
626  }, "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." )
627  .def ( "getSO3Coefficients",
628  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
629  {
630  //== Allocate memory for the numpy values
631  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.maxShellBand * ( ( self.maxShellBand * 2 ) + 1 ) * ( ( self.maxShellBand * 2 ) + 1 ) )];
632  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
633 
634  //== Allocate local variables
635  proshade_signed index = 0;
636 
637  //== Fill with zeroes as the matrix will be sparse
638  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 ); }
639 
640  //== Copy data to new memory
641  for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > ( self.maxShellBand ); bandIter++ )
642  {
643  for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
644  {
645  for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
646  {
647  index = order2 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * ( order1 + ( ( static_cast< proshade_signed > ( self.maxShellBand ) * 2 ) + 1 ) * bandIter );
648  npVals[index].real ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][0] );
649  npVals[index].imag ( self.so3Coeffs[self.so3CoeffsArrayIndex ( static_cast< proshade_signed > ( order1 - bandIter ), order2 - bandIter, bandIter )][1] );
650  }
651  }
652  }
653 
654  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
655  pybind11::capsule pyCapsuleSOCoeffs ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
656 
657  //== Create the output object
658  proshade_unsign orderSize = ( ( static_cast< proshade_unsign > ( self.maxShellBand ) * 2 ) + 1 );
659  pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
660  { self.maxShellBand, orderSize, orderSize },
661  { sizeof ( std::complex < proshade_double > ) * orderSize * orderSize,
662  sizeof ( std::complex < proshade_double > ) * orderSize,
663  sizeof ( std::complex < proshade_double > ) },
664  npVals,
665  pyCapsuleSOCoeffs );
666 
667  //== Done
668  return ( retArr );
669  }, "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." )
670  .def ( "getRotationFunctionMap",
671  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
672  {
673  //== Allocate memory for the numpy values
674  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) )];
675  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
676 
677  //== Copy the values
678  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] ); }
679 
680  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
681  pybind11::capsule pyCapsuleRotMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
682 
683  //== Create python readable object with the correct memory access
684  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
685  { ( self.maxShellBand * 2 ), ( self.maxShellBand * 2 ), ( self.maxShellBand * 2 ) }, // Shape
686  { ( self.maxShellBand * 2 ) * ( self.maxShellBand * 2 ) * sizeof(std::complex < proshade_double >),
687  ( self.maxShellBand * 2 ) * sizeof(std::complex < proshade_double >),
688  sizeof(std::complex < proshade_double >) }, // C-stype strides
689  npVals, // Data
690  pyCapsuleRotMap ); // Capsule (destructor)
691 
692  //== Done
693  return ( retArr );
694  }, "This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
695  .def ( "getRotationMatrixFromSOFTCoordinates",
696  [] ( ProSHADE_internal_data::ProSHADE_data &self, proshade_signed xPos, proshade_signed yPos, proshade_signed zPos ) -> pybind11::array_t < proshade_double >
697  {
698  //== Allocate memory for the numpy values
699  proshade_double* npVals = new proshade_double[9];
700  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
701 
702  //== Initialise local variables
703  proshade_double eulA, eulB, eulG;
704 
705  //== Compute the Euler angles from SOFT position
706  ProSHADE_internal_maths::getEulerZXZFromSOFTPosition ( static_cast< proshade_signed > ( self.maxShellBand ), xPos, yPos, zPos, &eulA, &eulB, &eulG );
707 
708  //== Compute the rotation matrix
710 
711  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
712  pybind11::capsule pyCapsuleRMSoft ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
713 
714  //== Create python readable object with the correct memory access
715  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
716  { 3 * sizeof(proshade_double), sizeof(proshade_double) },
717  npVals,
718  pyCapsuleRMSoft );
719 
720  //== Done
721  return ( retArr );
722  }, "This function converts a given rotation function map position onto the corresponding rotation matrix.", pybind11::arg ( "xPos" ), pybind11::arg ( "yPos" ), pybind11::arg ( "zPos" ) )
723  .def ( "getTranslationFunctionMap",
724  [] ( ProSHADE_internal_data::ProSHADE_data &self ) -> pybind11::array_t < std::complex < proshade_double > >
725  {
726  //== Allocate memory for the numpy values
727  std::complex<proshade_double>* npVals = new std::complex < proshade_double >[static_cast<proshade_unsign> ( self.getXDim() * self.getYDim() * self.getZDim() )];
728  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
729 
730  //== Copy the values
731  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] ); }
732 
733  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
734  pybind11::capsule pyCapsuleTrsMap ( npVals, []( void *f ) { std::complex<proshade_double>* foo = reinterpret_cast< std::complex<proshade_double>* > ( f ); delete foo; } );
735 
736  //== Create python readable object with the correct memory access
737  pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
738  { self.getXDim(), self.getYDim(), self.getZDim() }, // Shape
739  { self.getYDim() * self.getZDim() * sizeof(std::complex < proshade_double >),
740  self.getZDim() * sizeof(std::complex < proshade_double >),
741  sizeof(std::complex < proshade_double >) }, // C-stype strides
742  npVals, // Data
743  pyCapsuleTrsMap ); // Capsule (destructor)
744 
745  //== Done
746  return ( retArr );
747  }, "This function returns the translation function as a three-dimensional map of complex numbers." )
748 
749  //============================================ Member variables
750  .def_readwrite ( "fileName", &ProSHADE_internal_data::ProSHADE_data::fileName )
751  .def_readwrite ( "xDimSize", &ProSHADE_internal_data::ProSHADE_data::xDimSize )
752  .def_readwrite ( "yDimSize", &ProSHADE_internal_data::ProSHADE_data::yDimSize )
753  .def_readwrite ( "zDimSize", &ProSHADE_internal_data::ProSHADE_data::zDimSize )
754  .def_readwrite ( "aAngle", &ProSHADE_internal_data::ProSHADE_data::aAngle )
755  .def_readwrite ( "bAngle", &ProSHADE_internal_data::ProSHADE_data::bAngle )
756  .def_readwrite ( "cAngle", &ProSHADE_internal_data::ProSHADE_data::cAngle )
757  .def_readwrite ( "xDimIndices", &ProSHADE_internal_data::ProSHADE_data::xDimIndices )
758  .def_readwrite ( "yDimIndices", &ProSHADE_internal_data::ProSHADE_data::yDimIndices )
759  .def_readwrite ( "zDimIndices", &ProSHADE_internal_data::ProSHADE_data::zDimIndices )
760  .def_readwrite ( "xGridIndices", &ProSHADE_internal_data::ProSHADE_data::xGridIndices )
761  .def_readwrite ( "yGridIndices", &ProSHADE_internal_data::ProSHADE_data::yGridIndices )
762  .def_readwrite ( "zGridIndices", &ProSHADE_internal_data::ProSHADE_data::zGridIndices )
763  .def_readwrite ( "xAxisOrder", &ProSHADE_internal_data::ProSHADE_data::xAxisOrder )
764  .def_readwrite ( "yAxisOrder", &ProSHADE_internal_data::ProSHADE_data::yAxisOrder )
765  .def_readwrite ( "zAxisOrder", &ProSHADE_internal_data::ProSHADE_data::zAxisOrder )
766  .def_readwrite ( "xAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::xAxisOrigin )
767  .def_readwrite ( "yAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::yAxisOrigin )
768  .def_readwrite ( "zAxisOrigin", &ProSHADE_internal_data::ProSHADE_data::zAxisOrigin )
769  .def_readwrite ( "xCom", &ProSHADE_internal_data::ProSHADE_data::xCom )
770  .def_readwrite ( "yCom", &ProSHADE_internal_data::ProSHADE_data::yCom )
771  .def_readwrite ( "zCom", &ProSHADE_internal_data::ProSHADE_data::zCom )
772  .def_readwrite ( "xFrom", &ProSHADE_internal_data::ProSHADE_data::xFrom )
773  .def_readwrite ( "yFrom", &ProSHADE_internal_data::ProSHADE_data::yFrom )
774  .def_readwrite ( "zFrom", &ProSHADE_internal_data::ProSHADE_data::zFrom )
775  .def_readwrite ( "xTo", &ProSHADE_internal_data::ProSHADE_data::xTo )
776  .def_readwrite ( "yTo", &ProSHADE_internal_data::ProSHADE_data::yTo )
777  .def_readwrite ( "zTo", &ProSHADE_internal_data::ProSHADE_data::zTo )
778  .def_readwrite ( "spherePos", &ProSHADE_internal_data::ProSHADE_data::spherePos )
779  .def_readwrite ( "noSpheres", &ProSHADE_internal_data::ProSHADE_data::noSpheres )
780  .def_readwrite ( "maxShellBand", &ProSHADE_internal_data::ProSHADE_data::maxShellBand )
781  .def_readwrite ( "isEmpty", &ProSHADE_internal_data::ProSHADE_data::isEmpty )
782  .def_readwrite ( "inputOrder", &ProSHADE_internal_data::ProSHADE_data::inputOrder )
783 
784  //============================================ Description
785  .def ( "__repr__", [] ( ) { return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );
786 
787  //================================================ Export extra symmetry elements functions
788  pyProSHADE.def ( "computeGroupElementsForGroup",
789  [] ( proshade_double x, proshade_double y, proshade_double z, proshade_unsign fold ) -> pybind11::array_t < proshade_double >
790  {
791  //== Get the results
792  std::vector<std::vector< proshade_double > > retVec = ProSHADE_internal_data::computeGroupElementsForGroup ( x, y, z, static_cast< proshade_signed > ( fold ) );
793 
794  //== Allocate memory for the numpy values
795  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
796  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
797 
798  //== Copy the values
799  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); } }
800 
801  //== Create capsules to make sure memory is released properly from the allocating language (C++ in this case)
802  pybind11::capsule pyCapsuleGrEls ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
803 
804  //== Create python readable object with the correct memory access
805  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
806  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
807  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
808  npVals, // Data
809  pyCapsuleGrEls ); // Capsule (destructor)
810 
811  //== Done
812  return ( retArr );
813  }, "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" ) );
814  pyProSHADE.def ( "joinElementsFromDifferentGroups",
815  [] ( pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > first,
816  pybind11::array_t < proshade_double, pybind11::array::c_style | pybind11::array::forcecast > second,
817  proshade_double matrixTolerance,
818  bool combine ) -> pybind11::array_t < proshade_double >
819  {
820  //== Get array buffers
821  pybind11::buffer_info buf1 = first.request();
822  pybind11::buffer_info buf2 = second.request();
823 
824  //== Sanity check
825  if ( buf1.ndim != 3 || buf2.ndim != 3 )
826  {
827  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;
828  return ( pybind11::array_t < proshade_double > () );
829  }
830 
831  //== Convert input to C++ vectors
832  std::vector< std::vector< proshade_double > > fVec, sVec;
833 
834  proshade_double* dataPtr1 = reinterpret_cast < proshade_double* > ( buf1.ptr );
835  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf1.shape.at(0) ); elIt++ )
836  {
837  std::vector< proshade_double > rotMat;
838  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf1.shape.at(1) ); rowIt++ )
839  {
840  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf1.shape.at(2) ); colIt++ )
841  {
842  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 )] );
843  }
844  }
846  }
847 
848  proshade_double* dataPtr2 = reinterpret_cast < proshade_double* > ( buf2.ptr );
849  for ( proshade_unsign elIt = 0; elIt < static_cast<proshade_unsign> ( buf2.shape.at(0) ); elIt++ )
850  {
851  std::vector< proshade_double > rotMat;
852  for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf2.shape.at(1) ); rowIt++ )
853  {
854  for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf2.shape.at(2) ); colIt++ )
855  {
856  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 )] );
857  }
858  }
860  }
861 
862 
863  //== Get the results
864  std::vector< std::vector< proshade_double > > retVec = ProSHADE_internal_data::joinElementsFromDifferentGroups ( &fVec, &sVec, matrixTolerance, combine );
865 
866  //== Allocate memory for the numpy values
867  proshade_double* npVals = new proshade_double[static_cast<proshade_unsign> ( retVec.size() ) * 9];
868  ProSHADE_internal_misc::checkMemoryAllocation ( npVals, __FILE__, __LINE__, __func__ );
869 
870  //== Copy the values
871  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); } }
872 
873  //== Create capsule to make sure memory is released properly from the allocating language (C++ in this case)
874  pybind11::capsule pyCapsuleGrElsCombo ( npVals, []( void *f ) { proshade_double* foo = reinterpret_cast< proshade_double* > ( f ); delete foo; } );
875 
876  //== Create python readable object with the correct memory access
877  pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
878  { static_cast< int > ( retVec.size() ), 3, 3 }, // Shape
879  { 9 * sizeof(proshade_double), 3 * sizeof(proshade_double), sizeof(proshade_double) }, // C-stype strides
880  npVals, // Data
881  pyCapsuleGrElsCombo ); // Capsule (destructor)
882 
883  //== Done
884  return ( retArr );
885  }, "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" ) );
886 }
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:1843
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:3898
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:3102
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:4374
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_settings::allDetectedDAxes
std::vector< std::vector< proshade_unsign > > allDetectedDAxes
The vector of all detected dihedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:148
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:3878
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:3888
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:1742
ProSHADE_internal_data::ProSHADE_data::invertMirrorMap
void invertMirrorMap(ProSHADE_settings *settings)
Function for inverting the map to its mirror image.
Definition: ProSHADE_data.cpp:1188
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:1007
ProSHADE_settings::allDetectedTAxes
std::vector< proshade_unsign > allDetectedTAxes
The vector of all detected tetrahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:149
ProSHADE_settings::allDetectedOAxes
std::vector< proshade_unsign > allDetectedOAxes
The vector of all detected octahedral symmetry axes indices in allDetectedCAxes.
Definition: ProSHADE_settings.hpp:150
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:1525
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:3908
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:657
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:1799
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:1072
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:1005
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:954
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:1109
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:3918
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:1289
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:1699
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:891
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:3928
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: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:1453
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:151
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:146
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:2949
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:4363
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:1242
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:147
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:3645
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:519
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:1590