ProSHADE  0.7.6.1 (AUG 2021)
Protein Shape Detection
ProSHADE_internal_io Namespace Reference

This namespace contains the internal input/output functions. None of these should be directly accessed by the user. More...

Enumerations

enum  InputType { UNKNOWN, PDB, MAP, GEMMI }
 

Functions

InputType figureDataType (std::string fName)
 Function determining input data type. More...
 
bool isFilePDB (std::string fName)
 Function determining if the input data type is PDB. More...
 
bool isFileMAP (std::string fName)
 Function determining if the input data type is MAP. More...
 
void readInMapHeader (gemmi::Ccp4< float > *map, proshade_unsign *xDimInds, proshade_unsign *yDimInds, proshade_unsign *zDimInds, proshade_single *xDim, proshade_single *yDim, proshade_single *zDim, proshade_single *aAng, proshade_single *bAng, proshade_single *cAng, proshade_signed *xFrom, proshade_signed *yFrom, proshade_signed *zFrom, proshade_signed *xAxOrigin, proshade_signed *yAxOrigin, proshade_signed *zAxOrigin, proshade_unsign *xAxOrder, proshade_unsign *yAxOrder, proshade_unsign *zAxOrder, proshade_unsign *xGridInds, proshade_unsign *yGridInds, proshade_unsign *zGridInds)
 This function parses the CCP4 MAP file header as read in by gemmi. More...
 
void readInMapData (gemmi::Ccp4< float > *gemmiMap, proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder)
 This function converts the gemmi Ccp4 object data to ProSHADE internal map representation. More...
 
void readInMapData (gemmi::Ccp4< int8_t > *gemmiMap, proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder)
 This function converts the gemmi Ccp4 object data to ProSHADE mask representation. More...
 
void applyMask (proshade_double *&map, std::string maskFile, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_signed verbose, proshade_double *maskArray=nullptr, proshade_unsign maXInds=0, proshade_unsign maYInds=0, proshade_unsign maZInds=0)
 This function reads and applies the mask to the map. More...
 
void applyMaskFromArray (proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&mask, proshade_unsign xDimIndsMsk, proshade_unsign yDimIndsMsk, proshade_unsign zDimIndsMsk, proshade_signed verbose)
 This function applies the mask to the map. More...
 
void applyWeights (proshade_double *&map, std::string weightsFile, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_signed verbose, proshade_double *weightsArray=nullptr, proshade_unsign waXInds=0, proshade_unsign waYInds=0, proshade_unsign waZInds=0)
 This function reads and applies the Fourier weights to the map. More...
 
void applyWeightsFromArray (proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&weights, proshade_unsign xDimIndsWgh, proshade_unsign yDimIndsWgh, proshade_unsign zDimIndsWgh, proshade_signed verbose)
 This function applies the weights to the map. More...
 
void writeOutMapHeader (gemmi::Ccp4< float > *map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_single xDim, proshade_single yDim, proshade_single zDim, proshade_single aAng, proshade_single bAng, proshade_single cAng, proshade_signed xFrom, proshade_signed yFrom, proshade_signed zFrom, proshade_signed xAxOrigin, proshade_signed yAxOrigin, proshade_signed zAxOrigin, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder, proshade_unsign xGridInds, proshade_unsign yGridInds, proshade_unsign zGridInds, std::string title, int mode)
 This function parses the CCP4 MAP file header as read in by gemmi. More...
 
void writeRotationTranslationJSON (proshade_double trsX1, proshade_double trsY1, proshade_double trsZ1, proshade_double eulA, proshade_double eulB, proshade_double eulG, proshade_double trsX2, proshade_double trsY2, proshade_double trsZ2, std::string fileName)
 Function for writing out the optimal rotation and translation into a JSON file. More...
 

Detailed Description

This namespace contains the internal input/output functions. None of these should be directly accessed by the user.

The ProSHADE_internal_io namespace contains the helper functions for the data input and output. These should never be directly used by the user and these only serve to allow for self-documenting nature of the code. They are called internally by more advanced functions from the higher complexity classes.

Function Documentation

◆ applyMask()

void ProSHADE_internal_io::applyMask ( proshade_double *&  map,
std::string  maskFile,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_signed  verbose,
proshade_double *  maskArray = nullptr,
proshade_unsign  maXInds = 0,
proshade_unsign  maYInds = 0,
proshade_unsign  maZInds = 0 
)

This function reads and applies the mask to the map.

This function reads in the mask file and determines, if its dimensions are the same as the input map or not. If not, then the mask is re-sampled is such a way as to have the same number of indices as the input map. Once the map and mask dimensions are the same, the mask is applied by multiplying all map values by corresponding mask value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]maskFileFilename of where the mask is to be read from.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]verboseHow much std::out output would you like?
[in]maskArrayAn array of mask values (default nullptr) to be used instead of an input file.
[in]maXIndsThe size of maskArray x dimension in indices (defaults to 0).
[in]maYIndsThe size of maskArray y dimension in indices (defaults to 0).
[in]maZIndsThe size of maskArray z dimension in indices (defaults to 0).

Definition at line 294 of file ProSHADE_io.cpp.

295 {
296  //================================================ Report progress
297  std::stringstream hlpSS;
298  hlpSS << "Reading mask file " << maskFile;
299  ProSHADE_internal_messages::printProgressMessage ( verbose, 2, hlpSS.str() );
300 
301  //================================================ Are we reading from array or from file?
302  if ( ( maskArray != nullptr ) && ( maXInds != 0 ) && ( maYInds != 0 ) && ( maZInds != 0 ) )
303  {
304  //============================================ Array it is!
305  ProSHADE_internal_io::applyMaskFromArray ( map, xDimInds, yDimInds, zDimInds, maskArray, maXInds, maYInds, maZInds, verbose );
306  }
307  else
308  {
309  //============================================ Check if filename was given
310  if ( maskFile == "" ) { return; }
311 
312  //============================================ File it is! Open the mask
313  gemmi::Ccp4<float> mask;
314  mask.read_ccp4 ( gemmi::MaybeGzipped ( maskFile.c_str() ) );
315 
316  //============================================ Convert to XYZ and create complete mask, if need be
317  mask.setup ( gemmi::GridSetup::ReorderOnly, 0 );
318 
319  //============================================ Read in the rest of the mask file header
320  proshade_unsign xDI, yDI, zDI, xAOR, yAOR, zAOR, xGI, yGI, zGI;
321  proshade_single xDS, yDS, zDS, aA, bA, cA;
322  proshade_signed xF, yF, zF, xAO, yAO, zAO;
324  &xDI, &yDI, &zDI,
325  &xDS, &yDS, &zDS,
326  &aA, &bA, &cA,
327  &xF, &yF, &zF,
328  &xAO, &yAO, &zAO,
329  &xAOR, &yAOR, &zAOR,
330  &xGI, &yGI, &zGI );
331 
332  //============================================ Save the mask values to ProSHADE variable
333  proshade_double* internalMask = nullptr;
334  ProSHADE_internal_io::readInMapData ( &mask, internalMask, xDI, yDI, zDI, xAOR, yAOR, zAOR );
335 
336  //============================================ Apply mask from array
337  ProSHADE_internal_io::applyMaskFromArray ( map, xDimInds, yDimInds, zDimInds, internalMask, xDI, yDI, zDI, verbose );
338 
339  //============================================ Release the memory
340  delete[] internalMask;
341  }
342 
343  //================================================ Done
344  return ;
345 
346 }

◆ applyMaskFromArray()

void ProSHADE_internal_io::applyMaskFromArray ( proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_double *&  mask,
proshade_unsign  xDimIndsMsk,
proshade_unsign  yDimIndsMsk,
proshade_unsign  zDimIndsMsk,
proshade_signed  verbose 
)

This function applies the mask to the map.

This function firstlydetermines if its dimensions of mask are the same as the input map or not. If not, then the mask is re-sampled is such a way as to have the same number of indices as the input map. Once the map and mask dimensions are the same, the mask is applied by multiplying all map values by corresponding mask value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]xDimIndsThe size of the map x dimension in indices.
[in]yDimIndsThe size of the map y dimension in indices.
[in]zDimIndsThe size of the map z dimension in indices.
[in]maskA pointer to 1D array of the mask values.
[in]xDimIndsMskThe size of the mask x dimension in indices.
[in]yDimIndsMskThe size of the mask y dimension in indices.
[in]zDimIndsMskThe size of the mask z dimension in indices.
[in]verboseHow much std::out output would you like?

Definition at line 365 of file ProSHADE_io.cpp.

366 {
367  //================================================ Initialise local variables
368  size_t origVolume = xDimInds * yDimInds * zDimInds;
369  size_t newVolume = xDimIndsMsk * yDimIndsMsk * zDimIndsMsk;
370  proshade_double* maskFinal;
371 
372  //================================================ If mask has different number of indices than map, then re-sample mask
373  if ( ( xDimIndsMsk != xDimInds ) || ( yDimIndsMsk != yDimInds ) || ( zDimIndsMsk != zDimInds ) )
374  {
375  //============================================ Initialise variables
376  fftw_complex* origCoeffs = new fftw_complex [newVolume ];
377  fftw_complex* origCoeffsHKL = new fftw_complex [newVolume ];
378  fftw_complex* modifCoeffs = new fftw_complex [origVolume];
379  fftw_complex* modifCoeffsHKL = new fftw_complex [origVolume];
380  fftw_complex* inMap = new fftw_complex [newVolume ];
381  fftw_complex* outMap = new fftw_complex [origVolume];
382 
383  //============================================ Check memory allocation
384  ProSHADE_internal_misc::checkMemoryAllocation ( inMap, __FILE__, __LINE__, __func__ );
385  ProSHADE_internal_misc::checkMemoryAllocation ( outMap, __FILE__, __LINE__, __func__ );
386  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffs, __FILE__, __LINE__, __func__ );
387  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffs, __FILE__, __LINE__, __func__ );
388  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffsHKL, __FILE__, __LINE__, __func__ );
389  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffsHKL, __FILE__, __LINE__, __func__ );
390 
391  //============================================ Set array to zeroes
392  for ( size_t iter = 0; iter < origVolume; iter++ ) { modifCoeffsHKL[iter][0] = 0.0; modifCoeffsHKL[iter][1] = 0.0; }
393 
394  //============================================ Cope mask to Fourier input array
395  for ( size_t iter = 0; iter < newVolume; iter++ ) { inMap[iter][0] = mask[iter]; inMap[iter][1] = 0.0; }
396 
397  //============================================ Prepare Fourier transform plans
398  fftw_plan planForwardFourier = fftw_plan_dft_3d ( static_cast< int > ( xDimIndsMsk ), static_cast< int > ( yDimIndsMsk ), static_cast< int > ( zDimIndsMsk ), inMap, origCoeffs, FFTW_FORWARD, FFTW_ESTIMATE );
399  fftw_plan inverseFoourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), modifCoeffs, outMap, FFTW_BACKWARD, FFTW_ESTIMATE );
400 
401  //============================================ Compute pre and post changes
402  proshade_signed xPre, yPre, zPre;
403  xPre = std::abs ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsMsk ) ) / 2 );
404  yPre = std::abs ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsMsk ) ) / 2 );
405  zPre = std::abs ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsMsk ) ) / 2 );
406 
407  if ( ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsMsk ) ) % 2 ) == 1 ) { xPre -= 1; }
408  if ( ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsMsk ) ) % 2 ) == 1 ) { yPre -= 1; }
409  if ( ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsMsk ) ) % 2 ) == 1 ) { zPre -= 1; }
410 
411  //============================================ Run forward Fourier
412  fftw_execute ( planForwardFourier );
413 
414  //============================================ Initialise local variables
415  proshade_signed maskMapIndex = 0;
416  proshade_signed densMapIndex = 0;
417  proshade_signed xMaskPos, yMaskPos, zMaskPos, xDensPos, yDensPos, zDensPos;
418  proshade_signed maskH, maskK, maskL;
419 
420  //============================================ Convert mask to HKL for re-boxing
421  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimIndsMsk ); xIt++ )
422  {
423  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimIndsMsk ); yIt++ )
424  {
425  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimIndsMsk ); zIt++ )
426  {
427  //================================ Convert to HKL
428  maskH = xIt + static_cast< proshade_signed > ( (xDimIndsMsk+1) / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimIndsMsk ) ) { maskH -= xDimIndsMsk; }
429  maskK = yIt + static_cast< proshade_signed > ( (yDimIndsMsk+1) / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimIndsMsk ) ) { maskK -= yDimIndsMsk; }
430  maskL = zIt + static_cast< proshade_signed > ( (zDimIndsMsk+1) / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimIndsMsk ) ) { maskL -= zDimIndsMsk; }
431 
432  //================================ Find the positions
433  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimIndsMsk ) * ( yIt + static_cast< proshade_signed > ( yDimIndsMsk ) * xIt );
434  densMapIndex = maskL + static_cast< proshade_signed > ( zDimIndsMsk ) * ( maskK + static_cast< proshade_signed > ( yDimIndsMsk ) * maskH );
435 
436  //================================ Save the values
437  origCoeffsHKL[densMapIndex][0] = origCoeffs[maskMapIndex][0];
438  origCoeffsHKL[densMapIndex][1] = origCoeffs[maskMapIndex][1];
439  }
440  }
441  }
442 
443  //============================================ Rebox
444  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
445  {
446  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
447  {
448  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
449  {
450  //================================ Deal with X
451  if ( xDimIndsMsk >= xDimInds ) { xMaskPos = xIt + xPre; }
452  else { xMaskPos = xIt - xPre; }
453  xDensPos = xIt;
454 
455  //================================ Deal with Y
456  if ( yDimIndsMsk >= yDimInds ) { yMaskPos = yIt + yPre; }
457  else { yMaskPos = yIt - yPre; }
458  yDensPos = yIt;
459 
460  //================================ Deal with Z
461  if ( zDimIndsMsk >= zDimInds ) { zMaskPos = zIt + zPre; }
462  else { zMaskPos = zIt - zPre; }
463  zDensPos = zIt;
464 
465  //================================ Skip if mask value not available (because the modifCoeffsHKL array is zeroed, we do not need to do anything here)
466  if ( ( xMaskPos < 0 ) || ( xMaskPos >= static_cast< proshade_signed > ( xDimIndsMsk ) ) ) { continue; }
467  if ( ( yMaskPos < 0 ) || ( yMaskPos >= static_cast< proshade_signed > ( yDimIndsMsk ) ) ) { continue; }
468  if ( ( zMaskPos < 0 ) || ( zMaskPos >= static_cast< proshade_signed > ( zDimIndsMsk ) ) ) { continue; }
469 
470  //================================ Find the positions
471  maskMapIndex = zMaskPos + static_cast< proshade_signed > ( zDimIndsMsk ) * ( yMaskPos + static_cast< proshade_signed > ( yDimIndsMsk ) * xMaskPos );
472  densMapIndex = zDensPos + static_cast< proshade_signed > ( zDimInds ) * ( yDensPos + static_cast< proshade_signed > ( yDimInds ) * xDensPos );
473 
474  //================================ Copy values
475  modifCoeffsHKL[densMapIndex][0] = origCoeffsHKL[maskMapIndex][0];
476  modifCoeffsHKL[densMapIndex][1] = origCoeffsHKL[maskMapIndex][1];
477  }
478  }
479  }
480 
481  //============================================ Convert mask back to FFTW order
482  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
483  {
484  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
485  {
486  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
487  {
488  //================================ Convert to HKL
489  maskH = xIt + static_cast< proshade_signed > ( xDimInds / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimInds ) ) { maskH -= xDimInds; }
490  maskK = yIt + static_cast< proshade_signed > ( yDimInds / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimInds ) ) { maskK -= yDimInds; }
491  maskL = zIt + static_cast< proshade_signed > ( zDimInds / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimInds ) ) { maskL -= zDimInds; }
492 
493  //================================ Find the positions
494  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimInds ) * ( yIt + static_cast< proshade_signed > ( yDimInds ) * xIt );
495  densMapIndex = maskL + static_cast< proshade_signed > ( zDimInds ) * ( maskK + static_cast< proshade_signed > ( yDimInds ) * maskH );
496 
497  //================================ Save the values
498  modifCoeffs[densMapIndex][0] = modifCoeffsHKL[maskMapIndex][0];
499  modifCoeffs[densMapIndex][1] = modifCoeffsHKL[maskMapIndex][1];
500  }
501  }
502  }
503 
504  //============================================ Run inverse Fourier on the modified coefficients
505  fftw_execute ( inverseFoourier );
506 
507  //============================================ Delete old mask and allocate memory for the new, re-sampled mask
508  maskFinal = new proshade_double [origVolume];
509  ProSHADE_internal_misc::checkMemoryAllocation ( maskFinal, __FILE__, __LINE__, __func__ );
510 
511  //======================================== Copy results into a new, properly sampled mask
512  for ( size_t iter = 0; iter < origVolume; iter++ ) { maskFinal[iter] = outMap[iter][0]; }
513 
514  //======================================== Release remaining memory
515  fftw_destroy_plan ( planForwardFourier );
516  fftw_destroy_plan ( inverseFoourier );
517  delete[] origCoeffs;
518  delete[] modifCoeffs;
519  delete[] origCoeffsHKL;
520  delete[] modifCoeffsHKL;
521  delete[] inMap;
522  delete[] outMap;
523  }
524  else
525  {
526  maskFinal = new proshade_double [origVolume];
527  ProSHADE_internal_misc::checkMemoryAllocation ( maskFinal, __FILE__, __LINE__, __func__ );
528  for ( size_t iter = 0; iter < origVolume; iter++ ) { maskFinal[iter] = mask[iter]; }
529  }
530 
531  //================================================ Apply the mask to the map
532  for ( size_t iter = 0; iter < static_cast< size_t > ( xDimInds * yDimInds * zDimInds ); iter++ ) { map[iter] *= maskFinal[iter]; }
533 
534  //================================================ Release memory
535  delete[] maskFinal;
536 
537  //================================================ Report progress
538  ProSHADE_internal_messages::printProgressMessage ( verbose, 3, "Mask read in and applied successfully." );
539 
540  //================================================ Done
541  return ;
542 
543 }

◆ applyWeights()

void ProSHADE_internal_io::applyWeights ( proshade_double *&  map,
std::string  weightsFile,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_signed  verbose,
proshade_double *  weightsArray = nullptr,
proshade_unsign  waXInds = 0,
proshade_unsign  waYInds = 0,
proshade_unsign  waZInds = 0 
)

This function reads and applies the Fourier weights to the map.

This function reads in the weights file and determines, if its dimensions are the same as the input map or not. If not, then the weights are re-sampled is such a way as to have the same number of indices as the input map. Once the map and weights map dimensions are the same, the weights are applied by multiplying all map values by corresponding weight value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]weightsFileFilename of where the mask is to be read from.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]verboseHow much std::out output would you like?
[in]weightsArrayAn array of weights (default nullptr) to be used instead of input file.
[in]waXIndsThe size of weightsArray x dimension in indices (defaults to 0).
[in]waYIndsThe size of weightsArray y dimension in indices (defaults to 0).
[in]waZIndsThe size of weightsArray z dimension in indices (defaults to 0).

Definition at line 563 of file ProSHADE_io.cpp.

564 {
565  //================================================ Report progress
566  std::stringstream hlpSS;
567  ProSHADE_internal_messages::printProgressMessage ( verbose, 2, hlpSS.str() );
568 
569  //================================================ Are we reading from file, or from array?
570  if ( ( weightsArray != nullptr ) && ( waXInds != 0 ) && ( waYInds != 0 ) && ( waZInds != 0 ) )
571  {
572  //============================================ From array it is!
573  ProSHADE_internal_io::applyWeightsFromArray ( map, xDimInds, yDimInds, zDimInds, weightsArray, waXInds, waYInds, waZInds, verbose );
574  }
575  else
576  {
577  //============================================ Check if weights file was given
578  if ( weightsFile == "" ) { return; }
579 
580  //============================================ From file it is! Open the weights file
581  gemmi::Ccp4<float> weights;
582  weights.read_ccp4 ( gemmi::MaybeGzipped ( weightsFile.c_str() ) );
583 
584  //============================================ Convert to XYZ and create complete weights, if need be
585  weights.setup ( gemmi::GridSetup::ReorderOnly, 0 );
586 
587  //============================================ Read in the rest of the weights file header
588  proshade_unsign xDI, yDI, zDI, xAOR, yAOR, zAOR, xGI, yGI, zGI;
589  proshade_single xDS, yDS, zDS, aA, bA, cA;
590  proshade_signed xF, yF, zF, xAO, yAO, zAO;
592  &xDI, &yDI, &zDI,
593  &xDS, &yDS, &zDS,
594  &aA, &bA, &cA,
595  &xF, &yF, &zF,
596  &xAO, &yAO, &zAO,
597  &xAOR, &yAOR, &zAOR,
598  &xGI, &yGI, &zGI );
599 
600  //============================================ Save the weights values to local variable
601  proshade_double* internalWeights = nullptr;
602  ProSHADE_internal_io::readInMapData ( &weights, internalWeights, xDI, yDI, zDI, xAOR, yAOR, zAOR );
603 
604  //============================================ Apply weights from array
605  ProSHADE_internal_io::applyWeightsFromArray ( map, xDimInds, yDimInds, zDimInds, internalWeights, xDI, yDI, zDI, verbose );
606 
607  //============================================ Release the memory
608  delete[] internalWeights;
609  }
610 
611  //================================================ Done
612  return ;
613 
614 }

◆ applyWeightsFromArray()

void ProSHADE_internal_io::applyWeightsFromArray ( proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_double *&  weights,
proshade_unsign  xDimIndsWgh,
proshade_unsign  yDimIndsWgh,
proshade_unsign  zDimIndsWgh,
proshade_signed  verbose 
)

This function applies the weights to the map.

This function firstly determines if its dimensions of weights are the same as the input map or not. If not, then the weights are re-sampled is such a way as to have the same number of indices as the input map. Once the map and weights dimensions are the same, the weights are applied by multiplying all map Fourier coefficients by corresponding weights value. The map change is done in-place!

Parameters
[in]mapPointer to a the internal variable containing the map (it will be changed in place).
[in]xDimIndsThe size of the map x dimension in indices.
[in]yDimIndsThe size of the map y dimension in indices.
[in]zDimIndsThe size of the map z dimension in indices.
[in]weightsA pointer to 1D array of the weights values.
[in]xDimIndsWghThe size of the weights x dimension in indices.
[in]yDimIndsWghThe size of the weights y dimension in indices.
[in]zDimIndsWghThe size of the weights z dimension in indices.
[in]verboseHow much std::out output would you like?

Definition at line 633 of file ProSHADE_io.cpp.

634 {
635  //================================================ Initialise local variables
636  proshade_double* weightsFinal;
637  size_t origVolume = xDimInds * yDimInds * zDimInds;
638  size_t newVolume = xDimIndsWgh * yDimIndsWgh * zDimIndsWgh;
639 
640  //================================================ If weights have different number of indices than map, then re-sample weights in supplied space
641  if ( ( xDimIndsWgh != xDimInds ) || ( yDimIndsWgh != yDimInds ) || ( zDimIndsWgh != zDimInds ) )
642  {
643  //============================================ Initialise variables
644  fftw_complex* origCoeffs = new fftw_complex [newVolume ];
645  fftw_complex* origCoeffsHKL = new fftw_complex [newVolume ];
646  fftw_complex* modifCoeffs = new fftw_complex [origVolume];
647  fftw_complex* modifCoeffsHKL = new fftw_complex [origVolume];
648  fftw_complex* inMap = new fftw_complex [newVolume ];
649  fftw_complex* outMap = new fftw_complex [origVolume];
650 
651  //============================================ Check memory allocation
652  ProSHADE_internal_misc::checkMemoryAllocation ( inMap, __FILE__, __LINE__, __func__ );
653  ProSHADE_internal_misc::checkMemoryAllocation ( outMap, __FILE__, __LINE__, __func__ );
654  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffs, __FILE__, __LINE__, __func__ );
655  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffs, __FILE__, __LINE__, __func__ );
656  ProSHADE_internal_misc::checkMemoryAllocation ( origCoeffsHKL, __FILE__, __LINE__, __func__ );
657  ProSHADE_internal_misc::checkMemoryAllocation ( modifCoeffsHKL, __FILE__, __LINE__, __func__ );
658 
659  //============================================ Set array to zeroes
660  for ( size_t iter = 0; iter < origVolume; iter++ ) { modifCoeffsHKL[iter][0] = 0.0; modifCoeffsHKL[iter][1] = 0.0; }
661 
662  //============================================ Copy weights to Fourier input array
663  for ( size_t iter = 0; iter < newVolume; iter++ ) { inMap[iter][0] = weights[iter]; inMap[iter][1] = 0.0; }
664 
665  //============================================ Prepare Fourier transform plans
666  fftw_plan planForwardFourier = fftw_plan_dft_3d ( static_cast< int > ( xDimIndsWgh ), static_cast< int > ( yDimIndsWgh ), static_cast< int > ( zDimIndsWgh ), inMap, origCoeffs, FFTW_FORWARD, FFTW_ESTIMATE );
667  fftw_plan inverseFoourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), modifCoeffs, outMap, FFTW_BACKWARD, FFTW_ESTIMATE );
668 
669  //============================================ Compute pre and post changes
670  proshade_signed xPre, yPre, zPre;
671  xPre = std::abs ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsWgh ) ) / 2 );
672  yPre = std::abs ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsWgh ) ) / 2 );
673  zPre = std::abs ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsWgh ) ) / 2 );
674 
675  if ( ( ( static_cast< proshade_signed > ( xDimInds ) - static_cast< proshade_signed > ( xDimIndsWgh ) ) % 2 ) == 1 ) { xPre -= 1; }
676  if ( ( ( static_cast< proshade_signed > ( yDimInds ) - static_cast< proshade_signed > ( yDimIndsWgh ) ) % 2 ) == 1 ) { yPre -= 1; }
677  if ( ( ( static_cast< proshade_signed > ( zDimInds ) - static_cast< proshade_signed > ( zDimIndsWgh ) ) % 2 ) == 1 ) { zPre -= 1; }
678 
679  //============================================ Run forward Fourier
680  fftw_execute ( planForwardFourier );
681 
682  //============================================ Initialise local variables
683  proshade_signed maskMapIndex = 0;
684  proshade_signed densMapIndex = 0;
685  proshade_signed xMaskPos, yMaskPos, zMaskPos, xDensPos, yDensPos, zDensPos;
686  proshade_signed maskH, maskK, maskL;
687 
688  //============================================ Convert weights to HKL for re-boxing
689  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimIndsWgh ); xIt++ )
690  {
691  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimIndsWgh ); yIt++ )
692  {
693  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimIndsWgh ); zIt++ )
694  {
695  //================================ Convert to HKL
696  maskH = xIt + static_cast< proshade_signed > ( (xDimIndsWgh+1) / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimIndsWgh ) ) { maskH -= xDimIndsWgh; }
697  maskK = yIt + static_cast< proshade_signed > ( (yDimIndsWgh+1) / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimIndsWgh ) ) { maskK -= yDimIndsWgh; }
698  maskL = zIt + static_cast< proshade_signed > ( (zDimIndsWgh+1) / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimIndsWgh ) ) { maskL -= zDimIndsWgh; }
699 
700  //================================ Find the positions
701  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimIndsWgh ) * ( yIt + static_cast< proshade_signed > ( yDimIndsWgh ) * xIt );
702  densMapIndex = maskL + static_cast< proshade_signed > ( zDimIndsWgh ) * ( maskK + static_cast< proshade_signed > ( yDimIndsWgh ) * maskH );
703 
704  //================================ Save the values
705  origCoeffsHKL[densMapIndex][0] = origCoeffs[maskMapIndex][0];
706  origCoeffsHKL[densMapIndex][1] = origCoeffs[maskMapIndex][1];
707  }
708  }
709  }
710 
711  //============================================ Rebox
712  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
713  {
714  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
715  {
716  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
717  {
718  //================================ Deal with X
719  if ( xDimIndsWgh >= xDimInds ) { xMaskPos = xIt + xPre; }
720  else { xMaskPos = xIt - xPre; }
721  xDensPos = xIt;
722 
723  //================================ Deal with Y
724  if ( yDimIndsWgh >= yDimInds ) { yMaskPos = yIt + yPre; }
725  else { yMaskPos = yIt - yPre; }
726  yDensPos = yIt;
727 
728  //================================ Deal with Z
729  if ( zDimIndsWgh >= zDimInds ) { zMaskPos = zIt + zPre; }
730  else { zMaskPos = zIt - zPre; }
731  zDensPos = zIt;
732 
733  //================================ Skip if weights value not available (because the modifCoeffsHKL array is zeroed, we do not need to do anything here)
734  if ( ( xMaskPos < 0 ) || ( xMaskPos >= static_cast< proshade_signed > ( xDimIndsWgh ) ) ) { continue; }
735  if ( ( yMaskPos < 0 ) || ( yMaskPos >= static_cast< proshade_signed > ( yDimIndsWgh ) ) ) { continue; }
736  if ( ( zMaskPos < 0 ) || ( zMaskPos >= static_cast< proshade_signed > ( zDimIndsWgh ) ) ) { continue; }
737 
738  //================================ Find the positions
739  maskMapIndex = zMaskPos + static_cast< proshade_signed > ( zDimIndsWgh ) * ( yMaskPos + static_cast< proshade_signed > ( yDimIndsWgh ) * xMaskPos );
740  densMapIndex = zDensPos + static_cast< proshade_signed > ( zDimInds ) * ( yDensPos + static_cast< proshade_signed > ( yDimInds ) * xDensPos );
741 
742  //================================ Copy values
743  modifCoeffsHKL[densMapIndex][0] = origCoeffsHKL[maskMapIndex][0];
744  modifCoeffsHKL[densMapIndex][1] = origCoeffsHKL[maskMapIndex][1];
745  }
746  }
747  }
748 
749  //============================================ Convert weights back to FFTW order
750  for ( proshade_signed xIt = 0; xIt < static_cast< proshade_signed > ( xDimInds ); xIt++ )
751  {
752  for ( proshade_signed yIt = 0; yIt < static_cast< proshade_signed > ( yDimInds ); yIt++ )
753  {
754  for ( proshade_signed zIt = 0; zIt < static_cast< proshade_signed > ( zDimInds ); zIt++ )
755  {
756  //================================ Convert to HKL
757  maskH = xIt + static_cast< proshade_signed > ( xDimInds / 2 ); if ( maskH >= static_cast< proshade_signed > ( xDimInds ) ) { maskH -= xDimInds; }
758  maskK = yIt + static_cast< proshade_signed > ( yDimInds / 2 ); if ( maskK >= static_cast< proshade_signed > ( yDimInds ) ) { maskK -= yDimInds; }
759  maskL = zIt + static_cast< proshade_signed > ( zDimInds / 2 ); if ( maskL >= static_cast< proshade_signed > ( zDimInds ) ) { maskL -= zDimInds; }
760 
761  //================================ Find the positions
762  maskMapIndex = zIt + static_cast< proshade_signed > ( zDimInds ) * ( yIt + static_cast< proshade_signed > ( yDimInds ) * xIt );
763  densMapIndex = maskL + static_cast< proshade_signed > ( zDimInds ) * ( maskK + static_cast< proshade_signed > ( yDimInds ) * maskH );
764 
765  //================================ Save the values
766  modifCoeffs[densMapIndex][0] = modifCoeffsHKL[maskMapIndex][0];
767  modifCoeffs[densMapIndex][1] = modifCoeffsHKL[maskMapIndex][1];
768  }
769  }
770  }
771 
772  //============================================ Run inverse Fourier on the modified coefficients
773  fftw_execute ( inverseFoourier );
774 
775  //============================================ Delete old weights and allocate memory for the new, re-sampled weights
776  weightsFinal = new proshade_double [origVolume];
777  ProSHADE_internal_misc::checkMemoryAllocation ( weightsFinal, __FILE__, __LINE__, __func__ );
778 
779  //============================================ Copy results into a new, properly sampled weights
780  for ( size_t iter = 0; iter < origVolume; iter++ ) { weightsFinal[iter] = outMap[iter][0]; }
781 
782  //============================================ Release remaining memory
783  fftw_destroy_plan ( planForwardFourier );
784  fftw_destroy_plan ( inverseFoourier );
785  delete[] origCoeffs;
786  delete[] modifCoeffs;
787  delete[] origCoeffsHKL;
788  delete[] modifCoeffsHKL;
789  delete[] inMap;
790  delete[] outMap;
791  }
792  else
793  {
794  weightsFinal = new proshade_double [origVolume];
795  ProSHADE_internal_misc::checkMemoryAllocation ( weightsFinal, __FILE__, __LINE__, __func__ );
796  for ( size_t iter = 0; iter < origVolume; iter++ ) { weightsFinal[iter] = weights[iter]; }
797  }
798 
799  //================================================ Allocate memory for map Fourier transform
800  fftw_complex* inMap = new fftw_complex [origVolume];
801  fftw_complex* outMap = new fftw_complex [origVolume];
802  ProSHADE_internal_misc::checkMemoryAllocation ( inMap, __FILE__, __LINE__, __func__ );
803  ProSHADE_internal_misc::checkMemoryAllocation ( outMap, __FILE__, __LINE__, __func__ );
804  fftw_plan planForwardFourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), inMap, outMap, FFTW_FORWARD, FFTW_ESTIMATE );
805  fftw_plan inverseFoourier = fftw_plan_dft_3d ( static_cast< int > ( xDimInds ), static_cast< int > ( yDimInds ), static_cast< int > ( zDimInds ), outMap, inMap, FFTW_BACKWARD, FFTW_ESTIMATE );
806 
807  //================================================ Set data
808  for ( size_t iter = 0; iter < static_cast< size_t > ( origVolume ); iter++ ) { inMap[iter][0] = map[iter]; inMap[iter][1] = 0.0; }
809 
810  //================================================ Convert map to Fourier space
811  fftw_execute ( planForwardFourier );
812 
813  //================================================ Apply the weights to the map in Fourier space
814  proshade_double normFactor = static_cast<proshade_double> ( origVolume );
815  for ( size_t iter = 0; iter < static_cast< size_t > ( origVolume ); iter++ ) { outMap[iter][0] *= weightsFinal[iter] / normFactor; outMap[iter][1] *= weightsFinal[iter] / normFactor; }
816 
817  //================================================ Convert weighted map from Fourier space
818  fftw_execute ( inverseFoourier );
819 
820  //================================================ Copy results to map
821  for ( size_t iter = 0; iter < static_cast< size_t > ( xDimInds * yDimInds * zDimInds ); iter++ ) { map[iter] = inMap[iter][0]; }
822 
823  //================================================ Release memory
824  delete[] weightsFinal;
825  delete[] inMap;
826  delete[] outMap;
827  fftw_destroy_plan ( planForwardFourier );
828  fftw_destroy_plan ( inverseFoourier );
829 
830  //================================================ Report progress
831  ProSHADE_internal_messages::printProgressMessage ( verbose, 3, "Mask read in and applied successfully." );
832 
833  //================================================ Done
834  return ;
835 
836 }

◆ figureDataType()

ProSHADE_internal_io::InputType ProSHADE_internal_io::figureDataType ( std::string  fName)

Function determining input data type.

This function determines the type of the input structure. The possible outputs are MAP for MRC map files, PDB for mmCIF or PDB formatted data, or UNKNOWN if gemmi fail to read the file as co-ordinates as well as map.

Parameters
[in]fNameThe file name of the file for which the type should be determined.
[out]XProSHADE InputType variable with values UNKNOWN, MAP or PDB depending on the type of the input file.

Definition at line 918 of file ProSHADE_io.cpp.

919 {
920  //================================================ Try readin as PDB
921  if ( isFilePDB ( fName ) )
922  {
923  return ( PDB );
924  }
925 
926  //================================================ If not, try readin as MAP
927  if ( isFileMAP ( fName ) )
928  {
929  return ( MAP );
930  }
931 
932  //================================================ No luck? UNKNOWN it is ...
933  return ( UNKNOWN );
934 
935  //================================================ Done
936 
937 }

◆ isFileMAP()

bool ProSHADE_internal_io::isFileMAP ( std::string  fName)

Function determining if the input data type is MAP.

This function checks if the input file is a MAP file and can be read by the CMAP library.

Parameters
[in]fNameThe file name of the file for which the type should be determined.
[out]XBool value true if the file is a MAP file readable by CMAP and false otherwise.

Definition at line 60 of file ProSHADE_io.cpp.

61 {
62  gemmi::Ccp4<float> map;
63  try
64  {
65  map.read_ccp4 ( gemmi::MaybeGzipped (fName.c_str() ) );
66  }
67  catch ( std::runtime_error& e )
68  {
69  //============================================ Supress MSVC C4101 Unferenced variable warning / clang and gcc -Wunused-exception-parameter
70  (void)e;
71 
72  //============================================ Failed to read the map
73  return ( false );
74  }
75 
76  //================================================ Done
77  return ( true );
78 
79 }

◆ isFilePDB()

bool ProSHADE_internal_io::isFilePDB ( std::string  fName)

Function determining if the input data type is PDB.

This function checks if the input file is a PDB file and can be read by the gemmi library.

Parameters
[in]fNameThe file name of the file for which the type should be determined.
[out]XBool value true if the file is a PDB file readable by gemmi and false otherwise.

Definition at line 32 of file ProSHADE_io.cpp.

33 {
34  //================================================ Try reading the file using Gemmi
35  try
36  {
37  gemmi::Structure structure = gemmi::read_structure ( gemmi::MaybeGzipped ( fName ) );
38  }
39  catch ( std::runtime_error& e )
40  {
41  //============================================ Supress MSVC C4101 Unferenced variable warning / clang and gcc -Wunused-exception-parameter
42  (void)e;
43 
44  //============================================ Read failed. Done
45  return ( false );
46  }
47 
48  //================================================ Read successfull. Done
49  return ( true );
50 
51 }

◆ readInMapData() [1/2]

void ProSHADE_internal_io::readInMapData ( gemmi::Ccp4< float > *  gemmiMap,
proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_unsign  xAxOrder,
proshade_unsign  yAxOrder,
proshade_unsign  zAxOrder 
)

This function converts the gemmi Ccp4 object data to ProSHADE internal map representation.

This function firstly allocates the required memory for the ProSHADE internal map representation variable according to the grid size. Then, it iterates over the axes in such a way, so that the resulting ProSHADE variable would have XYZ axis order independently on the axis order of the Ccp4 gemmi object. This should not be necessary as the gemmi setup function should have been called by now, but one never knows.

Parameters
[in]gemmiMapPointer to a gemmi Ccp4 object containing the read in MAP file information.
[in]mapPointer reference to a variable to save the map data.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]xAxOrderThe order of the x-axis.
[in]yAxOrderThe order of the y-axis.
[in]zAxOrderThe order of the z-axis.

Definition at line 178 of file ProSHADE_io.cpp.

179 {
180  //================================================ Allocate internal variables
181  proshade_unsign *axOrdArr = new proshade_unsign[3];
182  proshade_unsign *axDimArr = new proshade_unsign[3];
183  proshade_unsign arrPos = 0;
184 
185  //================================================ Check memory allocation and fill in values
186  ProSHADE_internal_misc::checkMemoryAllocation ( axOrdArr, __FILE__, __LINE__, __func__ );
187  ProSHADE_internal_misc::checkMemoryAllocation ( axDimArr, __FILE__, __LINE__, __func__ );
188  axDimArr[0] = xDimInds;
189  axDimArr[1] = yDimInds;
190  axDimArr[2] = zDimInds;
191 
192  //================================================ Allocate the ProSHADE internal map variable memory
193  map = new proshade_double [xDimInds * yDimInds * zDimInds];
194  ProSHADE_internal_misc::checkMemoryAllocation ( map, __FILE__, __LINE__, __func__ );
195 
196  //================================================ Copy read in data to internal map variable
197  for ( axOrdArr[0] = 0; axOrdArr[0] < axDimArr[xAxOrder-1]; axOrdArr[0]++ )
198  {
199  for ( axOrdArr[1] = 0; axOrdArr[1] < axDimArr[yAxOrder-1]; axOrdArr[1]++ )
200  {
201  for ( axOrdArr[2] = 0; axOrdArr[2] < axDimArr[zAxOrder-1]; axOrdArr[2]++ )
202  {
203  arrPos = axOrdArr[2] + axDimArr[zAxOrder-1] * ( axOrdArr[1] + axDimArr[yAxOrder-1] * axOrdArr[0] );
204  map[arrPos] = static_cast< proshade_double > ( gemmiMap->grid.get_value_q( static_cast< int > ( axOrdArr[xAxOrder-1] ),
205  static_cast< int > ( axOrdArr[yAxOrder-1] ),
206  static_cast< int > ( axOrdArr[zAxOrder-1] ) ) );
207  }
208  }
209  }
210 
211  //================================================ Release internal variables memory
212  delete[] axDimArr;
213  delete[] axOrdArr;
214 
215  //================================================ Done
216  return ;
217 
218 }

◆ readInMapData() [2/2]

void ProSHADE_internal_io::readInMapData ( gemmi::Ccp4< int8_t > *  gemmiMap,
proshade_double *&  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_unsign  xAxOrder,
proshade_unsign  yAxOrder,
proshade_unsign  zAxOrder 
)

This function converts the gemmi Ccp4 object data to ProSHADE mask representation.

This function firstly allocates the required memory for the ProSHADE mask representation variable according to the grid size. Then, it iterates over the axes in such a way, so that the resulting ProSHADE variable would have XYZ axis order independently on the axis order of the Ccp4 gemmi object. This should not be necessary as the gemmi setup function should have been called by now, but one never knows.

Parameters
[in]gemmiMapPointer to a gemmi Ccp4 object containing the read in mask file information.
[in]mapPointer reference to a variable to save the map data.
[in]xDimIndsThe size of x dimension in indices.
[in]yDimIndsThe size of y dimension in indices.
[in]zDimIndsThe size of z dimension in indices.
[in]xAxOrderThe order of the x-axis.
[in]yAxOrderThe order of the y-axis.
[in]zAxOrderThe order of the z-axis.

Definition at line 234 of file ProSHADE_io.cpp.

235 {
236  //================================================ Allocate internal variables
237  proshade_unsign *axOrdArr = new proshade_unsign[3];
238  proshade_unsign *axDimArr = new proshade_unsign[3];
239  proshade_unsign arrPos = 0;
240 
241  //================================================ Check memory allocation and fill in values
242  ProSHADE_internal_misc::checkMemoryAllocation ( axOrdArr, __FILE__, __LINE__, __func__ );
243  ProSHADE_internal_misc::checkMemoryAllocation ( axDimArr, __FILE__, __LINE__, __func__ );
244  axDimArr[0] = xDimInds;
245  axDimArr[1] = yDimInds;
246  axDimArr[2] = zDimInds;
247 
248  //================================================ Allocate the ProSHADE internal map variable memory
249  map = new proshade_double [xDimInds * yDimInds * zDimInds];
250  ProSHADE_internal_misc::checkMemoryAllocation ( map, __FILE__, __LINE__, __func__ );
251 
252  //================================================ Copy read in data to internal map variable
253  for ( axOrdArr[0] = 0; axOrdArr[0] < axDimArr[xAxOrder-1]; axOrdArr[0]++ )
254  {
255  for ( axOrdArr[1] = 0; axOrdArr[1] < axDimArr[yAxOrder-1]; axOrdArr[1]++ )
256  {
257  for ( axOrdArr[2] = 0; axOrdArr[2] < axDimArr[zAxOrder-1]; axOrdArr[2]++ )
258  {
259  arrPos = axOrdArr[2] + axDimArr[zAxOrder-1] * ( axOrdArr[1] + axDimArr[yAxOrder-1] * axOrdArr[0] );
260  map[arrPos] = static_cast< proshade_double > ( gemmiMap->grid.get_value_q( static_cast< int > ( axOrdArr[xAxOrder-1] ),
261  static_cast< int > ( axOrdArr[yAxOrder-1] ),
262  static_cast< int > ( axOrdArr[zAxOrder-1] ) ) );
263  }
264  }
265  }
266 
267  //================================================ Release internal variables memory
268  delete[] axDimArr;
269  delete[] axOrdArr;
270 
271  //================================================ Done
272  return ;
273 
274 }

◆ readInMapHeader()

void ProSHADE_internal_io::readInMapHeader ( gemmi::Ccp4< float > *  map,
proshade_unsign *  xDimInds,
proshade_unsign *  yDimInds,
proshade_unsign *  zDimInds,
proshade_single *  xDim,
proshade_single *  yDim,
proshade_single *  zDim,
proshade_single *  aAng,
proshade_single *  bAng,
proshade_single *  cAng,
proshade_signed *  xFrom,
proshade_signed *  yFrom,
proshade_signed *  zFrom,
proshade_signed *  xAxOrigin,
proshade_signed *  yAxOrigin,
proshade_signed *  zAxOrigin,
proshade_unsign *  xAxOrder,
proshade_unsign *  yAxOrder,
proshade_unsign *  zAxOrder,
proshade_unsign *  xGridInds,
proshade_unsign *  yGridInds,
proshade_unsign *  zGridInds 
)

This function parses the CCP4 MAP file header as read in by gemmi.

This function uses the gemmi Ccp4 object, which contains all the information read in from a MAP file (including the header), to parse out the ProSHADE required information from the header and saving it to the supplied variables.

Parameters
[in]mapA gemmi Ccp4 objecct containing all the data read in from a MAP file.
[in]xDimIndsAddress to a variable to save the x-axis size in indices.
[in]yDimIndsAddress to a variable to save the y-axis size in indices.
[in]zDimIndsAddress to a variable to save the z-axis size in indices.
[in]xDimAddress to a variable to save the x dimension size in angstroms.
[in]yDimAddress to a variable to save the y dimension size in angstroms.
[in]zDimAddress to a variable to save the z dimension size in angstroms.
[in]aAngAddress to a variable to save the a angle in degrees.
[in]bAngAddress to a variable to save the b angle in degrees.
[in]cAngAddress to a variable to save the c angle in degrees.
[in]xFromAddress to a variable to save the starting index along the x-axis.
[in]yFromAddress to a variable to save the starting index along the y-axis.
[in]zFromAddress to a variable to save the starting index along the z-axis.
[in]xAxOriginAddress to a variable to save the map origin positon along the x-axis.
[in]yAxOriginAddress to a variable to save the map origin positon along the y-axis.
[in]zAxOriginAddress to a variable to save the map origin positon along the z-axis.
[in]xAxOrderAddress to a variable to save the order of x axis.
[in]yAxOrderAddress to a variable to save the order of y axis.
[in]zAxOrderAddress to a variable to save the order of z axis.
[in]xGridIndsAddress to a variable to save the grid indices count along the x-axis.
[in]yGridIndsAddress to a variable to save the grid indices count along the y-axis.
[in]zGridIndsAddress to a variable to save the grid indices count along the z-axis.

Definition at line 109 of file ProSHADE_io.cpp.

110 {
111  //================================================ Read in the map file header
112  *xDimInds = static_cast<proshade_unsign> ( map->header_i32 ( 1 ) );
113  *yDimInds = static_cast<proshade_unsign> ( map->header_i32 ( 2 ) );
114  *zDimInds = static_cast<proshade_unsign> ( map->header_i32 ( 3 ) );
115 
116  *xFrom = static_cast<proshade_signed> ( map->header_i32 ( 5 ) );
117  *yFrom = static_cast<proshade_signed> ( map->header_i32 ( 6 ) );
118  *zFrom = static_cast<proshade_signed> ( map->header_i32 ( 7 ) );
119 
120  *xDim = static_cast<proshade_single> ( map->header_float ( 11 ) );
121  *yDim = static_cast<proshade_single> ( map->header_float ( 12 ) );
122  *zDim = static_cast<proshade_single> ( map->header_float ( 13 ) );
123 
124  *aAng = static_cast<proshade_single> ( map->header_float ( 14 ) );
125  *bAng = static_cast<proshade_single> ( map->header_float ( 15 ) );
126  *cAng = static_cast<proshade_single> ( map->header_float ( 16 ) );
127 
128  *xAxOrigin = static_cast<proshade_signed> ( map->header_i32 ( 50 ) ) + (*xFrom);
129  *yAxOrigin = static_cast<proshade_signed> ( map->header_i32 ( 51 ) ) + (*yFrom);
130  *zAxOrigin = static_cast<proshade_signed> ( map->header_i32 ( 52 ) ) + (*zFrom);
131 
132  *xAxOrder = static_cast<proshade_unsign> ( map->header_i32 ( 17 ) );
133  *yAxOrder = static_cast<proshade_unsign> ( map->header_i32 ( 18 ) );
134  *zAxOrder = static_cast<proshade_unsign> ( map->header_i32 ( 19 ) );
135 
136  *xGridInds = static_cast<proshade_unsign> ( map->header_i32 ( 8 ) );
137  *yGridInds = static_cast<proshade_unsign> ( map->header_i32 ( 9 ) );
138  *zGridInds = static_cast<proshade_unsign> ( map->header_i32 ( 10 ) );
139 
140  //================================================ Deal with sampling being different from cell size
141  if ( *xGridInds != *xDimInds )
142  {
143  *xDim = *xDim * ( static_cast<proshade_single> ( *xDimInds ) / static_cast<proshade_single> ( *xGridInds ) );
144  *xGridInds = *xDimInds;
145  }
146 
147  if ( *yGridInds != *yDimInds )
148  {
149  *yDim = *yDim * ( static_cast<proshade_single> ( *yDimInds ) / static_cast<proshade_single> ( *yGridInds ) );
150  *yGridInds = *yDimInds;
151  }
152 
153  if ( *zGridInds != *zDimInds )
154  {
155  *zDim = *zDim * ( static_cast<proshade_single> ( *zDimInds ) / static_cast<proshade_single> ( *zGridInds ) );
156  *zGridInds = *zDimInds;
157  }
158 
159  //================================================ Done
160  return ;
161 
162 }

◆ writeOutMapHeader()

void ProSHADE_internal_io::writeOutMapHeader ( gemmi::Ccp4< float > *  map,
proshade_unsign  xDimInds,
proshade_unsign  yDimInds,
proshade_unsign  zDimInds,
proshade_single  xDim,
proshade_single  yDim,
proshade_single  zDim,
proshade_single  aAng,
proshade_single  bAng,
proshade_single  cAng,
proshade_signed  xFrom,
proshade_signed  yFrom,
proshade_signed  zFrom,
proshade_signed  xAxOrigin,
proshade_signed  yAxOrigin,
proshade_signed  zAxOrigin,
proshade_unsign  xAxOrder,
proshade_unsign  yAxOrder,
proshade_unsign  zAxOrder,
proshade_unsign  xGridInds,
proshade_unsign  yGridInds,
proshade_unsign  zGridInds,
std::string  title,
int  mode 
)

This function parses the CCP4 MAP file header as read in by gemmi.

This function uses the gemmi Ccp4 object, which contains all the information read in from a MAP file (including the header), to parse out the ProSHADE required information from the header and saving it to the supplied variables.

Parameters
[in]mapA gemmi Ccp4 objecct containing all the data read in from a MAP file.
[in]xDimIndsVariable holding the x-axis size in indices.
[in]yDimIndsVariable holding the y-axis size in indices.
[in]zDimIndsVariable holding the z-axis size in indices.
[in]xDimVariable holding the x dimension size in angstroms.
[in]yDimVariable holding the y dimension size in angstroms.
[in]zDimVariable holding the z dimension size in angstroms.
[in]aAngVariable holding the a angle in degrees.
[in]bAngVariable holding the b angle in degrees.
[in]cAngVariable holding the c angle in degrees.
[in]xFromVariable holding the starting index along the x-axis.
[in]yFromVariable holding the starting index along the y-axis.
[in]zFromVariable holding the starting index along the z-axis.
[in]xAxOriginVariable holding the map origin positon along the x-axis.
[in]yAxOriginVariable holding the map origin positon along the y-axis.
[in]zAxOriginVariable holding the map origin positon along the z-axis.
[in]xAxOrderVariable holding the order of x axis.
[in]yAxOrderVariable holding the order of y axis.
[in]zAxOrderVariable holding the order of z axis.
[in]xGridIndsVariable holding the grid indices count along the x-axis.
[in]yGridIndsVariable holding the grid indices count along the y-axis.
[in]zGridIndsVariable holding the grid indices count along the z-axis.
[in]titleThe title to be written into the MAP file.
[in]modeThe variable type of the data, please leave two (float) unless you require any specific other mode.

Definition at line 868 of file ProSHADE_io.cpp.

869 {
870  //================================================ Fill in the map file header
871  map->set_header_i32 ( 1 , static_cast<int32_t> ( xDimInds ) ); // Number of columns in 3D data array (fast axis)
872  map->set_header_i32 ( 2 , static_cast<int32_t> ( yDimInds ) ); // Number of columns in 3D data array (medium axis)
873  map->set_header_i32 ( 3 , static_cast<int32_t> ( zDimInds ) ); // Number of columns in 3D data array (slow axis)
874  map->set_header_i32 ( 4 , static_cast<int32_t> ( mode ) ); // Map mode
875  map->set_header_i32 ( 5 , static_cast<int32_t> ( xFrom ) ); // Starting index (fast axis)
876  map->set_header_i32 ( 6 , static_cast<int32_t> ( yFrom ) ); // Starting index (medium axis)
877  map->set_header_i32 ( 7 , static_cast<int32_t> ( zFrom ) ); // Starting index (slow axis)
878  map->set_header_i32 ( 8 , static_cast<int32_t> ( xGridInds ) ); // Grid sampling (fast axis)
879  map->set_header_i32 ( 9 , static_cast<int32_t> ( yGridInds ) ); // Grid sampling (medium axis)
880  map->set_header_i32 ( 10, static_cast<int32_t> ( zGridInds ) ); // Grid sampling (slow axis)
881  map->set_header_float ( 11, static_cast<float> ( xDim ) ); // Grid dimension in Angstrom (fast axis)
882  map->set_header_float ( 12, static_cast<float> ( yDim ) ); // Grid dimension in Angstrom (medium axis)
883  map->set_header_float ( 13, static_cast<float> ( zDim ) ); // Grid dimension in Angstrom (slow axis)
884  map->set_header_float ( 14, static_cast<float> ( aAng ) ); // Alpha angle in degrees
885  map->set_header_float ( 15, static_cast<float> ( bAng ) ); // Beta angle in degrees
886  map->set_header_float ( 16, static_cast<float> ( cAng ) ); // Gamma angle in degrees
887  map->set_header_i32 ( 17, static_cast<int32_t> ( xAxOrder ) ); // MAPC
888  map->set_header_i32 ( 18, static_cast<int32_t> ( yAxOrder ) ); // MAPR
889  map->set_header_i32 ( 19, static_cast<int32_t> ( zAxOrder ) ); // MAPS
890  if ( map->grid.spacegroup ) { map->set_header_i32 ( 23, static_cast<int32_t> ( map->grid.spacegroup->ccp4 ) ); } // Space group
891  else { map->set_header_i32 ( 23, static_cast<int32_t> ( 1 ) ); }
892  map->set_header_i32 ( 24, static_cast<int32_t> ( map->grid.spacegroup->operations().order() * 80 ) ); // NSYMBT - size of extended header (which follows main header) in bytes
893  map->set_header_str ( 27, "CCP4" ); // Code for the type of extended header
894  map->set_header_i32 ( 28, static_cast<int32_t> ( 20140 ) ); // Version
895  map->set_header_i32 ( 50, static_cast<int32_t> ( xAxOrigin ) ); // Origin of the map (fast axis)
896  map->set_header_i32 ( 51, static_cast<int32_t> ( yAxOrigin ) ); // Origin of the map (medium axis)
897  map->set_header_i32 ( 52, static_cast<int32_t> ( zAxOrigin ) ); // Origin of the map (slow axis)
898  map->set_header_str ( 53, "MAP" ); // File format
899  if ( gemmi::is_little_endian() ) { map->set_header_i32 ( 54, static_cast<int32_t> ( 0x00004144 ) ); } // Machine stamp encoding byte ordering of data
900  else { map->set_header_i32 ( 54, static_cast<int32_t> ( 0x11110000 ) ); }
901  map->set_header_i32 ( 56, static_cast<int32_t> ( 1 ) ); // Number of labels used
902  std::memset ( reinterpret_cast<void*> ( &(map->ccp4_header.at( 56 )) ), ' ', static_cast< size_t > ( 800 + map->grid.spacegroup->operations().order() * 80 ) ); // 56 is used because the vector is indexed from 0
903  map->set_header_str ( 57, title ); // Title
904 
905  //================================================ Done
906  return ;
907 
908 }

◆ writeRotationTranslationJSON()

void ProSHADE_internal_io::writeRotationTranslationJSON ( proshade_double  trsX1,
proshade_double  trsY1,
proshade_double  trsZ1,
proshade_double  eulA,
proshade_double  eulB,
proshade_double  eulG,
proshade_double  trsX2,
proshade_double  trsY2,
proshade_double  trsZ2,
std::string  fileName 
)

Function for writing out the optimal rotation and translation into a JSON file.

This function takes the initial translation (assuming that the centre of rotation is not at the centre of indices), the rotation (in terms of the Euler angles) and the translation detected by the overlay task and proceeds to write a JSON file containing all of these information in the order in which they should be applied with the exception that depending on around which point the rotation should be done, either the translation to origin or to the map centre should be applied.

This function assumes that the second translation includes the reverse of the first translation already.

Parameters
[in]trsX1The translation required to get the rotation centre to origin along the x-axis.
[in]trsY1The translation required to get the rotation centre to origin along the y-axis.
[in]trsZ1The translation required to get the rotation centre to origin along the z-axis.
[in]eulAThe optimal rotation Euler angle alpha.
[in]eulBThe optimal rotation Euler angle beta.
[in]eulGThe optimal rotation Euler angle gamma.
[in]trsX2The optimal translation along the x-axis + reverse of the trsX1.
[in]trsY2The optimal translation along the y-axis + reverse of the trsY1.
[in]trsZ2The optimal translation along the z-axis + reverse of the trsZ1.
[in]fileNameThe file name of the file for which the information should be written into.

Definition at line 960 of file ProSHADE_io.cpp.

961 {
962  //================================================ Open file for writing
963  std::ofstream jsonFile;
964  jsonFile.open ( fileName );
965 
966  //================================================ Check file opening success
967  if ( !jsonFile.is_open( ) )
968  {
969  throw ProSHADE_exception ( "Failed to open JSON output file.", "E000056", __FILE__, __LINE__, __func__, "Failed to open json file to which the rotation and\n : translation would be written into. Most likely cause is\n : lack of rights to write in the current folder." );
970  }
971 
972  //================================================ Get rotation matrix from Euler angles
973  proshade_double* rotMat = new proshade_double[9];
974  ProSHADE_internal_misc::checkMemoryAllocation ( rotMat, __FILE__, __LINE__, __func__ );
976 
977  //================================================ Write the info
978  jsonFile << "{\n";
979  jsonFile << " \"translationToOrigin\" : [ " << trsX1 << ", " << trsY1 << ", " << trsZ1 << " ], \n";
980 
981  jsonFile << " \"rotationMatrix:\" : [ " << rotMat[0] << ", " << rotMat[1] << ", " << rotMat[2] << ", \n";
982  jsonFile << " " << rotMat[3] << ", " << rotMat[4] << ", " << rotMat[5] << ", \n";
983  jsonFile << " " << rotMat[6] << ", " << rotMat[7] << ", " << rotMat[8] << "], \n";
984 
985  jsonFile << " \"translationFromRotCenToOverlay\" : [ " << trsX2 << ", " << trsY2 << ", " << trsZ2 << " ] \n";
986  jsonFile << "}\n";
987 
988  //================================================ Close file
989  jsonFile.close ( );
990 
991  //================================================ Release memory
992  delete[] rotMat;
993 
994  //================================================ Done
995  return ;
996 
997 }
ProSHADE_internal_io::isFilePDB
bool isFilePDB(std::string fName)
Function determining if the input data type is PDB.
Definition: ProSHADE_io.cpp:32
ProSHADE_exception
This class is the representation of ProSHADE exception.
Definition: ProSHADE_exceptions.hpp:37
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_internal_io::isFileMAP
bool isFileMAP(std::string fName)
Function determining if the input data type is MAP.
Definition: ProSHADE_io.cpp:60
ProSHADE_internal_io::readInMapData
void readInMapData(gemmi::Ccp4< float > *gemmiMap, proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_unsign xAxOrder, proshade_unsign yAxOrder, proshade_unsign zAxOrder)
This function converts the gemmi Ccp4 object data to ProSHADE internal map representation.
Definition: ProSHADE_io.cpp:178
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_io::applyWeightsFromArray
void applyWeightsFromArray(proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&weights, proshade_unsign xDimIndsWgh, proshade_unsign yDimIndsWgh, proshade_unsign zDimIndsWgh, proshade_signed verbose)
This function applies the weights to the map.
Definition: ProSHADE_io.cpp:633
ProSHADE_internal_io::applyMaskFromArray
void applyMaskFromArray(proshade_double *&map, proshade_unsign xDimInds, proshade_unsign yDimInds, proshade_unsign zDimInds, proshade_double *&mask, proshade_unsign xDimIndsMsk, proshade_unsign yDimIndsMsk, proshade_unsign zDimIndsMsk, proshade_signed verbose)
This function applies the mask to the map.
Definition: ProSHADE_io.cpp:365
ProSHADE_internal_messages::printProgressMessage
void printProgressMessage(proshade_signed verbose, proshade_signed messageLevel, std::string message)
General stdout message printing.
Definition: ProSHADE_messages.cpp:70
ProSHADE_internal_io::readInMapHeader
void readInMapHeader(gemmi::Ccp4< float > *map, proshade_unsign *xDimInds, proshade_unsign *yDimInds, proshade_unsign *zDimInds, proshade_single *xDim, proshade_single *yDim, proshade_single *zDim, proshade_single *aAng, proshade_single *bAng, proshade_single *cAng, proshade_signed *xFrom, proshade_signed *yFrom, proshade_signed *zFrom, proshade_signed *xAxOrigin, proshade_signed *yAxOrigin, proshade_signed *zAxOrigin, proshade_unsign *xAxOrder, proshade_unsign *yAxOrder, proshade_unsign *zAxOrder, proshade_unsign *xGridInds, proshade_unsign *yGridInds, proshade_unsign *zGridInds)
This function parses the CCP4 MAP file header as read in by gemmi.
Definition: ProSHADE_io.cpp:109