23 #include <pybind11/pybind11.h>
24 #include <pybind11/stl.h>
25 #include <pybind11/numpy.h>
28 void add_dataClass ( pybind11::module& pyProSHADE )
31 pybind11::class_ < ProSHADE_internal_data::ProSHADE_data > ( pyProSHADE,
"ProSHADE_data" )
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 )
38 pybind11::buffer_info buf = mapData.request();
39 proshade_unsign len =
static_cast< proshade_unsign
> ( buf.size );
42 double* npVals =
new double[
static_cast<proshade_unsign
> ( len )];
48 for ( proshade_unsign iter = 0; iter < len; iter++ ) { npVals[iter] =
static_cast < double > ( mapData.at(iter) ); }
50 else if ( buf.ndim == 3 )
52 float* dataPtr =
reinterpret_cast < float*
> ( buf.ptr );
53 for ( proshade_unsign xIt = 0; xIt < static_cast<proshade_unsign> ( buf.shape.at(0) ); xIt++ )
55 for ( proshade_unsign yIt = 0; yIt < static_cast<proshade_unsign> ( buf.shape.at(1) ); yIt++ )
57 for ( proshade_unsign zIt = 0; zIt < static_cast<proshade_unsign> ( buf.shape.at(2) ); zIt++ )
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 )] );
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 );
73 static_cast<int> ( len ),
92 .def (
"readInStructure",
96 pybind11::buffer_info maskArr_buf = maskArr.request();
97 pybind11::buffer_info weightsArr_buf = weightsArr.request();
100 if ( ( maskArr_buf.shape.at(0) != 0 ) && ( maskArr_buf.shape.at(1) != 0 ) && ( maskArr_buf.shape.at(2) != 0 ) )
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 ); }
106 proshade_double* mskArr =
new proshade_double[maskArr_buf.shape.at(0)*maskArr_buf.shape.at(1)*maskArr_buf.shape.at(2)];
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]; }
114 if ( ( weightsArr_buf.shape.at(0) != 0 ) && ( weightsArr_buf.shape.at(1) != 0 ) && ( weightsArr_buf.shape.at(2) != 0 ) )
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 ); }
120 proshade_double* wghArr =
new proshade_double[weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2)];
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]; }
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) ) );
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 );
145 if ( ( weightsArr_buf.shape.at(0) != 0 ) && ( weightsArr_buf.shape.at(1) != 0 ) && ( weightsArr_buf.shape.at(2) != 0 ) )
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 ); }
151 proshade_double* wghArr =
new proshade_double[weightsArr_buf.shape.at(0)*weightsArr_buf.shape.at(1)*weightsArr_buf.shape.at(2)];
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]; }
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) ) );
167 self.readInStructure ( fName, inputO, settings,
nullptr, 0, 0, 0,
nullptr, 0, 0, 0 );
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 )
181 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( {
self.xDimIndices,
self.yDimIndices,
self.zDimIndices },
182 {
self.yDimIndices *
self.zDimIndices *
sizeof(proshade_double),
183 self.zDimIndices *
sizeof(proshade_double),
184 sizeof(proshade_double) },
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" ) )
200 .def (
"getReBoxBoundaries",
204 proshade_signed* retVals =
new proshade_signed[6];
210 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->
forceBounds[iter]; }
217 static_cast< proshade_signed
> (
self.xDimIndices ),
218 static_cast< proshade_signed
> (
self.yDimIndices ),
219 static_cast< proshade_signed
> (
self.zDimIndices ),
224 self.xDimSize,
self.yDimSize,
self.zDimSize, retVals, settings->
boundsExtraSpace );
230 if ( settings->
useSameBounds && (
self.inputOrder == 0 ) ) {
for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->
forceBounds[iter] = retVals[iter]; } }
234 pybind11::capsule pyCapsuleReBox2 ( retVals, [](
void *f ) { proshade_signed* foo =
reinterpret_cast< proshade_signed*
> ( f );
delete foo; } );
237 pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 },
238 {
sizeof(proshade_signed) },
244 },
"This function finds the boundaries enclosing positive map values and adds some extra space." )
245 .def (
"createNewMapFromBounds",
249 proshade_signed* newBounds =
new proshade_signed[6];
253 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
256 self.createNewMapFromBounds ( settings, newStr, newBounds );
260 },
"This function creates a new structure from the calling structure and new bounds values." )
277 .def (
"detectSymmetryInStructure",
282 },
"This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg (
"settings" ) )
285 .def (
"getRecommendedSymmetryAxes",
289 float* npVals =
new float[
static_cast<unsigned int> ( settings->
detectedSymmetry.size() * 7 )];
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] ); } }
296 pybind11::capsule pyCapsuleStrRecSym ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
299 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
detectedSymmetry.size() ),
static_cast<int> ( 7 ) },
300 { 7 *
sizeof(float),
sizeof(
float) },
302 pyCapsuleStrRecSym );
306 },
"This function returns the recommended symmetry axes as a 2D numpy array." )
307 .def (
"getAllCSyms",
311 float* npVals =
new float[
static_cast<unsigned int> ( settings->
allDetectedCAxes.size() * 7 )];
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) ); } }
318 pybind11::capsule pyCapsuleStrSymList ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
321 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
allDetectedCAxes.size() ), 7 },
322 { 7 *
sizeof(float),
sizeof(
float) },
324 pyCapsuleStrSymList );
328 },
"This function returns all symmetry axes as a 2D numpy array." )
329 .def (
"getNonCSymmetryAxesIndices",
333 pybind11::dict retDict;
334 pybind11::list dList, tList, oList, iList;
337 for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.size() ); dIt++ )
339 pybind11::list memList;
340 for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.at(dIt).size() ); memIt++ )
344 dList.append ( memList );
348 for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( settings->
allDetectedTAxes.size() ); tIt++ )
354 for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( settings->
allDetectedOAxes.size() ); oIt++ )
360 for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( settings->
allDetectedIAxes.size() ); iIt++ )
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;
373 },
"This function returns array of non-C axes indices." )
374 .def (
"getAllGroupElements",
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 ); }
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) ) ); }
386 std::vector < std::vector< proshade_double > > vals =
self.getAllGroupElements ( settings, axesList, groupType, matrixTolerance );
389 pybind11::list retList;
392 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
395 float* npVals =
new float[
static_cast<unsigned int> ( 9 )];
399 for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] =
static_cast< float > ( vals.at(iter).at(it) ); }
402 pybind11::capsule pyCapsuleGrpEl ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
405 pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 },
406 { 3 *
sizeof(float),
sizeof(
float) },
411 retList.append ( retArr );
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 )
418 .def (
"getMapCOMProcessChange",
422 std::vector< proshade_double > vals =
self.getMapCOMProcessChange ();
425 float* npVals =
new float[
static_cast<unsigned int> ( 3 )];
429 for ( proshade_unsign iter = 0; iter < 3; iter++ ) { npVals[iter] =
static_cast< float > ( vals.at(iter) ); }
432 pybind11::capsule pyCapsuleSymShiftDat ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
435 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( vals.size() ) },
438 pyCapsuleSymShiftDat );
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." )
446 .def (
"getBestRotationMapPeaksEulerAngles",
450 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
453 float* npVals =
new float[
static_cast<unsigned int> (vals.size())];
457 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] =
static_cast< float > ( vals.at(iter) ); }
460 pybind11::capsule pyCapsuleEuPeak ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
463 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<unsigned int> (vals.size()) },
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",
475 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
478 proshade_double* retMat =
new proshade_double[9];
483 pybind11::capsule pyCapsuleRMPeak ( retMat, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
486 pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 },
487 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
493 },
"This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg (
"settings" ) )
498 .def (
"getOverlayTranslations",
502 std::vector< proshade_double > vals =
self.getBestTranslationMapPeaksAngstrom ( staticStructure );
505 pybind11::dict retDict;
506 pybind11::list rotCen, toOverlay;
509 rotCen.append (
self.originalPdbRotCenX );
510 rotCen.append (
self.originalPdbRotCenY );
511 rotCen.append (
self.originalPdbRotCenZ );
513 toOverlay.append (
self.originalPdbTransX );
514 toOverlay.append (
self.originalPdbTransY );
515 toOverlay.append (
self.originalPdbTransZ );
518 retDict[ pybind11::handle ( pybind11::str (
"centreOfRotation" ).ptr ( ) ) ] = rotCen;
519 retDict[ pybind11::handle ( pybind11::str (
"rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
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" ) )
527 .def (
"findSHIndex",
531 proshade_signed index = seanindex (
static_cast< int > ( order ),
532 static_cast< int > ( band ),
533 static_cast< int > (
self.spheres[shell]->getLocalBandwidth() ) );
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",
542 std::complex<proshade_double>* npVals =
new std::complex<proshade_double>[
static_cast<proshade_unsign
> (
self.noSpheres * pow(
self.maxShellBand, 2.0 ) )];
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 ); }
549 proshade_signed pyPosSH;
550 proshade_signed pyPos;
553 for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> (
self.noSpheres ); shIt++ )
555 for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> (
self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
557 for ( proshade_signed order = -bnd; order <= bnd; order++ )
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] );
570 pybind11::capsule pyCapsuleSHs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
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 > ) },
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." )
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 ) )];
590 proshade_signed index = 0;
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 ); }
596 proshade_double emReal, emImag;
597 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
599 for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
601 for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
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 );
612 pybind11::capsule pyCapsuleEMs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
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 > ) },
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",
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 ) )];
635 proshade_signed index = 0;
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 ); }
641 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
643 for ( proshade_signed order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
645 for ( proshade_signed order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
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] );
655 pybind11::capsule pyCapsuleSOCoeffs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
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 > ) },
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",
674 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> ( (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) )];
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] ); }
681 pybind11::capsule pyCapsuleRotMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
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 ) },
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 >) },
694 },
"This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
695 .def (
"getRotationMatrixFromSOFTCoordinates",
699 proshade_double* npVals =
new proshade_double[9];
703 proshade_double eulA, eulB, eulG;
712 pybind11::capsule pyCapsuleRMSoft ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
715 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
716 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
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",
727 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.getXDim() *
self.getYDim() *
self.getZDim() )];
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] ); }
734 pybind11::capsule pyCapsuleTrsMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
737 pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
738 {
self.getXDim(),
self.getYDim(),
self.getZDim() },
739 {
self.getYDim() *
self.getZDim() *
sizeof(std::complex < proshade_double >),
740 self.getZDim() *
sizeof(std::complex < proshade_double >),
741 sizeof(std::complex < proshade_double >) },
747 },
"This function returns the translation function as a three-dimensional map of complex numbers." )
785 .def (
"__repr__", [] ( ) {
return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );
788 pyProSHADE.def (
"computeGroupElementsForGroup",
789 [] ( proshade_double x, proshade_double y, proshade_double z, proshade_unsign fold ) -> pybind11::array_t < proshade_double >
795 proshade_double* npVals =
new proshade_double[
static_cast<proshade_unsign
> ( retVec.size() ) * 9];
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); } }
802 pybind11::capsule pyCapsuleGrEls ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
805 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
806 {
static_cast< int > ( retVec.size() ), 3, 3 },
807 { 9 *
sizeof(proshade_double), 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
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 >
821 pybind11::buffer_info buf1 = first.request();
822 pybind11::buffer_info buf2 = second.request();
825 if ( buf1.ndim != 3 || buf2.ndim != 3 )
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 > () );
832 std::vector< std::vector< proshade_double > > fVec, sVec;
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++ )
837 std::vector< proshade_double > rotMat;
838 for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf1.shape.at(1) ); rowIt++ )
840 for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf1.shape.at(2) ); colIt++ )
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++ )
851 std::vector< proshade_double > rotMat;
852 for ( proshade_unsign rowIt = 0; rowIt < static_cast<proshade_unsign> ( buf2.shape.at(1) ); rowIt++ )
854 for ( proshade_unsign colIt = 0; colIt < static_cast<proshade_unsign> ( buf2.shape.at(2) ); colIt++ )
867 proshade_double* npVals =
new proshade_double[
static_cast<proshade_unsign
> ( retVec.size() ) * 9];
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); } }
874 pybind11::capsule pyCapsuleGrElsCombo ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
877 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > (
878 {
static_cast< int > ( retVec.size() ), 3, 3 },
879 { 9 *
sizeof(proshade_double), 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
881 pyCapsuleGrElsCombo );
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" ) );