22 #include <pybind11/pybind11.h>
23 #include <pybind11/stl.h>
24 #include <pybind11/numpy.h>
27 void add_dataClass ( pybind11::module& pyProSHADE )
30 pybind11::class_ < ProSHADE_internal_data::ProSHADE_data > ( pyProSHADE,
"ProSHADE_data" )
33 .def ( pybind11::init < ProSHADE_settings* > ( ), pybind11::arg (
"settings" ) )
34 .def ( pybind11::init ( [] (
ProSHADE_settings* settings, 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 )
37 pybind11::buffer_info buf = mapData.request();
38 proshade_unsign len = buf.size;
41 double* npVals =
new double[
static_cast<proshade_unsign
> ( len )];
47 for ( proshade_unsign iter = 0; iter < len; iter++ ) { npVals[iter] =
static_cast < double > ( mapData.at(iter) ); }
49 else if ( buf.ndim == 3 )
51 float* dataPtr =
reinterpret_cast < float*
> ( buf.ptr );
52 for ( proshade_unsign xIt = 0; xIt < static_cast<proshade_unsign> ( buf.shape.at(0) ); xIt++ )
54 for ( proshade_unsign yIt = 0; yIt < static_cast<proshade_unsign> ( buf.shape.at(1) ); yIt++ )
56 for ( proshade_unsign zIt = 0; zIt < static_cast<proshade_unsign> ( buf.shape.at(2) ); zIt++ )
58 npVals[zIt + buf.shape.at(2) * ( yIt + buf.shape.at(1) * xIt )] =
static_cast < double > ( dataPtr[zIt + buf.shape.at(2) * ( yIt + buf.shape.at(1) * xIt )] );
65 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;
66 exit ( EXIT_FAILURE );
73 static_cast<int> ( len ),
91 .def (
"writeMap", &
ProSHADE_internal_data::ProSHADE_data::writeMap,
"Function for writing out the internal structure representation in MRC MAP format.", pybind11::arg (
"fname" ), pybind11::arg (
"title" ) =
"Created by ProSHADE and written by GEMMI", pybind11::arg (
"mode" ) = 2 )
92 .def (
"writePdb", &
ProSHADE_internal_data::ProSHADE_data::writePdb,
"This function writes out the PDB formatted file coresponding to the structure.", pybind11::arg (
"fname" ), pybind11::arg (
"euA" ) = 0.0, pybind11::arg (
"euB" ) = 0.0, pybind11::arg (
"euG" ) = 0.0, pybind11::arg (
"trsX" ) = 0.0, pybind11::arg (
"trsY" ) = 0.0, pybind11::arg (
"trsZ" ) = 0.0, pybind11::arg (
"firstModel" ) = true )
97 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( {
self.xDimIndices,
self.yDimIndices,
self.zDimIndices },
98 {
self.yDimIndices *
self.zDimIndices *
sizeof(proshade_double),
99 self.zDimIndices *
sizeof(proshade_double),
100 sizeof(proshade_double) },
108 .def (
"processInternalMap", &
ProSHADE_internal_data::ProSHADE_data::processInternalMap,
"This function simply clusters several map manipulating functions which should be called together. These include centering, phase removal, normalisation, adding extra space, etc.", pybind11::arg (
"settings" ) )
116 .def (
"getReBoxBoundaries",
120 proshade_signed* retVals =
new proshade_signed[6];
126 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { retVals[iter] = settings->
forceBounds[iter]; }
133 self.xDimSize,
self.yDimSize,
self.zDimSize, retVals );
137 self.xDimSize,
self.yDimSize,
self.zDimSize, retVals, settings->
boundsExtraSpace );
143 if ( settings->
useSameBounds && (
self.inputOrder == 0 ) ) {
for ( proshade_unsign iter = 0; iter < 6; iter++ ) { settings->
forceBounds[iter] = retVals[iter]; } }
147 pybind11::capsule pyCapsuleReBox2 ( retVals, [](
void *f ) { proshade_signed* foo =
reinterpret_cast< proshade_signed*
> ( f );
delete foo; } );
150 pybind11::array_t < proshade_signed > retArr = pybind11::array_t < proshade_signed > ( { 6 },
151 {
sizeof(proshade_signed) },
157 },
"This function finds the boundaries enclosing positive map values and adds some extra space." )
158 .def (
"createNewMapFromBounds",
162 proshade_signed* newBounds =
new proshade_signed[6];
166 for ( proshade_unsign iter = 0; iter < 6; iter++ ) { newBounds[iter] = bounds.at(iter); }
169 self.createNewMapFromBounds ( settings, newStr, newBounds );
173 },
"This function creates a new structure from the calling structure and new bounds values." )
190 .def (
"detectSymmetryInStructure",
204 },
"This function runs the symmetry detection algorithms on this structure and saves the results in the settings object.", pybind11::arg (
"settings" ) )
207 .def (
"getRecommendedSymmetryAxes",
211 float* npVals =
new float[
static_cast<unsigned int> ( settings->
detectedSymmetry.size() * 6 )];
215 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->
detectedSymmetry.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 6; it++ ) { npVals[(iter*6)+it] = settings->
detectedSymmetry.at(iter)[it]; } }
218 pybind11::capsule pyCapsuleStrRecSym ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
221 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
detectedSymmetry.size() ),
static_cast<int> ( 6 ) },
222 { 6 *
sizeof(float),
sizeof(
float) },
224 pyCapsuleStrRecSym );
228 },
"This function returns the recommended symmetry axes as a 2D numpy array." )
229 .def (
"getAllCSyms",
233 float* npVals =
new float[
static_cast<unsigned int> ( settings->
allDetectedCAxes.size() * 6 )];
237 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( settings->
allDetectedCAxes.size() ); iter++ ) {
for ( proshade_unsign it = 0; it < 6; it++ ) { npVals[(iter*6)+it] = settings->
allDetectedCAxes.at(iter).at(it); } }
240 pybind11::capsule pyCapsuleStrSymList ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
243 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<int> ( settings->
allDetectedCAxes.size() ), 6 },
244 { 6 *
sizeof(float),
sizeof(
float) },
246 pyCapsuleStrSymList );
250 },
"This function returns all symmetry axes as a 2D numpy array." )
251 .def (
"getNonCSymmetryAxesIndices",
255 pybind11::dict retDict;
256 pybind11::list dList, tList, oList, iList;
259 for ( proshade_unsign dIt = 0; dIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.size() ); dIt++ )
261 pybind11::list memList;
262 for ( proshade_unsign memIt = 0; memIt < static_cast<proshade_unsign> ( settings->
allDetectedDAxes.at(dIt).size() ); memIt++ )
266 dList.append ( memList );
270 for ( proshade_unsign tIt = 0; tIt < static_cast<proshade_unsign> ( settings->
allDetectedTAxes.size() ); tIt++ )
276 for ( proshade_unsign oIt = 0; oIt < static_cast<proshade_unsign> ( settings->
allDetectedOAxes.size() ); oIt++ )
282 for ( proshade_unsign iIt = 0; iIt < static_cast<proshade_unsign> ( settings->
allDetectedIAxes.size() ); iIt++ )
288 retDict[ pybind11::handle ( pybind11::str (
"D" ).ptr ( ) ) ] = dList;
289 retDict[ pybind11::handle ( pybind11::str (
"T" ).ptr ( ) ) ] = tList;
290 retDict[ pybind11::handle ( pybind11::str (
"O" ).ptr ( ) ) ] = oList;
291 retDict[ pybind11::handle ( pybind11::str (
"I" ).ptr ( ) ) ] = iList;
295 },
"This function returns array of non-C axes indices." )
296 .def (
"getAllGroupElements",
300 pybind11::buffer_info buf = axList.request();
301 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 ); }
304 std::vector< proshade_unsign > axesList;
308 std::vector < std::vector< proshade_double > > vals =
self.getAllGroupElements ( settings, axesList, groupType, matrixTolerance );
311 pybind11::list retList;
314 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ )
317 float* npVals =
new float[
static_cast<unsigned int> ( 9 )];
321 for ( proshade_unsign it = 0; it < 9; it++ ) { npVals[it] = vals.at(iter).at(it); }
324 pybind11::capsule pyCapsuleGrpEl ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
327 pybind11::array_t < float > retArr = pybind11::array_t<float> ( { 3, 3 },
328 { 3 *
sizeof(float),
sizeof(
float) },
333 retList.append ( retArr );
338 },
"This function returns the group elements as rotation matrices of any point group described by the detected axes.", pybind11::arg (
"settings" ), pybind11::arg (
"axList" ), pybind11::arg (
"groupType" ) =
"", pybind11::arg(
"matrixTolerance" ) = 0.05 )
342 .def (
"getBestRotationMapPeaksEulerAngles",
346 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
349 float* npVals =
new float[
static_cast<unsigned int> (vals.size())];
353 for ( proshade_unsign iter = 0; iter < static_cast<proshade_unsign> ( vals.size() ); iter++ ) { npVals[iter] = vals.at(iter); }
356 pybind11::capsule pyCapsuleEuPeak ( npVals, [](
void *f ) {
float* foo =
reinterpret_cast< float*
> ( f );
delete foo; } );
359 pybind11::array_t < float > retArr = pybind11::array_t<float> ( {
static_cast<unsigned int> (vals.size()) },
366 },
"This function returns a vector of three floats, the three Euler angles of the best peak in the rotation map.", pybind11::arg (
"settings" ) )
367 .def (
"getBestRotationMapPeaksRotationMatrix",
371 std::vector< proshade_double > vals =
self.getBestRotationMapPeaksEulerAngles ( settings );
374 proshade_double* retMat =
new proshade_double[9];
379 pybind11::capsule pyCapsuleRMPeak ( retMat, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
382 pybind11::array_t < proshade_double > retArr = pybind11::array_t<proshade_double> ( { 3, 3 },
383 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
389 },
"This function returns a rotation matrix representing the best peak in the rotation map.", pybind11::arg (
"settings" ) )
390 .def (
"rotateMap", &
ProSHADE_internal_data::ProSHADE_data::rotateMap,
"This function rotates a map based on the given Euler angles.", pybind11::arg (
"settings" ), pybind11::arg (
"eulerAlpha" ), pybind11::arg (
"eulerBeta" ), pybind11::arg (
"eulerGamma" ) )
393 .def (
"getOverlayTranslations",
397 std::vector< proshade_double > vals =
self.getBestTranslationMapPeaksAngstrom ( staticStructure, eulA, eulB, eulG );
400 pybind11::dict retDict;
401 pybind11::list rotCen, toOverlay;
404 rotCen.append (
self.originalPdbRotCenX );
405 rotCen.append (
self.originalPdbRotCenY );
406 rotCen.append (
self.originalPdbRotCenZ );
408 toOverlay.append (
self.originalPdbTransX );
409 toOverlay.append (
self.originalPdbTransY );
410 toOverlay.append (
self.originalPdbTransZ );
413 retDict[ pybind11::handle ( pybind11::str (
"centreOfRotation" ).ptr ( ) ) ] = rotCen;
414 retDict[ pybind11::handle ( pybind11::str (
"rotCenToOverlay" ).ptr ( ) ) ] = toOverlay;
418 },
"This function returns the vector from optimal rotation centre to origin and the optimal overlay translation vector. These two vectors allow overlaying the inputs (see documentation for details on how the two vectors should be used).", pybind11::arg (
"staticStructure" ), pybind11::arg (
"eulA" ), pybind11::arg (
"eulB" ), pybind11::arg (
"eulG" ) )
419 .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 (
"settings" ), pybind11::arg (
"trsX" ), pybind11::arg (
"trsY" ), pybind11::arg (
"trsZ" ) )
422 .def (
"findSHIndex",
426 proshade_signed index = seanindex ( order, band,
self.spheres[shell]->getLocalBandwidth() );
430 },
"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." )
431 .def (
"getSphericalHarmonics",
435 std::complex<proshade_double>* npVals =
new std::complex<proshade_double>[
static_cast<proshade_unsign
> (
self.noSpheres * pow(
self.maxShellBand, 2.0 ) )];
439 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 ); }
442 proshade_signed pyPosSH;
443 proshade_signed pyPos;
446 for ( proshade_signed shIt = 0; shIt < static_cast<proshade_signed> (
self.noSpheres ); shIt++ )
448 for ( proshade_signed bnd = 0; bnd < static_cast<proshade_signed> (
self.spheres[shIt]->getLocalBandwidth() ); bnd++ )
450 for ( proshade_signed order = -bnd; order <= bnd; order++ )
452 pyPosSH = ( shIt * pow(
self.maxShellBand, 2.0 ) );
453 pyPos = seanindex ( order, bnd,
self.spheres[shIt]->getLocalBandwidth() );
454 npVals[pyPosSH+pyPos].real (
self.sphericalHarmonics[shIt][pyPos][0] );
455 npVals[pyPosSH+pyPos].imag (
self.sphericalHarmonics[shIt][pyPos][1] );
461 pybind11::capsule pyCapsuleSHs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
464 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
465 {
static_cast<int> (
self.noSpheres ),
static_cast<int> ( pow (
self.maxShellBand, 2.0 ) ) },
466 {
sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( pow (
self.maxShellBand, 2.0 ) ),
sizeof ( std::complex < proshade_double > ) },
472 },
"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." )
477 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) )];
481 proshade_signed index = 0;
484 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 ); }
487 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
489 for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
491 for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
493 index = order2 + ( (
self.maxShellBand * 2 ) + 1 ) * ( order1 + ( (
self.maxShellBand * 2 ) + 1 ) * bandIter );
494 npVals[index].real (
self.eMatrices[bandIter][order1][order2][0] );
495 npVals[index].imag (
self.eMatrices[bandIter][order1][order2][1] );
501 pybind11::capsule pyCapsuleEMs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
504 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
505 {
self.maxShellBand, ( (
self.maxShellBand * 2 ) + 1 ), ( (
self.maxShellBand * 2 ) + 1 ) },
506 {
sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) ),
507 sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( (
self.maxShellBand * 2 ) + 1 ),
508 sizeof ( std::complex < proshade_double > ) },
514 },
"This function returns the E matrix values (these are the integral over all spheres of ( c1^(l,m) * c2^(l,m') ) values) obtained when rotation function or self-rotation function are computed. The returned array has three dimensions, first being the band, second being the order1 and third being the order2. Please note that as negative indexing does not simply work, the order indexing starts from 0 - i.e. array[1][0][0] means band 1 ; order1 = -1 and order2 = -1." )
515 .def (
"getSO3Coefficients",
519 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.maxShellBand * ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) )];
523 proshade_signed index = 0;
526 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 ); }
529 for ( proshade_signed bandIter = 0; bandIter < static_cast< proshade_signed > (
self.maxShellBand ); bandIter++ )
531 for ( proshade_unsign order1 = 0; order1 < ( ( bandIter * 2 ) + 1 ); order1++ )
533 for ( proshade_unsign order2 = 0; order2 < ( ( bandIter * 2 ) + 1 ); order2++ )
535 index = order2 + ( (
self.maxShellBand * 2 ) + 1 ) * ( order1 + ( (
self.maxShellBand * 2 ) + 1 ) * bandIter );
536 npVals[index].real (
self.so3Coeffs[
self.so3CoeffsArrayIndex ( order1 - bandIter, order2 - bandIter, bandIter )][0] );
537 npVals[index].imag (
self.so3Coeffs[
self.so3CoeffsArrayIndex ( order1 - bandIter, order2 - bandIter, bandIter )][1] );
543 pybind11::capsule pyCapsuleSOCoeffs ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
546 pybind11::array_t < std::complex<proshade_double> > retArr = pybind11::array_t < std::complex<proshade_double > > (
547 {
self.maxShellBand, ( (
self.maxShellBand * 2 ) + 1 ), ( (
self.maxShellBand * 2 ) + 1 ) },
548 {
sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( ( (
self.maxShellBand * 2 ) + 1 ) * ( (
self.maxShellBand * 2 ) + 1 ) ),
549 sizeof ( std::complex < proshade_double > ) *
static_cast<int> ( (
self.maxShellBand * 2 ) + 1 ),
550 sizeof ( std::complex < proshade_double > ) },
556 },
"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." )
557 .def (
"getRotationFunctionMap",
561 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> ( (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) )];
565 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] ); }
568 pybind11::capsule pyCapsuleRotMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
571 pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
572 { (
self.maxShellBand * 2 ), (
self.maxShellBand * 2 ), (
self.maxShellBand * 2 ) },
573 { (
self.maxShellBand * 2 ) * (
self.maxShellBand * 2 ) *
sizeof(std::complex < proshade_double >),
574 (
self.maxShellBand * 2 ) *
sizeof(std::complex < proshade_double >),
575 sizeof(std::complex < proshade_double >) },
581 },
"This function returns the (self) rotation function as a three-dimensional map of complex numbers." )
582 .def (
"getRotationMatrixFromSOFTCoordinates",
586 proshade_double* npVals =
new proshade_double[9];
590 proshade_double eulA, eulB, eulG;
599 pybind11::capsule pyCapsuleRMSoft ( npVals, [](
void *f ) { proshade_double* foo =
reinterpret_cast< proshade_double*
> ( f );
delete foo; } );
602 pybind11::array_t < proshade_double > retArr = pybind11::array_t < proshade_double > ( { 3, 3 },
603 { 3 *
sizeof(proshade_double),
sizeof(proshade_double) },
609 },
"This function converts a given rotation function map position onto the corresponding rotation matrix.", pybind11::arg (
"xPos" ), pybind11::arg (
"yPos" ), pybind11::arg (
"zPos" ) )
610 .def (
"getTranslationFunctionMap",
614 std::complex<proshade_double>* npVals =
new std::complex < proshade_double >[
static_cast<proshade_unsign
> (
self.getXDim() *
self.getYDim() *
self.getZDim() )];
618 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] ); }
621 pybind11::capsule pyCapsuleTrsMap ( npVals, [](
void *f ) { std::complex<proshade_double>* foo =
reinterpret_cast< std::complex<proshade_double>*
> ( f );
delete foo; } );
624 pybind11::array_t < std::complex < proshade_double > > retArr = pybind11::array_t < std::complex < proshade_double > > (
625 {
self.getXDim(),
self.getYDim(),
self.getZDim() },
626 {
self.getYDim() *
self.getZDim() *
sizeof(std::complex < proshade_double >),
627 self.getZDim() *
sizeof(std::complex < proshade_double >),
628 sizeof(std::complex < proshade_double >) },
634 },
"This function returns the translation function as a three-dimensional map of complex numbers." )
672 .def (
"__repr__", [] (
const ProSHADE_internal_data::ProSHADE_data &a ) {
return "<ProSHADE_data class object> (This class contains all information, results and available functionalities for a structure)"; } );