snapshotgadget.h 12 KB
Newer Older
jcl's avatar
jcl committed
1
// ============================================================================
2
// Copyright Jean-Charles LAMBERT - 2008-2016
3 4 5 6 7 8 9 10
//           Centre de donneeS Astrophysiques de Marseille (CeSAM)
// e-mail:   Jean-Charles.Lambert@lam.fr
// address:  Aix Marseille Universite, CNRS, LAM
//           Laboratoire d'Astrophysique de Marseille
//           Pole de l'Etoile, site de Chateau-Gombert
//           38, rue Frederic Joliot-Curie
//           13388 Marseille cedex 13 France
//           CNRS UMR 7326
jcl's avatar
jcl committed
11 12 13
// ============================================================================

/**
jclamber's avatar
jclamber committed
14
	@author Jean-Charles Lambert <Jean-Charles.Lambert@lam.fr>
jcl's avatar
jcl committed
15 16 17 18 19 20 21 22 23
 */

#ifndef SNAPSHOTGADGET_H
#define SNAPSHOTGADGET_H

#include <string>
#include <assert.h>
#include <fstream>
#include <map>
24 25
#include <iostream>
#include <sstream>
jcl's avatar
jcl committed
26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43
#include "snapshotinterface.h"

namespace uns {

typedef struct io_header_1
{
  int      npart[6];
  double   mass[6];
  double   time;
  double   redshift;
  int      flag_sfr;
  int      flag_feedback;
  int      npartTotal[6];
  int      flag_cooling;
  int      num_files;
  double   BoxSize;
  double   Omega0;
  double   OmegaLambda;
44
  double   HubbleParam;
jcl's avatar
jcl committed
45 46 47
  char     fill[256- 6*4- 6*8- 2*8- 2*4- 6*4- 2*4 - 4*8];  /* fills to 256 Bytes */
} t_io_header_1;

48 49 50 51 52 53 54 55 56 57
template <class T> class CG2_T_Header {
public:
  T   time;
  T   redshift;
  T   BoxSize;
  T   Omega0;
  T   OmegaLambda;
  T   HubbleParam;
};

jcl's avatar
jcl committed
58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83
typedef struct s_data_type_header {
  int size_bytes;
  int items;
} t_data_type_header;

enum DATA_TYPE {
  INT   =4,
  FLOAT =4,
  DOUBLE=8,
  CHAR  =1
};

const int NB_DATA_HEADER = 14;
const t_data_type_header dth[NB_DATA_HEADER] = {
  {INT      ,6}, // int      npart[6];
  {DOUBLE   ,6}, // double   mass[6];
  {DOUBLE   ,1}, // double   time;
  {DOUBLE   ,1}, // double   redshift;
  {INT      ,1}, // int      flag_sfr;
  {INT      ,1}, // int      flag_feedback;
  {INT      ,6}, // int      npartTotal[6];
  {INT      ,1}, // int      flag_cooling;
  {INT      ,1}, // int      num_files;
  {DOUBLE   ,1}, // double   BoxSize;
  {DOUBLE   ,1}, // double   Omega0;
  {DOUBLE   ,1}, // double   OmegaLambda;
84
  {DOUBLE   ,1}, // double   HubbleParam;
jcl's avatar
jcl committed
85 86 87 88 89 90 91 92 93
  {CHAR     ,96} // char     fill[256- 6*4- 6*8- 2*8- 2*4- 6*4- 2*4 - 4*8];
};


  // ---------------------------------------------------
  // CSnapshotGadgetIn
  // READING class
  // ---------------------------------------------------

94
  template <class T> class CSnapshotGadgetIn : public CSnapshotInterfaceIn<T> {
95

jcl's avatar
jcl committed
96 97 98 99 100 101 102 103 104
  public:
    // READING constrcuctor
    CSnapshotGadgetIn(const std::string, const std::string, const std::string, const bool verb=false);
    ~CSnapshotGadgetIn();
    int nextFrame(uns::UserSelection &);
//    int nextFrame();
    ComponentRangeVector * getSnapshotRange();
    int getNbody() { return getNtotal();}
    // virtual function implemented
105 106
    bool getData(const std::string,int *n,T **);
    bool getData(const std::string,       T * );
jcl's avatar
jcl committed
107 108
    bool getData(const std::string,int *n,int   **);
    bool getData(const std::string,       int   * );
109
    bool getData(const std::string, const std::string ,int *,T **);
jcl's avatar
jcl committed
110 111 112 113 114
    bool getData(const std::string, const std::string ,int *,int   **);
    bool isNewFrame() { return first_loc;}
   int close();

  private:
115

116
   std::map<std::string, std::vector<T> > data_vector; // map container to track vector storage
117 118


jcl's avatar
jcl committed
119 120 121 122 123 124 125 126
   bool first_loc;
   int open(const std::string);

   int read (uns::UserSelection &);
   int getVersion() const { return version;}
   const uns::ComponentRangeVector getCRV() const { return crv;}
   int   getNtotal() const { return npartTotal;}
     std::string filename,file0;
127 128
  std::ifstream in; // input descriptor
  std::ios::pos_type in_start_block;
jcl's avatar
jcl committed
129 130 131

  int multiplefiles;
  bool lonely_file;
132

jcl's avatar
jcl committed
133
  //data
LAMBERT Jean-charles's avatar
LAMBERT Jean-charles committed
134 135
  T * mass, * pos, * vel, * acc, *pot, * rho,
    * hsml, * age, * metal, * intenerg, * temp, * nh, *sfr;
jcl's avatar
jcl committed
136 137
  int * id;
  // new data for Sergey
138
  T * zs, * zsmt, * im, * cm, * ssl;
jcl's avatar
jcl committed
139 140
  int czs, czsmt; // constant for zs and zsmt
  int bits; // to store the bits components
141
  T tframe,redshift;
142
  int ntotmasses;
jcl's avatar
jcl committed
143
  t_io_header_1 header;
144
  CG2_T_Header<T> t_header;
jcl's avatar
jcl committed
145 146 147 148 149 150 151 152 153 154 155 156 157 158
  int npartTotal, npart, npart_total_local;
  int array_vs_file_size; // array vs file 0:same  1:smaler(half) 2:bigger(double)
  inline void checkFileVsArray(const int bytes_to_read, const int size_data, const int npart) {
    int bytes_array = size_data * npart;
    if (bytes_array == bytes_to_read) {
      array_vs_file_size = 0; // same size
    } else
      if (bytes_array < bytes_to_read) {
        array_vs_file_size = 1; // file data type bigger than array data size
        //assert(2*bytes_array == bytes_to_read);
      } else {
        array_vs_file_size = 2; // file data type smaller than array data size
        //assert(bytes_to_read*2 == bytes_array);
      }
159
    if (this->verbose) {
jcl's avatar
jcl committed
160 161 162
      std::cerr << "file_vs_array_size ="<<array_vs_file_size<<" bytes_to_read="<<bytes_to_read<<" bytes_array ="<<bytes_array<<"\n";
    }
  }
163 164 165 166 167 168 169 170 171 172 173 174
  int sizeRealOnFile() {
    int ret;
    switch (array_vs_file_size) {
    case 0 : ret=sizeof(T); break;
    case 1 : ret=sizeof(double); break;
    case 2 : ret=sizeof(float); break;
    default: std::cerr << "Wrong array_vs_file_size ["<<array_vs_file_size<<"]\nabort...";
      std::exit(1);
      break;
    }
    return ret;
  }
jcl's avatar
jcl committed
175 176 177 178 179

  bool isLittleEndian();
  bool swap;
  uns::ComponentRangeVector  crv;
  void storeComponents();
180

jcl's avatar
jcl committed
181
  // member data
182 183 184 185 186 187 188 189 190 191 192 193
  T * getMass()   { return mass; }
  T   getTime()   { return tframe;}
  T   getRedshift() { return redshift;}
  T * getPos()    { return pos; }
  T * getVel()    { return vel; }
  T * getAcc()    { return acc; }
  T * getPot()    { return pot; }
  T * getAge(int & n)   { n=header.npartTotal[4]; return age;}
  T * getMetal(int & n) { n=header.npartTotal[0]+header.npartTotal[4]; return metal;}
  T * getMetalGas(int & n) { n=header.npartTotal[0]; return metal;}
  T * getMetalStars(int & n) { n=header.npartTotal[4]; return metal+header.npartTotal[0];}
  T * getTemp(int & n) { n=header.npartTotal[0]; return temp;}
LAMBERT Jean-charles's avatar
LAMBERT Jean-charles committed
194 195
  T * getNh(int & n) { n=header.npartTotal[0]; return nh;}
  T * getSfr(int & n) { n=header.npartTotal[0]; return sfr;}
196 197 198 199 200 201 202 203 204 205 206 207 208 209
  T * getU(int & n) { n=header.npartTotal[0]; return intenerg;}
  T * getRho(int & n) { n=header.npartTotal[0]; return rho;}
  T * getHsml(int & n) { n=header.npartTotal[0]; return hsml;}
  T * getZs(int & n) { n=czs*(header.npartTotal[0]+header.npartTotal[4]); return zs;}
  T * getZsGas(int & n) { n=czs*(header.npartTotal[0]); return zs;}
  T * getZsStars(int & n) { n=czs*(header.npartTotal[4]); return zs+header.npartTotal[0]*czs;}
  T * getZsmt(int & n) { n=czsmt*(header.npartTotal[0]+header.npartTotal[4]); return zsmt;}
  T * getZsmtGas(int & n) { n=czsmt*(header.npartTotal[0]); return zsmt;}
  T * getZsmtStars(int & n) { n=czsmt*(header.npartTotal[4]); return zsmt+header.npartTotal[0]*czsmt;}
  T * getIm(int & n) { n=header.npartTotal[4]; return im;}
  T * getSsl(int & n) { n=header.npartTotal[4]; return ssl;}
  T * getCm(int & n) { n=1*(header.npartTotal[0]+header.npartTotal[4]); return cm;}
  T * getCmGas(int & n) { n=1*(header.npartTotal[0]); return cm;}
  T * getCmStars(int & n) { n=1*(header.npartTotal[4]); return ((cm==NULL)?cm:cm+header.npartTotal[0]*1);}
210
  bool getHeader(std::string name, T * data);
jcl's avatar
jcl committed
211 212 213 214 215 216 217
  //fortran offset record length
  int frecord_offset;
  //control
  bool is_open, is_read;
  bool status;
  int bytes_counter;
  // method
218
  template <class U> int readCompData(U ** data, const int * index2, const int * npartOffset,
jcl's avatar
jcl committed
219
                                      const int dim, const int nsel);
220 221
  template <class U> int readGasStarsUnknownArray(U ** data, int * n,const int * compOffset);
  template <class U> int readOneArray(U ** data, const int compid,const int * compOffset);
222

223
  template <class U> int readStreamBlock(const std::string blockname, std::vector<U> &data);
jcl's avatar
jcl committed
224 225 226 227 228 229 230 231 232 233 234 235 236
  bool readBlockName();
  std::string block_name;
  int readHeader(const int);

  int readData(char * ptr,const size_t size_bytes,const  int items);
  bool guessVersion();
  int version;
  void unitConversion();
  // skip Block
  inline void skipBlock() {
    int len1 = readFRecord();
    in.seekg(len1,std::ios::cur);
    int len2 = readFRecord();
237
    if (this->verbose) std::cerr << "skipping block name ["<<block_name<<"]\n";
238 239
    if (len2==len1) ; // remove warning....
    assert(len1==len2 && in.good());
jcl's avatar
jcl committed
240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266
    if (block_name == "AGE" || block_name == "Z" ) {
      //std::cerr << "len1 = " << len1 << "\nlen2 = " << len2 << "\n";
    }
  }
  inline void skipData(int len) {
    bytes_counter += (len);
    in.seekg(len,std::ios::cur);
    assert(in.good());
  }
  // swap bytes
  inline void swapBytes(void * x,const int size) {
    char * p=(char *) x;
    for (int i=0;i<size/2;i++) {
     int t=*(p+i); *(p+i)=*(p+size-i-1); *(p+size-i-1)=t;
    }
  }
  // read fortran record
  inline int readFRecord() {
    int len; in.read((char *) &len,sizeof(int));
    // We SWAP data
    if (swap) { // swapping requested
      swapBytes(&len,sizeof(int));
    }
    assert(in.good());
    return len;
  }

267
  }; // end of class CSnapshotGadgetIn
jcl's avatar
jcl committed
268 269 270 271 272 273 274 275 276 277

  //
  //
  //

  // ---------------------------------------------------
  // CSnapshotGadgetout
  // WRITING class
  // ---------------------------------------------------

278
  template <class T> class CSnapshotGadgetOut : public CSnapshotInterfaceOut<T> {
jcl's avatar
jcl committed
279 280 281 282 283 284 285

  public:
    // WRITING constructor
    CSnapshotGadgetOut(const std::string, const std::string, const bool);
    ~CSnapshotGadgetOut();
    int setHeader(void * );
    int setNbody(const int _n);
286 287
    int setData(std::string, T);
    int setData(std::string, const int , T *,const bool _addr=false);
jcl's avatar
jcl committed
288
    // array by double keys
289
    int setData(std::string, std::string, const int , T *,const bool _addr=false);
jcl's avatar
jcl committed
290
    int setData(std::string, std::string, const int , int   *,const bool _addr=false);
291

jcl's avatar
jcl committed
292
    int setData(std::string, const int , int *,const bool _addr=false);
293
    int setData(std::string, const int ,
294
        T *, T *, T *, const bool _addr=false);
295

jcl's avatar
jcl committed
296 297

    int save();
298
    std::vector<double> moveToCom();
jcl's avatar
jcl committed
299 300 301

  private:
    //data
LAMBERT Jean-charles's avatar
LAMBERT Jean-charles committed
302 303
    T * mass[6], * pos[6], *acc[6], * vel[6], * pot[6], * rho, * hsml,
       * age, * metal[6], * intenerg, * temp, * nh, * sfr;
jcl's avatar
jcl committed
304 305 306
    int * id[6];
    int ntot_withmasses;
    std::ofstream out;
307
    // Map to associate EXTRA tags with vectors
308
    std::map<std::string, std::vector <T>  > s_mapStringVector;
jcl's avatar
jcl committed
309 310 311 312 313 314 315 316
    // Map to associate the strings with the bool values
    std::map<std::string, bool> ptrIsAlloc[6];
    t_io_header_1 header;
    int bits, npartTotal;
    int bytes_counter;
    int version;
    int setHeader(t_io_header_1 *);
    int writeData(char * ptr,const size_t size_bytes,const int items);
317
    int writeDataValue(T value, const size_t size_bytes,const int items);
jcl's avatar
jcl committed
318 319 320 321 322 323 324 325 326 327
    void saveFile();
    void setupHeader(bool check=false);
    int writeHeader();
    int write();
    bool writeBlockName(std::string, int);
    inline void writeFRecord(const int len) {
      out.write((char *) &len,sizeof(int));
      assert(out.good());
    }
    // array by keys
328 329 330
    int setMass(std::string, const int _n, T * _mass, const bool addr);
    int setPos (std::string, const int _n, T * _pos , const bool addr);
    int setVel (std::string, const int _n, T * _vel , const bool addr);
jcl's avatar
jcl committed
331
    int setId  (std::string, const int _n, int   * _id  , const bool addr);
332 333 334 335 336 337 338 339
    int setPot (std::string, const int _n, T * _pot , const bool addr);
    int setAcc (std::string, const int _n, T * _acc , const bool addr);

    int setRho (const int _n, T * _rho , const bool addr);
    int setHsml(const int _n, T * _hsml, const bool addr);
    int setU   (const int _n, T * _U   , const bool addr);
    int setAge (const int _n, T * _age   , const bool addr);
    int setTemp(const int _n, T * _temp, const bool addr);
LAMBERT Jean-charles's avatar
LAMBERT Jean-charles committed
340 341
    int setNh  (const int _n, T * _nh, const bool addr);
    int setSfr (const int _n, T * _sfr, const bool addr);
342 343
    int setMetalGas(const int _n, T * _mg   , const bool addr);
    int setMetalStars(const int _n, T * _ms   , const bool addr);
344
    // extra TAGS
345
    int setExtra(std::string tag, const int _n, T * _data, const bool addr);
346
    int setHeader(std::string tag,  T  _data);
347

348
  }; // end of class CSnapshotGadgetOut
jcl's avatar
jcl committed
349 350 351
} // namespace

#endif