pdb.h

// -*-c++-*-
#ifndef STORMM_PDB_H
#define STORMM_PDB_H
 
#include <string>
#include <vector>
#include "copyright.h"
#include "Accelerator/hybrid.h"
#include "Accelerator/gpu_enumerators.h"
#include "DataTypes/stormm_vector_types.h"
#include "FileManagement/file_enumerators.h"
#include "Parsing/textfile.h"
#include "Synthesis/condensate.h"
#include "Synthesis/phasespace_synthesis.h"
#include "Topology/atomgraph.h"
#include "Topology/atomgraph_enumerators.h"
#include "Trajectory/coordinateframe.h"
#include "Trajectory/coordinate_series.h"
#include "Trajectory/phasespace.h"
#include "molecule_format_enumerators.h"
 
namespace stormm {
namespace structure {
 
using card::default_hpc_format;
using card::HybridFormat;
using diskutil::PrintSituation;
using parse::TextFile;
using synthesis::Condensate;
using synthesis::PhaseSpaceSynthesis;
using topology::AtomGraph;
using topology::UnitCellType;
using trajectory::CoordinateFrame;
using trajectory::CoordinateSeries;
using trajectory::PhaseSpace;
 
class Pdb {
public:

  Pdb(int atom_count_in = 0, int model_count_in = 1);
 
  Pdb(const std::string &file_name_in, int model_index = -1,
      const std::vector<char> &alternate_prefs = { ' ', 'A' });
 
  Pdb(const std::string &file_name_in, const std::vector<int> &model_indices,
      const std::vector<char> &alternate_prefs = { ' ', 'A' });
 
  Pdb(int atom_count_in, int model_count_in, const std::vector<double> &x_coordinates_in,
      const std::vector<double> &y_coordinates_in, const std::vector<double> &z_coordinates_in,
      const std::vector<char4> &atom_names_in, const std::vector<char4> &residue_names_in,
      const std::vector<char> &chain_names_in, const std::vector<double> &occupancies_in,
      const std::vector<double> &b_factors_in, const std::vector<char2> &elements_in = {},
      const std::vector<int> &atom_numbers_in = {},
      const std::vector<int> &residue_numbers_in = {},
      const std::vector<int> &anisotropy_matrices_in = {},
      const std::vector<PdbAtomKind> &atom_classes_in = {},
      const std::vector<std::string> &header_in = {},
      const std::vector<std::string> &remark_in = {},
      const std::vector<std::string> &extras_in = {});
 
  Pdb(const AtomGraph *ag, const CoordinateFrame *cf);
 
  Pdb(const AtomGraph &ag, const CoordinateFrame &cf);
 
  Pdb(const AtomGraph *ag, const PhaseSpace *ps);  
 
  Pdb(const AtomGraph &ag, const PhaseSpace &ps);
 
  template <typename T>
  Pdb(const AtomGraph *ag, const CoordinateSeries<T> *cf, const std::vector<int> &frames = {});
 
  template <typename T>
  Pdb(const AtomGraph &ag, const CoordinateSeries<T> &cf, const std::vector<int> &frames = {});
 
  template <typename T> Pdb(const AtomGraph *ag, const CoordinateSeries<T> *cs, int frame_number);
 
  template <typename T> Pdb(const AtomGraph &ag, const CoordinateSeries<T> &cs, int frame_number);

  int getAtomCount() const;

  int getModelCount() const;

  bool getFormalChargePrinting() const;
  
  CoordinateFrame exportCoordinateFrame(int model_index,
                                        HybridFormat layout = default_hpc_format) const;

  PhaseSpace exportPhaseSpace(int model_index, HybridFormat layout = default_hpc_format) const;

  template <typename T>
  CoordinateSeries<T> getCoordinateSeries(const std::vector<int> model_sequence = {}) const;

  void writeToFile(const std::string &file_name, PrintSituation expectation,
                   const std::vector<int> &model_indices = {}) const;
  
  void loadTopologicalData(const AtomGraph *ag);
  void loadTopologicalData(const AtomGraph &ag);

  template <typename T>
  void loadCoordinates(const T* xcrd_in, const T* ycrd_in, const T* zcrd_in, int model_index = 0,
                       double gpos_scale = 1.0);
 
  void loadCoordinates(const llint* xcrd_in, const llint* ycrd_in, const llint* zcrd_in,
                       const int* xcrd_ovrf, const int* ycrd_ovrf, const int* zcrd_ovrf,
                       double gpos_scale = 1.0, int model_index = 0);
  
  void loadCoordinates(const CoordinateFrame *cf, int model_index = 0);
 
  void loadCoordinates(const CoordinateFrame &cf, int model_index = 0);
 
  void loadCoordinates(const PhaseSpace *ps, int model_index = 0);
 
  void loadCoordinates(const PhaseSpace &ps, int model_index = 0);
 
  template <typename T>
  void loadCoordinates(const CoordinateSeries<T> *cs, const std::vector<int2> &frame_mapping = {});
 
  template <typename T>
  void loadCoordinates(const CoordinateSeries<T> &cs, const std::vector<int2> &frame_mapping = {});
 
  void loadCoordinates(const PhaseSpaceSynthesis *poly_ps, int system_index, int model_index = 0);
 
  void loadCoordinates(const PhaseSpaceSynthesis &poly_ps, int system_index, int model_index = 0);
 
  void loadCoordinates(const Condensate *cdns, int system_index, int model_index = 0);
 
  void loadCoordinates(const Condensate &cdns, int system_index, int model_index = 0);

  void setFormalChargePrinting(bool print_formal_charges_in = true);
  
private:
  std::string file_name;                  
  int atom_count;                         
  int model_count;                        
  int padded_atom_count;                  
  int ter_card_count;                     
  UnitCellType unit_cell;                 
  std::vector<PdbAtomKind> atom_classes;  
  std::vector<int> atom_numbers;          
  std::vector<char4> atom_names;          
  std::vector<char4> residue_names;       
  std::vector<char> chain_names;          
  std::vector<int> residue_numbers;       
  std::vector<double> x_coordinates;      
  std::vector<double> y_coordinates;      
  std::vector<double> z_coordinates;      
  std::vector<double> box_dimensions;     
  std::vector<double> occupancies;        
  std::vector<double> b_factors;          
  std::vector<char2> atom_symbols;        
  std::vector<int> atom_formal_charges;   
  std::vector<int> ter_card_locations;    
 
  // The following are anisotropic temperature factors for each atom, providing the six unique
  // elements of each atom's matrix.
  std::vector<int> anisotropy_xx;  
  std::vector<int> anisotropy_xy;  
  std::vector<int> anisotropy_xz;  
  std::vector<int> anisotropy_yy;  
  std::vector<int> anisotropy_yz;  
  std::vector<int> anisotropy_zz;  
 
  // The following are output settings which will be called upon when printing the Pdb object's
  // contents to a file.
  bool print_formal_charges;
  
  void allocate();
  
  void computeTransforms(const int model_index, std::vector<double> *umat,
                         std::vector<double> *invu) const;

  int countAtomRecords(const TextFile &tf, int line_start = 0, int line_end = -1) const;

  void readAtomRecords(const TextFile &tf, int line_start, int line_end, int atom_start,
                       const std::vector<char> &alternate_prefs = { ' ', 'A' });

  void seekUnitCellDimensions(const TextFile &tf, int line_start, int line_end, int box_start);
  
  void validateModelIndex(int model_index, const char* caller = nullptr) const;

  void writeModel(std::ofstream *foutp, const int model_index) const;
};
 
} // namespace structure
} // namespace stormm
 
#include "pdb.tpp"
 
#endif