Changeset 3261 in Sophya for trunk/Cosmo/SimLSS

Timestamp:

Jun 4, 2007, 7:02:03 PM (18 years ago)

Author:

cmv

Message:

add sigma computation with bad pixel cmv 04/06/2007

Location:

trunk/Cosmo/SimLSS

Files:

• cmvdefsurv.cc (modified) (3 diffs)
• cmvobserv3d.cc (modified) (4 diffs)
• genefluct3d.cc (modified) (4 diffs)
• genefluct3d.h (modified) (1 diff)

trunk/Cosmo/SimLSS/cmvdefsurv.cc

@@ -301 +301 @@
      <<","<<hifactor*Msol2FluxHI(mhiref,dlumlim[1])<<"] W/m^2"<<endl;
  double sfhi = fhi / (dnuhiz*1e9) / Jansky2Watt_cst;
- cout<<"If spread over "<<dnuhiz<<" GHz, flux density = "<<sfhi<<" Jy"<<endl;
+ cout<<"If spread over pixel depth ("<<dnuhiz<<" GHz), flux density = "<<sfhi<<" Jy"<<endl;
 
  // --- Signal analysis
@@ -337 +337 @@
  cout<<"Synchrotron: T="<<Tsynch408<<" K ("<<nuhaslam<<" GHz), "
      <<tsynch<<" K ("<<nuhiz<<" GHz)"<<endl
-     <<"             P="<<psynch<<" W"<<endl;
- cout<<"             flux density = "<<ssynch<<" Jy"<<endl;
+     <<"             P="<<psynch<<" W for pixel"<<endl;
+ cout<<"             flux density = "<<ssynch<<" Jy for pixel solid angle"<<endl;
 
  // CMB
@@ -345 +345 @@
  double pcmb = facpolar * planck(nuhiz*1.e+9) * surfeff * angsol * (dnuhiz*1e9);
  double scmb = pcmb / surfeff / (dnuhiz*1.e+9) / Jansky2Watt_cst;
- cout<<"CMB: T="<<tcmb<<" K -> P="<<pcmb<<" W"<<endl;
- cout<<"                       flux density = "<<scmb<<" Jy"<<endl;
+ cout<<"CMB: T="<<tcmb<<" K -> P="<<pcmb<<" W for pixel"<<endl;
+ cout<<"                       flux density = "<<scmb<<" Jy for pixel solid angle"<<endl;
 
  // AGN

trunk/Cosmo/SimLSS/cmvobserv3d.cc

@@ -435 +435 @@
  cout<<"1+rgen: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
      <<rs2<<" -> "<<sqrt(rs2)<<endl;
+ nm = fluct3d.MeanSigma2(rm,rs2,0.,1e200,true,0.);
+ cout<<"1+rgen with_neg_a_zero: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
+     <<rs2<<" -> "<<sqrt(rs2)<<endl;
  PrtTim(">>>> End Converting fluctuations into mass");
 
@@ -441 +444 @@
  cout<<rm<<" galaxies put into survey"<<endl;
  nm = fluct3d.MeanSigma2(rm,rs2,0.);
- cout<<"galaxy: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
+ cout<<"galaxy mass: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
+     <<rs2<<" -> "<<sqrt(rs2)<<endl;
+ nm = fluct3d.MeanSigma2(rm,rs2,0.,1e200,true,0.);
+ cout<<"galaxy mass with_neg_a_zero: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
      <<rs2<<" -> "<<sqrt(rs2)<<endl;
  PrtTim(">>>> End Converting mass into galaxy number");
 
- cout<<"\n--- Set negative pixels to BAD"<<endl;
+ cout<<"\n--- Set negative and null pixels to BAD"<<endl;
  nm = fluct3d.SetToVal(0.,1e+200,-999.);
  cout<<nm<<" negative in survey set to BAD"<<endl;
  nm = fluct3d.MeanSigma2(rm,rs2,-998.);
- cout<<"galaxy: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
+ cout<<"galaxy mass: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
      <<rs2<<" -> "<<sqrt(rs2)<<endl;
  PrtTim(">>>> End Set negative pixels to BAD etc...");
@@ -457 +463 @@
  cout<<nm<<" galaxies into survey after poisson"<<endl;
  nm = fluct3d.MeanSigma2(rm,rs2,-998.);
- cout<<"galaxy: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
+ cout<<"galaxy number: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
+     <<rs2<<" -> "<<sqrt(rs2)<<endl;
+ nm = fluct3d.MeanSigma2(rm,rs2,0.,1e200,true,0.);
+ cout<<"galaxy number with_neg_a_zero: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
      <<rs2<<" -> "<<sqrt(rs2)<<endl;
  PrtTim(">>>> End Apply poisson on galaxy number");
@@ -474 +483 @@
  double mhi = fluct3d.TurnNGal2Mass(tirhmdndm,true);
  cout<<mhi<<" MSol in survey / "<<mass_by_mpc3*fluct3d.GetVol()<<endl;
- nm = fluct3d.MeanSigma2(rm,rs2,0.);
+ nm = fluct3d.MeanSigma2(rm,rs2,-998.);
  cout<<"HI mass: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
      <<rs2<<" -> "<<sqrt(rs2)<<endl;
  cout<<"Equivalent: "<<rm*nm/fluct3d.NPix()<<" Msol / pixels"<<endl;
+ nm = fluct3d.MeanSigma2(rm,rs2,0.,1e200,true,0.);
+ cout<<"HI mass with_neg_a_zero: ("<<nm<<") Mean = "<<rm<<", Sigma^2 = "
+     <<rs2<<" -> "<<sqrt(rs2)<<endl;
  PrtTim(">>>> End Convert Galaxy number to HI mass");
 

trunk/Cosmo/SimLSS/genefluct3d.cc

@@ -854 +854 @@
 }
 
-int_8 GeneFluct3D::MeanSigma2(double& rm,double& rs2,double vmin,double vmax)
-// mean,sigma^2 pour pixels avec valeurs ]vmin,vmax[ extremites exclues
-//   -> mean and sigma^2 are NOT computed with pixels values vmin and vmax
-{
+int_8 GeneFluct3D::MeanSigma2(double& rm,double& rs2,double vmin,double vmax
+                             ,bool useout,double vout)
+// Calcul de mean,sigma2 dans le cube reel avec valeurs ]vmin,vmax[ extremites exclues
+// useout = false: ne pas utiliser les pixels hors limites pour calculer mean,sigma2
+//          true : utiliser les pixels hors limites pour calculer mean,sigma2
+//                 en remplacant leurs valeurs par "vout"
+{
+ bool tstval = (vmax>vmin)? true: false;
+ if(lp_>0) {
+   cout<<"--- MeanSigma2: ";
+   if(tstval) cout<<"range=]"<<vmin<<","<<vmax<<"[";
+   if(useout) cout<<", useout="<<useout<<" vout="<<vout;
+   cout<<endl;
+ }
  check_array_alloc();
 
@@ -866 +876 @@
    int_8 ip = IndexR(i,j,l);
    double v = data_[ip];
-   if(v<=vmin || v>=vmax) continue;
+   if(tstval) {
+     if(v<=vmin || v>=vmax) {if(useout) v=vout; else continue;}
+   }
    rm += v;
    rs2 += v*v;
@@ -877 +889 @@
  }
 
+ if(lp_>0) cout<<"  n="<<n<<" m="<<rm<<" s2="<<rs2<<" s="<<sqrt(fabs(rs2))<<endl;
+
  return n;
 }
@@ -882 +896 @@
 int_8 GeneFluct3D::SetToVal(double vmin, double vmax,double val0)
 // set to "val0" if out of range ]vmin,vmax[ extremites exclues
-//     ->  vmin and vmax are set to val0
-{
+// cad set to "val0" if in [vmin,vmax] -> vmin and vmax are set to val0
+{
+ if(lp_>0) cout<<"--- SetToVal set to="<<val0
+               <<" when in range=["<<vmin<<","<<vmax<<"]"<<endl;
  check_array_alloc();
 

trunk/Cosmo/SimLSS/genefluct3d.h

@@ -74 +74 @@
 
   int_8 VarianceFrReal(double R,double& var);
-  int_8 MeanSigma2(double& rm,double& rs2,double vmin=-1.e+150,double vmax=1.e+150);
+  int_8 MeanSigma2(double& rm,double& rs2,double vmin=1.,double vmax=-1.
+                  ,bool useout=false,double vout=0.);
   int_8 NumberOfBad(double vmin=-1.e+150,double vmax=1.e+150);
   int_8 SetToVal(double vmin, double vmax,double val0=0.);