Changeset 3318 in Sophya for trunk/Cosmo/SimLSS/cmvtransf.cc

Timestamp:

Aug 28, 2007, 11:41:49 AM (18 years ago)

Author:

cmv

Message:

TransfertTabulate pour avoir les fct de transfert de CMBFast cmv 28/08/2007

File:

• trunk/Cosmo/SimLSS/cmvtransf.cc (modified) (5 diffs)

trunk/Cosmo/SimLSS/cmvtransf.cc

@@ -14 +14 @@
 
 void usage(void);
-void usage(void) {cout<<"cmvtuniv [k1,k2,npt] [h100,Om0,Ob0,tcmb]"<<endl;}
+void usage(void) {
+cout
+<<"cmvtransf [options]"<<endl
+<<"  -k k1,k2,npt : k range in Mpc^-1"<<endl
+<<"  -U h100,Om0,Ob0,tcmb : cosmology"<<endl
+<<"  -F filename : read also CMBFast file"<<endl
+<<endl;
+}
 
 int main(int narg,char *arg[])
@@ -20 +27 @@
  double k1 = 1e-6, k2 = 10.; int_4 npt = 100;
  double h100 = 0.71, Om0 = 0.267804, Ob0 = 0.0444356, tcmb = T_CMB_Par;
+ string fcmbfast = "";
 
- if(narg>1) sscanf(arg[1],"%lf,%lf,%d",&k1,&k2,&npt);
+ char c;
+  while((c = getopt(narg,arg,"hk:U:F:")) != -1) {
+  switch (c) {
+  case 'k' :
+    sscanf(optarg,"%d,%lf,%lf",&npt,&k1,&k2);
+    break;
+  case 'U' :
+    sscanf(optarg,"%lf,%lf,%lf,%lf",&h100,&Om0,&Ob0,&tcmb);
+    break;
+  case 'F' :
+    fcmbfast = optarg;
+    break;
+  case 'h' :
+  default :
+    usage(); return -1;
+  }
+ }
+
  if(k1<=0.) {cout<<"k1 must be >0"<<endl; return -1;}
- if(npt<=0) npt = 100;
+ if(npt<=0) npt = 1000;
+ if(k2<=k1) {k1=1.e-4; k2=1.e+2;}
  cout<<"k1="<<k1<<"  k2="<<k2<<" npt="<<npt<<endl;
- if(narg>2) sscanf(arg[2],"%lf,%lf,%lf,%lf",&h100,&Om0,&Ob0,&tcmb);
  cout<<"h100="<<h100<<" Om0="<<Om0<<" Ob0="<<Ob0<<" Tcmb="<<tcmb<<endl;
+ cout<<"fcmbfast="<<fcmbfast<<endl;
 
  cout<<"TransfertEisenstein with baryon"<<endl;
@@ -43 +69 @@
  TransfertEisenstein tfnob(h100,Om0-Ob0,Ob0,tcmb,true);  // No baryons
 
- const int n = 5;
- char *vname[n] = {"k","t","tnosc1","tnosc2","tnob"};
+ TransfertTabulate tfcmbfast(h100,Om0-Ob0,Ob0);
+ bool cmbfastOK = false;
+ if(fcmbfast.size()>0) {
+   cout<<endl<<"TransfertTabulate for CMBfast"<<endl;
+   tfcmbfast.ReadCMBFast(fcmbfast);
+   tfcmbfast.SetInterpTyp(1);
+   if(tfcmbfast.NPoints()>0) cmbfastOK = true;
+ }
+
+ const int n = 6;
+ char *vname[n] = {"k","t","tnosc1","tnosc2","tnob","tcmbf"};
  NTuple nt(n,vname);
  double xnt[n];
+ for(int i=0;i<n;i++) xnt[i]=0.;
 
  double lnk1 = log10(k1), lnk2=log10(k2), dlnk=(lnk2-lnk1)/npt;
@@ -56 +92 @@
    xnt[3]  = tfnosc2(k);
    xnt[4]  = tfnob(k);
+   if(cmbfastOK) xnt[5]  = tfcmbfast(k);
    nt.Fill(xnt);
  }
@@ -73 +110 @@
 set k log10(k)
 
-n/plot nt.t%$k ! ! "nsta crossmarker3"
-n/plot nt.tnob%$k ! ! "nsta red crossmarker3 same"
-n/plot nt.tnosc1%$k ! ! "nsta blue circlemarker3 same"
-n/plot nt.tnosc2%$k ! ! "nsta green trianglemarker3 same"
+# Eisenstein
+n/plot nt.t%$k ! ! "nsta connectpoints"
+n/plot nt.tnob%$k ! ! "nsta red same connectpoints"
+n/plot nt.tnosc2%$k ! ! "nsta blue same connectpoints"
+n/plot nt.tnosc1%$k ! ! "nsta green same connectpoints"
 
-n/plot nt.t/tnosc1%$k ! ! "nsta crossmarker3"
-n/plot nt.t/tnosc2%$k ! ! "nsta crossmarker3 red same"
+n/plot nt.t/tnob%$k   tnob>0   ! "nsta red connectpoints"
+n/plot nt.t/tnosc2%$k tnosc2>0 ! "nsta blue same connectpoints"
+n/plot nt.t/tnosc1%$k tnosc1>0 ! "nsta green same connectpoints"
+
+# CMBFast
+n/plot nt.t%$k ! ! "nsta connectpoints"
+n/plot nt.tcmbf%$k ! ! "nsta red same connectpoints"
+
+n/plot nt.(tcmbf-t)%$k tcmbf>0. ! "nsta connectpoints"
+n/plot nt.(tcmbf-t)/tcmbf%$k tcmbf>0. ! "nsta connectpoints"
 */