Changeset 3787 in Sophya for trunk/Cosmo/RadioBeam/calcpk2.cc

Timestamp:

Jun 25, 2010, 12:00:30 PM (15 years ago)

Author:

ansari

Message:

Ajout classes / programmes de calcul d'effet de lobe sur les radio-sources, Reza 25/06/2010

File:

• trunk/Cosmo/RadioBeam/calcpk2.cc (modified) (8 diffs)

trunk/Cosmo/RadioBeam/calcpk2.cc

@@ -7 +7 @@
     R. Ansari , C. Magneville - Juin 2010 
 
-Usage: calcpk2 InMapLSS convFacLSS InMapSync convFacSync OutPk 
+Usage: calcpk2 InMapLSS convFacLSS InMapSync convFacSync InMapRadioSource convFacRsc OutPk [PixNoiseLevel]
 ---------------------------------------------------------------  */
 
@@ -22 +22 @@
 
 #include "qhist.h" 
+#include "lobe.h" 
+#include "cubedef.h" 
+
 #include "histinit.h"
 #include "fftwserver.h"
@@ -30 +33 @@
 typedef ThSDR48RandGen RandomGenerator ;
 //--- Declaration des fonctions 
-TArray<r_4> CleanForeground(TArray<r_4>& maplss, TArray<r_4>& mapsync, TArray<r_8>& synctemp, TArray<r_8>& specidx);
+TArray<r_4> CleanForeground(TArray<r_4>& maplss, TArray<r_4>& mapsync, double freq0, double dfreq, 
+                            TArray<r_8>& synctemp, TArray<r_8>& specidx);
 
 //-------------------------------------------------------------------------
@@ -39 +43 @@
 {
   if (narg<6) {
-    cout << " Usage: calcpk2 InMapLSS convFacLSS InMapSync convFacSync OutPk  " << endl;
+    cout << " Usage: calcpk2 InMapLSS convFacLSS InMapSync convFacSync  OutPk [PixNoiseLevel] " << endl;
     return 1;
   }
@@ -50 +54 @@
     r_4 rfacsync = atof(arg[4]);
     string outname = arg[5];
+
+    double pixsignoise = 0.;
+    bool fgaddnoise=false;
+    if (narg>6) {
+      pixsignoise=atof(arg[6]);
+      fgaddnoise=true;
+    }
+
+
     TArray<r_4> inmaplss, inmapsync;
     const char * tits[2]={"LSS", "Sync"};
@@ -68 +81 @@
       cout << " ... Mean=" << mean << " Sigma=" << sigma << endl;
     }
+
+    if (fgaddnoise) { 
+      cout << " calcpk2: adding noise to LSS input cube ... " << endl;
+      BeamEffect::AddNoise(inmaplss, pixsignoise);
+    }
+
     tm.Split(" After read ");
     TArray<r_8> synctemp, specidx; 
     cout << "calcpk2[3] : calling CleanForeground(...) " << endl;
-    TArray<r_4> inmap = CleanForeground(inmaplss, inmapsync, synctemp, specidx);
+    double freq0 = Freq0MHz;
+    double dfreq = FreqSizeMHz/(double)NFreq;
+    TArray<r_4> inmap = CleanForeground(inmaplss, inmapsync, freq0, dfreq, synctemp, specidx);
     double mean, sigma;
     MeanSigma(inmap, mean, sigma);
@@ -117 +138 @@
 
 /* --Fonction-- */
-TArray<r_4> CleanForeground(TArray<r_4>& maplss, TArray<r_4>& mapsync, TArray<r_8>& synctemp, TArray<r_8>& specidx)
+TArray<r_4> CleanForeground(TArray<r_4>& maplss, TArray<r_4>& mapsync, double freq0, double dfreq, 
+                            TArray<r_8>& synctemp, TArray<r_8>& specidx)
+// Inputs : maplss, mapsyc, freq0, dfreq
+// Outputs : synctemp, specidx  (reconstructed foreground temperature and spectral index
+// Return_Array : foreground subtracted LSS signal 
 {
   bool smo;
@@ -131 +156 @@
   Vector vlnf(maplss.SizeZ());
   int nprt = 0;
-  double freq0=840.;  // Frequence premier index en k (MHz)
-  double dfreq=1.;   // largeur en frequence de chaque plan (Mhz)  
+  // double freq0 : Frequence premier index en k (MHz)
+  // double dfreq :   // largeur en frequence de chaque plan (Mhz)  
   for(sa_size_t i=0; i<maplss.SizeX(); i++) 
     for(sa_size_t j=0; j<maplss.SizeY(); j++)  {