Changeset 3141 in Sophya for trunk/Cosmo/SimLSS/genefluct3d.h


Ignore:
Timestamp:
Jan 17, 2007, 6:44:04 PM (19 years ago)
Author:
cmv
Message:

chgt HProf->HistoErr + spectre 2D cmv 17/01/2007

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Cosmo/SimLSS/genefluct3d.h

    r3134 r3141  
    55#include "genericfunc.h"
    66#include "tarray.h"
    7 #include "hisprof.h"
     7#include "histerr.h"
     8#include "hist2err.h"
    89#include "perandom.h"
     10
     11#include "FFTW/fftw3.h"
     12#include "FitsIO/fitsio.h"
    913
    1014#include <vector>
     
    1216#include "pkspectrum.h"
    1317
    14 #include "FFTW/fftw3.h"
    1518
    1619namespace SOPHYA {
     
    1922class GeneFluct3D {
    2023public:
    21   GeneFluct3D(TArray< complex<r_8 > >& T,PkSpectrumZ& pkz);
     24  GeneFluct3D(TArray< complex<r_8 > >& T);
    2225  virtual ~GeneFluct3D(void);
    2326
     
    2528  void SetSize(long nx,long ny,long nz,double dx,double dy,double dz);  // Mpc
    2629
    27   inline void SetZ(double z) {pkz_.SetZ(z);}
    28   double GetZ(void) {return pkz_.GetZ();}
     30  TArray< complex<r_8> >& GetComplexArray(void) {return T_;}
     31  fftw_complex* GetComplexPointer(void) {return fdata_;}
     32  TArray<r_8>& GetRealArray(void) {return R_;}
     33  r_8* GetRealPointer(void) {return data_;}
     34
     35  // Pour adressage data_[ip]
     36  inline int_8 IndexR(long i,long j,long k) {return (int_8)(k+NTz_*(j+Ny_*i));}
     37  // Pour adressage fdata_[ip][0-1]
     38  inline int_8 IndexC(long i,long j,long k) {return (int_8)(k+NCz_*(j+Ny_*i));}
     39  // - On peut aussi adresser:
     40  // TArray< complex<r_8> >& pk = gf3d.GetComplexArray(); pk(i,j,k) = ...;
     41  // TArray<r_8>& rgen = gf3d.GetRealArray(); rgen(i,j,k) = ...;
     42
     43  vector<long> GetNpix(void) {return N_;}
     44  int_8 NPix(void) {return NRtot_;}
     45
     46  vector<r_8> GetDinc(void) {return D_;}
     47  double GetDVol(void) {return dVol_;}
    2948  double GetVol(void) {return Vol_;}
    30   double GetDVol(void) {return dVol_;}
    31   int_8 NPix(void) {return NRtot_;}
     49
     50  vector<r_8> GetKinc(void) {return Dk_;}
     51  vector<r_8> GetKnyq(void) {return Knyq_;}
    3252  double GetKmax(void) {return sqrt(Knyqx_*Knyqx_+Knyqy_*Knyqy_+Knyqz_*Knyqz_);}
    33   vector<r_8> GetKinc(void)
    34     {vector<r_8> dk; dk.push_back(Dkx_); dk.push_back(Dky_); dk.push_back(Dkz_); return dk;}
    35   vector<r_8> GetKnyq(void)
    36     {vector<r_8> kn; kn.push_back(Knyqx_); kn.push_back(Knyqy_); kn.push_back(Knyqz_); return kn;}
    37   vector<r_8> GetDinc(void)
    38     {vector<r_8> d; d.push_back(Dx_); d.push_back(Dy_); d.push_back(Dz_); return d;}
    39   vector<long> GetNpix(void)
    40     {vector<long> d; d.push_back(Nx_); d.push_back(Ny_); d.push_back(Nz_); return d;}
     53  double GetKTmax(void) {return sqrt(Knyqx_*Knyqx_+Knyqy_*Knyqy_);}
     54  double GetKincMin(void)
     55    {vector<r_8>::const_iterator it = min_element(Dk_.begin(), Dk_.end()); return *it;}
     56  double GetKTincMin(void) {return min(Dk_[0],Dk_[1]);}
    4157
    42   void ComputeFourier0(void);
    43   void ComputeFourier(void);
     58  void ComputeFourier0(GenericFunc& pk_at_z);
     59  void ComputeFourier(GenericFunc& pk_at_z);
    4460  void FilterByPixel(void);
    45   int  ComputeSpectrum(HProf& hp);
     61
    4662  void ComputeReal(void);
     63
    4764  void ReComputeFourier(void);
     65
     66  int  ComputeSpectrum(HistoErr& herr);
     67  int  ComputeSpectrum2D(Histo2DErr& herr);
    4868
    4969  int_8 VarianceFrReal(double R,double& var);
     
    5979  void AddNoise2Real(double snoise);
    6080
     81  void WriteFits(string cfname,int bitpix=FLOAT_IMG);
     82  void ReadFits(string cfname);
     83
     84  void WritePPF(string cfname,bool write_real=true);
     85  void ReadPPF(string cfname);
     86
    6187  void Print(void);
    6288
    6389protected:
     90  void setsize(long nx,long ny,long nz,double dx,double dy,double dz);
     91  void setalloc(void);
     92  void setpointers(bool from_real);
    6493  long manage_coefficients(void);
    6594  double compute_power_carte(void);
     95  void check_array_alloc(void);
    6696  inline double pixelfilter(double x)
    6797    {return (x<0.025) ? 1.-x*x/6.*(1.-x*x/20.): sin(x)/x;}
    6898
    6999
    70   long Nx_,Ny_,Nz_;
    71   sa_size_t SzK_[3];
     100  long Nx_,Ny_,Nz_;  vector<long> N_;
    72101  long NCz_,NTz_;
    73102  int_8 NRtot_;
    74   double Dx_,Dy_,Dz_;
     103
     104  double Dx_,Dy_,Dz_;  vector<double> D_;
     105
     106  double Dkx_,Dky_,Dkz_;  vector<double> Dk_;
     107  double Knyqx_,Knyqy_,Knyqz_;  vector<double> Knyq_;
     108  double Dk3_;
    75109  double dVol_, Vol_;
    76   double Dkx_,Dky_,Dkz_;
    77   double Knyqx_,Knyqy_,Knyqz_;
    78   double Dk3_;
    79110
    80111  fftw_plan pf_,pb_;
    81112  unsigned short nthread_;
    82113
    83   //   = 0 : T_ vide
    84   //   = 1 : T_ contient une realisation du spectre (espace fourier)
    85   //   = 2 : T_ contient une realisation dans l'espace reel
    86   //   = 3 : T_ contient la TF d'une realisation dans l'espace reel (d_rho/rho)
    87   //   = 4 : T_ contient (1 + d_rho/rho)
    88   //   = ... etc ...
    89   unsigned short Tcontent_;
     114  bool array_allocated_;  // true if array has been allocated
    90115
    91   TArray< complex<r_8 > >& T_;
    92   PkSpectrumZ& pkz_;
     116  TArray< complex<r_8> >& T_;
     117  fftw_complex *fdata_;
     118  TArray<r_8> R_;
     119  double *data_;
    93120};
    94121
Note: See TracChangeset for help on using the changeset viewer.