Changeset 3289 in Sophya for trunk/Cosmo/SimLSS

Timestamp:

Aug 6, 2007, 6:20:25 PM (18 years ago)

Author:

cmv

Message:

petites mises en formes cmv 06/08/2007

Location:

trunk/Cosmo/SimLSS

Files:

• cmvconcherr.cc (modified) (1 diff)
• cmvobserv3d.cc (modified) (12 diffs)
• cmvtstagn.cc (modified) (1 diff)
• genefluct3d.cc (modified) (3 diffs)
• genefluct3d.h (modified) (2 diffs)

trunk/Cosmo/SimLSS/cmvconcherr.cc

@@ -307 +307 @@
 openppf cmvconcherr.ppf
 
+c++exec int i = 700; \
+cout<<i<<" "<<herrconc.BinCenter(i)<<" "<<herrconc(i)<<" "<<herrconc.Error2(i)<<" "<<herrconc.NEntBin(i)<<endl;
+
 zone 2 2
 disp herrconc "hbincont err"

trunk/Cosmo/SimLSS/cmvobserv3d.cc

@@ -25 +25 @@
      <<" -G : compute Pk(z=0) and apply growth factor in real space"<<endl
      <<"      (default: no, spectrum Pk(z=z_median) for all cube)"<<endl
+     <<" -F : filter spectrum by pixel shape (0=no 1=yes(default)"<<endl
      <<" -x nx,dx : size along x axis (npix,Mpc)"<<endl
      <<" -y ny,dy : size along y axis (npix,Mpc)"<<endl
@@ -94 +95 @@
  bool use_growth_factor = false;
  unsigned short nthread=0;
+ int filter_by_pixel = 1;
 
  // *** What to do
@@ -109 +111 @@
 
  char c;
- while((c = getopt(narg,arg,"ha0PWSV2Gx:y:z:s:Z:M:A:T:")) != -1) {
+ while((c = getopt(narg,arg,"ha0PWSV2GF:x:y:z:s:Z:M:A:T:")) != -1) {
   int nth = 0;
   switch (c) {
@@ -121 +123 @@
     use_growth_factor = true;
     break;
+  case 'F' :
+    sscanf(optarg,"%d",&filter_by_pixel);
+    break;
   case 'x' :
     sscanf(optarg,"%ld,%lf",&nx,&dx);
@@ -177 +182 @@
  cout<<"kmin="<<kmin<<" ("<<lkmin<<"), kmax="<<kmax<<" ("<<lkmax<<") Mpc^-1"
      <<", npt="<<npt<<endl;
+ cout<<"Filter by pixel = "<<filter_by_pixel<<endl;
  cout<<"R="<<R<<" Rg="<<Rg<<" Mpc, sigmaR="<<sigmaR<<endl;
  cout<<"nstar= "<<nstar<<"  mstar="<<mstar<<"  alpha="<<alpha<<endl;
@@ -337 +343 @@
  PrtTim(">>>> End Computing a realization in Fourier space");
 
- if(1) {
-   cout<<"\n--- Checking realization spectra"<<endl;
-   HistoErr hpkgen(0.,knyqmax,nherr);
-   hpkgen.ReCenterBin(); hpkgen.Zero();
-   hpkgen.Show();
-   fluct3d.ComputeSpectrum(hpkgen);
-   {
-   tagobs = "hpkgen"; posobs.PutObject(hpkgen,tagobs);
-   }
-   PrtTim(">>>> End Checking realization spectra");
- }
+ cout<<"\n--- Checking realization spectra"<<endl;
+ HistoErr hpkgen(0.,knyqmax,nherr);
+ hpkgen.ReCenterBin(); hpkgen.Zero();
+ hpkgen.Show();
+ fluct3d.ComputeSpectrum(hpkgen);
+ {
+ tagobs = "hpkgen"; posobs.PutObject(hpkgen,tagobs);
+ }
+ PrtTim(">>>> End Checking realization spectra");
 
  if(comp2dspec) {
@@ -361 +365 @@
  }
 
- if(1) {
+ if(filter_by_pixel!=0) {
+   // A_FAIRE (???) ne faut-il pas dans le cas ou on demande
+   // use_growth_factor filtrer uniquement apres avoir appliquer
+   // le facteur de croissance  ??????
    cout<<"\n--- Computing convolution by pixel shape"<<endl;
    fluct3d.FilterByPixel();
    PrtTim(">>>> End Computing convolution by pixel shape");
- }
-
- if(wfits) {
-   fluct3d.WriteFits("!cmvobserv3d_k0.fits");
-   PrtTim(">>>> End WriteFits");
- }
- if(wppf) {
-   fluct3d.WritePPF("cmvobserv3d_k0.ppf",false);
-   PrtTim(">>>> End WritePPF");
- }
-
- if(1) {
+
    cout<<"\n--- Checking realization spectra after pixel shape convol."<<endl;
    HistoErr hpkgenf(0.,knyqmax,nherr);
@@ -386 +382 @@
    }
    PrtTim(">>>> End Checking realization spectra");
- }
-
- if(comp2dspec) {
-   cout<<"\n--- Checking realization 2D spectra after pixel shape convol."<<endl;
-   Histo2DErr hpkgenf2(0.,ktnyqmax,nherrt,0.,kznyqmax,nherrz);
-   hpkgenf2.ReCenterBin(); hpkgenf2.Zero();
-   hpkgenf2.Show();
-   fluct3d.ComputeSpectrum2D(hpkgenf2);
-   {
-   tagobs = "hpkgenf2"; posobs.PutObject(hpkgenf2,tagobs);
+
+   if(comp2dspec) {
+     cout<<"\n--- Checking realization 2D spectra after pixel shape convol."<<endl;
+     Histo2DErr hpkgenf2(0.,ktnyqmax,nherrt,0.,kznyqmax,nherrz);
+     hpkgenf2.ReCenterBin(); hpkgenf2.Zero();
+     hpkgenf2.Show();
+     fluct3d.ComputeSpectrum2D(hpkgenf2);
+     {
+     tagobs = "hpkgenf2"; posobs.PutObject(hpkgenf2,tagobs);
+     }
+     PrtTim(">>>> End Checking realization 2D spectra");
    }
-   PrtTim(">>>> End Checking realization 2D spectra");
+ }
+
+ if(wfits) {
+   fluct3d.WriteFits("!cmvobserv3d_k0.fits");
+   PrtTim(">>>> End WriteFits");
+ }
+ if(wppf) {
+   fluct3d.WritePPF("cmvobserv3d_k0.ppf",false);
+   PrtTim(">>>> End WritePPF");
  }
 
@@ -414 +419 @@
    cout<<"rgen.Min = "<<rmin<<" , Max="<<rmax<<endl;
    PrtTim(">>>> End Applying growth factor");
- }
+   // A_FAIRE (???) ne faut-il pas maintenant retourner dans fourier
+   // appliquer les filtre des pixels et re-retourner dans reel ?
+ }
+
 
  if(wfits) {
@@ -431 +439 @@
  int_8 nm;
  double rm,rs2;
- if(1) {
-   cout<<"\n--- Check mean and variance in real space"<<endl;
-   fluct3d.NumberOfBad(-1.,1e+200);
-   nm = fluct3d.MeanSigma2(rm,rs2);
-   PrtTim(">>>> End Check mean and variance in real space");
- }
+ cout<<"\n--- Check mean and variance in real space"<<endl;
+ fluct3d.NumberOfBad(-1.,1e+200);
+ nm = fluct3d.MeanSigma2(rm,rs2);
+ PrtTim(">>>> End Check mean and variance in real space");
 
  if(compvarreal) {
@@ -549 +555 @@
  // -- NE PAS FAIRE CA SI ON VEUT CONTINUER LA SIMULATION -> d_rho/rho ecrase
  
- if(1) {
-   cout<<endl<<"\n--- ReComputing spectrum from real space"<<endl;
-   fluct3d.ReComputeFourier();
-   PrtTim(">>>> End ReComputing spectrum");
- }
+ cout<<endl<<"\n--- ReComputing spectrum from real space"<<endl;
+ fluct3d.ReComputeFourier();
+ PrtTim(">>>> End ReComputing spectrum");
 
  if(wfits) {
@@ -564 +568 @@
  }
 
- if(1) {
-   cout<<endl<<"\n--- Computing final spectrum"<<endl;
-   HistoErr hpkrec(0.,knyqmax,nherr);
-   hpkrec.ReCenterBin();
-   hpkrec.Show();
-   fluct3d.ComputeSpectrum(hpkrec);
-   tagobs = "hpkrec"; posobs.PutObject(hpkrec,tagobs);
-   PrtTim(">>>> End Computing final spectrum");
- }
+ cout<<endl<<"\n--- Computing final spectrum"<<endl;
+ HistoErr hpkrec(0.,knyqmax,nherr);
+ hpkrec.ReCenterBin();
+ hpkrec.Show();
+ fluct3d.ComputeSpectrum(hpkrec);
+ tagobs = "hpkrec"; posobs.PutObject(hpkrec,tagobs);
+ PrtTim(">>>> End Computing final spectrum");
 
  if(comp2dspec) {

trunk/Cosmo/SimLSS/cmvtstagn.cc

@@ -104 +104 @@
      sumn += n; sumn2 += n*n;
      double fsum = 0.; unsigned long nf = 0;
-     for(long j=0;j<n;j++) {
+     for(unsigned long j=0;j<n;j++) {
        double f = agn.FluxJY();
        hlflux.Add(log10(f));

trunk/Cosmo/SimLSS/genefluct3d.cc

@@ -462 +462 @@
      S.ReSize(Nx_,Ny_);
      if(s==0) l=0; else if(s==1) l=(Nz_+1)/2; else  l=Nz_-1;
-     sprintf(str,"z%d",l);
+     sprintf(str,"z%ld",l);
      for(i=0;i<Nx_;i++) for(j=0;j<Ny_;j++) S(i,j)=data_[IndexR(i,j,l)];
      pos<<PPFNameTag(str)<<S; S.RenewObjId();
@@ -471 +471 @@
      S.ReSize(Nz_,Nx_);
      if(s==0) j=0; else if(s==1) j=(Ny_+1)/2; else  j=Ny_-1;
-     sprintf(str,"y%d",j);
+     sprintf(str,"y%ld",j);
      for(i=0;i<Nx_;i++) for(l=0;l<Nz_;l++) S(l,i)=data_[IndexR(i,j,l)];
      pos<<PPFNameTag(str)<<S; S.RenewObjId();
@@ -480 +480 @@
      S.ReSize(Nz_,Ny_);
      if(s==0) i=0; else if(s==1) i=(Nx_+1)/2; else  i=Nx_-1;
-     sprintf(str,"x%d",i);
+     sprintf(str,"x%ld",i);
      for(j=0;j<Ny_;j++) for(l=0;l<Nz_;l++) S(l,j)=data_[IndexR(i,j,l)];
      pos<<PPFNameTag(str)<<S; S.RenewObjId();

trunk/Cosmo/SimLSS/genefluct3d.h

@@ -3 +3 @@
 
 #include "machdefs.h"
+#include <math.h>
 #include "genericfunc.h"
 #include "tarray.h"
@@ -67 +68 @@
   double GetKincMin(void)
     {vector<r_8>::const_iterator it = min_element(Dk_.begin(), Dk_.end()); return *it;}
+  double GetKincMax(void)
+    {vector<r_8>::const_iterator it = max_element(Dk_.begin(), Dk_.end()); return *it;}
   double GetKTincMin(void) {return min(Dk_[0],Dk_[1]);}
+  double GetKTincMax(void) {return max(Dk_[0],Dk_[1]);}
 
   void ComputeFourier0(GenericFunc& pk_at_z);