Changeset 3367 in Sophya for trunk/Cosmo/SimLSS/cmvdefsurv.cc

Timestamp:

Oct 30, 2007, 6:06:49 PM (18 years ago)

Author:

cmv

Message:

gestion correct du lobe / taille pixel, evolution bruit Msol en fct de z, cmv 30/10/2007

File:

• trunk/Cosmo/SimLSS/cmvdefsurv.cc (modified) (8 diffs)

trunk/Cosmo/SimLSS/cmvdefsurv.cc

@@ -28 +28 @@
 inline double min2rad(double tmin) {return tmin*M_PI/(180.*60.);}
 inline double sec2rad(double tsec) {return tsec*M_PI/(180.*3600.);}
+
+inline double sr2deg2(double asr) {return asr*pow(rad2deg(1.),2.);}
+inline double sr2min2(double asr) {return asr*pow(rad2min(1.),2.);}
+inline double deg22sr(double adeg) {return adeg*pow(deg2rad(1.),2.);}
+inline double min22sr(double amin) {return amin*pow(min2rad(1.),2.);}
 
 double LargeurDoppler(double v /* km/s */, double nu);
@@ -54 +59 @@
       <<" -O surf,tobs,eta : surface effective (m^2), temps d\'observation (s), efficacite d\'ouverture"<<endl
       <<" -N Tsys : temperature du system (K)"<<endl
-      <<" -L lobewidth,freqlob : taille du lobe d\'observation (FWHM) en arcmin (def= rien)"<<endl
-      <<"                        pour la frequence freqlob en MHz"<<endl
-      <<"            Si freqlob<=0 : la frequence de reference est celle du redshift etudie"<<endl
-      <<"            Si freqlob absent : la frequence de reference 1.4 GHz"<<endl
-      <<"            lobewidth est pris comme la FWHM du lobe gaussien equivalent"<<endl
+      <<" -L lobearea,freqlob : angle solide du lobe d\'observation en arcmin^2 (def= celle du pixel)"<<endl
+      <<"                       pour la frequence freqlob en MHz"<<endl
+      <<"            Si freqlob<=0 : la frequence de reference est celle du redshift etudie (def)"<<endl
       <<" -2 : two polarisations measured"<<endl
       <<" -M  : masse de HI de reference (MSol), si <=0 mean schechter in pixel"<<endl
@@ -107 +110 @@
  double tobs = 6000., surftot = 400000., eta_eff  = 1.;
  // taille du lobe d'observation en arcmin pour la frequence
- double lobewidth0 = -1., lobefreq0 = Fr_HyperFin_Par*1.e3;
+ double lobearea0 = -1., lobefreq0 = -1.;
 
  //-- Variance cosmique (default = standard SDSSII)
@@ -138 +141 @@
     break;
   case 'L' :
-    sscanf(optarg,"%lf,%lf",&lobewidth0,&lobefreq0);
+    sscanf(optarg,"%lf,%lf",&lobearea0,&lobefreq0);
     break;
   case 'N' :
@@ -363 +366 @@
  cout<<"Pixel solid angle = "<<angsol_pix[0]<<" sr ("<<angsol_pix[0]/(4.*M_PI)<<" *4Pi)"
      <<" in ["<<angsol_pix[1]<<","<<angsol_pix[2]<<"]"<<endl;
+ cout<<"                  = "<<sr2min2(angsol_pix[0])<<"\'^2"
+     <<" in ["<<sr2min2(angsol_pix[1])<<","<<sr2min2(angsol_pix[2])<<"]"<<endl;
  cout<<"Total solid angle = "<<angsol_tot[0]<<" sr ("<<angsol_tot[0]/(4.*M_PI)<<" *4Pi)"
      <<" in ["<<angsol_tot[1]<<","<<angsol_tot[2]<<"]"<<endl;
@@ -452 +457 @@
  double angSEq[4], angEq[4];
  for(int i=0;i<4;i++) {
-   double nu =Fr_HyperFin_Par;
+   double nu = Fr_HyperFin_Par;
    if(i<3) nu = nuhiz[i];
    angSEq[i] = AngsolEqTelescope(nu,surftot);
@@ -464 +469 @@
      <<"        angular diameter = "<<angEq[3]<<" rd = "<<rad2min(angEq[3])<<"\'"<<endl;
 
- // Pour une observation ou le lobe est + petit que le pixel de simulation
+ // Pour une observation ou le lobe est + petit ou grand que le pixel de simulation
  // La taille angulaire du lobe change avec la frequence donc avec le redshift
- // La taille du lobe "lobewidth" est donnee pour une frequence de reference "lobefreq0"
+ // La taille du lobe "lobearea0" est donnee pour une frequence de reference "lobefreq0" en MHz
  double nlobes[3] = {1.,1.,1.};
- if(lobewidth0>0.) {
-   if(lobefreq0<=0.) lobefreq0 = nuhiz[0]*1.e3; // MHz
-   cout<<"\n--- Lobe: width="<<lobewidth0<<"\' pour "<<lobefreq0<<" MHz"<<endl;
-   double lobewidth[3], slobe[3], lobecyl[3], lobearea[3];
-   for(int i=0;i<3;i++) {
-     lobewidth[i] = lobewidth0*lobefreq0/(nuhiz[i]*1.e3); // arcmin
-     slobe[i] = lobewidth[i]/2.35482; // sigma du lobe en arcmin
-     lobecyl[i] = sqrt(8.)*slobe[i]; // diametre du lobe cylindrique equiv en arcmin
-     lobearea[i] = M_PI*lobecyl[i]*lobecyl[i]/4.; // en arcmin^2 (hypothese lobe gaussien)
-     nlobes[i] = rad2min(adtx[i])*rad2min(adty[i])/lobearea[i];
-   }
-   cout<<"Lobe width changed to "<<lobewidth[0]<<"\' pour "<<nuhiz[0]*1.e3<<" MHz"
-       <<" in ["<<lobewidth[1]<<","<<lobewidth[2]<<"]"<<endl;
-   cout<<"Lobe sigma "<<slobe[0]<<"\' in ["<<slobe[1]<<","<<slobe[2]<<"]"<<endl;
-   cout<<"Lobe cylindrique diametre "<<lobecyl[0]<<"\' in ["<<lobecyl[1]<<","<<lobecyl[2]<<"]"<<endl;
-   cout<<"Lobe cylindrique area "<<lobearea[0]<<"\'^2 in ["<<lobearea[1]<<","<<lobearea[2]<<"]"<<endl;
-   cout<<"Number of beams in one transversal pixel "<<nlobes[0]<<" in ["<<nlobes[1]<<","<<nlobes[2]<<"]"<<endl;
- }
+ // Si "lobefreq0" negatif c'est la frequence du centre du cube
+ if(lobefreq0<=0.) lobefreq0 = nuhiz[0]*1.e3;
+ // Si "lobearea0" negatif c'est la taille du pixel ramenee a lobefreq0
+ if(lobearea0<=0.) lobearea0 = sr2min2(angsol_pix[0])*pow(nuhiz[0]*1.e3/lobefreq0,2.);
+ cout<<"\n--- Lobe: angsol="<<lobearea0<<"\'^2 pour "<<lobefreq0<<" MHz"<<endl;
+ double lobearea[3];
+ for(int i=0;i<3;i++) {
+   lobearea[i] = lobearea0*pow(lobefreq0/(nuhiz[0]*1.e3),2.);
+   nlobes[i] = sr2min2(angsol_pix[i])/lobearea[i];
+ }
+ cout<<"Lobe cylindrique area "<<lobearea[0]<<"\'^2 in ["<<lobearea[1]<<","<<lobearea[2]<<"]"<<endl;
+ cout<<"Number of beams in one transversal pixel "<<nlobes[0]<<" in ["<<nlobes[1]<<","<<nlobes[2]<<"]"<<endl;
 
  // --- 
@@ -765 +764 @@
        <<"   Pk      = "<<pkbruit[0]<<" Mpc^3 in ["<<pkbruit[1]<<","<<pkbruit[2]<<"]"<<endl;
 
-   if(fabs(nlobes[0]-1.)>1e-6) {
-     double slim_nl[3], SsN_nl[3], smass_nl[3], pkbruit_nl[3];
-     for(int i=0;i<3;i++) {
-       slim_nl[i] = slim[i] * sqrt(nlobes[i]);
-       SsN_nl[i] = ssig_2polar[i] / slim_nl[i];
-       smass_nl[i] = mhiref / ssig_2polar[i] * slim_nl[i];
-       pkbruit_nl[i] = pow(smass_nl[i]/mhiref,2.)*vol_pixel;
-     }
-     cout<<"for "<<nlobes[0]<<" lobes[i] in ["<<nlobes[1]<<","<<nlobes[2]<<"] :"<<endl
+   double slim_nl[3], SsN_nl[3], smass_nl[3], pkbruit_nl[3];
+   for(int i=0;i<3;i++) {
+     slim_nl[i] = slim[i] * sqrt(nlobes[i]);
+     SsN_nl[i] = ssig_2polar[i] / slim_nl[i];
+     smass_nl[i] = mhiref / ssig_2polar[i] * slim_nl[i];
+     pkbruit_nl[i] = pow(smass_nl[i]/mhiref,2.)*vol_pixel;
+   }
+   cout<<"for "<<nlobes[0]<<" lobes[i] in ["<<nlobes[1]<<","<<nlobes[2]<<"] :"<<endl
        <<"   Slim    = "<<slim_nl[0]*1.e6<<" mu_Jy in ["<<slim_nl[1]*1.e6<<","<<slim_nl[2]*1.e6<<"]"<<endl
        <<"   S/N     = "<<SsN_nl[0]<<" in ["<<SsN_nl[1]<<","<<SsN_nl[2]<<"]"<<endl
        <<"   Mass HI = "<<smass_nl[0]<<" Msol in ["<<smass_nl[1]<<","<<smass_nl[2]<<"]"<<endl
        <<"   Pk      = "<<pkbruit_nl[0]<<" Mpc^3 in ["<<pkbruit_nl[1]<<","<<pkbruit_nl[2]<<"]"<<endl;
-   }
 
  }