Changeset 607 in Sophya for trunk/SophyaLib/SkyT/skymixer.cc

Timestamp:

Nov 20, 1999, 10:00:54 PM (26 years ago)

Author:

ansari

Message:

Modifs preparatoire pour Garching MAP , Reza 20/11/99

File:

• trunk/SophyaLib/SkyT/skymixer.cc (modified) (8 diffs)

trunk/SophyaLib/SkyT/skymixer.cc

@@ -32 +32 @@
                   PixelMap<double>& mapToAdd, RadSpectra& rs, double K=1.);
 
+int CheckCards(DataCards & dc, string & msg);
+char * BuildFITSFileName(string const & fname);
+void RadSpec2Nt(RadSpectra & rs, POutPersist & so, string name);
+void SpectralResponse2Nt(SpectralResponse& sr, POutPersist & so, string name);
+
+//  ----- Variable globale ------------
+static char mapPath[256];     // Path for input maps
+static int  hp_nside = 32;    // HealPix NSide
+static int  nskycomp = 0;     // Number of sky components
+static int  debuglev = 0;     // debuglevel
+
 // -------------------------------------------------------------------------
 //                             main program 
@@ -37 +48 @@
 int main(int narg, char * arg[])
 {
-  if ((narg < 2) || (strcmp(arg[1], "-h") == 0) ) {
-    cout << " skymixer / Error args \n Usage: skymixer parameterFile " << endl;
+  if ((narg < 3) || ((narg > 1) && (strcmp(arg[1], "-h") == 0) )) {
+    cout << "  Usage: skymixer parameterFile outputfitsname [outppfname]" << endl;
     exit(0);
   }
@@ -46 +57 @@
 string msg;
 int rc = 0;
-  
+POutPersist * so = NULL;
+
 try {
   string dcard = arg[1];
@@ -52 +64 @@
   DataCards dc(dcard);
 
-  //  Cheking datacards
-  if ( (!dc.HasKey("SKYMIX")) || (!dc.HasKey("FILTER")) ) {
-     rc = 71; 
-     msg = "Invalid parameters  - NO @SKYMIX or @FILTER card ";
-     goto problem;
-  }
-
-  // Decoding number of component and pixelisation parameter
-  int mg = 32;
-  int ncomp = 0;
-  ncomp = dc.IParam("SKYMIX", 0, 0);
-  mg = dc.IParam("SKYMIX", 1, 32);
-  if (ncomp < 1)  {
-    msg = "Invalid parameters  - Check datacards @SKYMIX ";
-    rc = 72;
-    goto problem;
-  }
-
-  int kc;
-  string key;
-  char buff[256];
-  bool pb = false;
-  for(kc=0; kc<ncomp; kc++) {
-    sprintf(buff, "SCMFILE%d", kc+1);
-    key = buff;
-    if (dc.NbParam(key) < 1) {   
-      msg = "Missing or invalid card : " + key;
-      pb = true;  break;
-    }
-    sprintf(buff, "SCSPEC%d", kc+1);
-    key = buff;
-    if (dc.NbParam(key) < 1)  {   
-      msg = "Missing or invalid card : " + key;
-      pb = true;  break;
-    }
-
-    }
-
-  if (pb)  {
-    rc = 72;
-    goto problem;
-  }
-
-  cout << " skymix/Info : NComp = " <<  ncomp << " M_SphereGorski= " << mg << endl;
-
-  SphereGorski<float> outgs(mg);
-  cout << " Creating Output Gorski Map  NbPixels= " << outgs.NbPixels() << endl;
+  rc = CheckCards(dc, msg);
+  if (rc) goto Fin;
+
+  cout << " skymix/Info : NComp = " <<  nskycomp << " SphereGorski_NSide= " << hp_nside << endl;
+  cout << "  ... MapPath = " << (string)mapPath << "  DbgLev= " << debuglev << endl;
+
+// We create an output persist file for writing debug objects
+  if (debuglev > 0) so = new POutPersist("skymixdbg.ppf");
+
+  SphereGorski<float> outgs(hp_nside);
+  cout << " Output Gorski Map  created - NbPixels= " << outgs.NbPixels() << endl;
   outgs.SetPixels(0.);
 
@@ -111 +87 @@
   cout << filt << endl;
 
+// FOR debug
+  if (debuglev > 0)  SpectralResponse2Nt(filt, *so, "filter");
+
   PrtTim(" After FilterCreation ");
-  }
+
+  SphereGorski<float> * ings = NULL;  // Our input map  
+  FitsIoServer fios; // Our FITS IO Server
+  char * flnm, buff[64];
+  string key;
+
+  // Loop over sky component 
+  int sk;
+  for(sk = 0; sk<nskycomp; sk++) {
+    cout << " Processing sky component No " << sk+1 << endl;
+    if (ings) { delete ings;  ings = NULL; }
+    ings = new SphereGorski<float>(hp_nside);
+    sprintf(buff, "MAPFITSFILE%d", sk+1);
+    key = buff;
+    flnm = BuildFITSFileName(dc.SParam(key, 0));
+    cout << " Reading Input FITS map " << (string)flnm << endl;
+    fios.load(*ings, flnm, 1);
+
+    if (debuglev > 4) {  // Writing tne input map to the outppf
+      FIO_SphereGorski<float> fiog(ings);
+      fiog.Write(*so, key);
+      }  
+    }
+  }   // End of try block
 
 catch (PException exc) {
   msg = exc.Msg();
   cerr << " !!!! skymixer - Catched exception - Msg= " << exc.Msg() << endl;
-  rc = 70;
+  rc = 90;
   }   
 
 
-problem:
+Fin:
+if (so) delete so;  //  Closing the debug ppf file 
 if (rc == 0)  return(0);
 cerr << " Error condition -> Rc= " << rc << endl;
@@ -128 +131 @@
 }
 
+/* Nouvelle-Fonction */
+int CheckCards(DataCards & dc, string & msg)
+//   Function to check datacards
+{
+mapPath[0] = '\0';
+hp_nside = 32;
+nskycomp = 0;
+debuglev  = 0;
+
+int rc = 0;
+  //  Cheking datacards
+  if ( (!dc.HasKey("SKYMIX")) || (!dc.HasKey("FILTER")) ) {
+     rc = 71; 
+     msg = "Invalid parameters  - NO @SKYMIX or @FILTER card ";
+     return(rc);
+  }
+
+  // Decoding number of component and pixelisation parameter
+  int mg = 32;
+  int ncomp = 0;
+  ncomp = dc.IParam("SKYMIX", 0, 0);
+  mg = dc.IParam("SKYMIX", 1, 32);
+  if (ncomp < 1)  {
+    msg = "Invalid parameters  - Check datacards @SKYMIX ";
+    rc = 72;
+    return(rc);
+  }
+
+//  Checking input FITS file specifications 
+  int kc;
+  string key, key2;
+  char buff[256];
+  bool pb = false;
+  for(kc=0; kc<ncomp; kc++) {
+    sprintf(buff, "MAPFITSFILE%d", kc+1);
+    key = buff;
+    if (dc.NbParam(key) < 1) {   
+      msg = "Missing or invalid card : " + key;
+      pb = true;  break;
+    }
+    sprintf(buff, "SPECTRAFITSFILE%d", kc+1);
+    key = buff;
+    sprintf(buff, "SPECTRAFUNC%d", kc+1);
+    key2 = buff;
+    if ( (dc.NbParam(key) < 1) && (dc.NbParam(key2) < 2) )  {   
+      msg = "Missing card or invalid parameters : " + key + " or " + key2;
+      pb = true;  break;
+    }
+
+    }
+
+  if (pb)  {
+    rc = 72;
+    return(72);
+  }
+
+//  Checking detection filter specification 
+  key = "FILTER";
+  if (dc.NbParam(key) < 3) { 
+    msg = "Missing card or invalid parameters : " + key;
+    rc = 73;  return(rc); 
+    }
+
+// Initialiazing parameters
+  rc = 0;
+  msg = "OK";
+  nskycomp = ncomp;
+  hp_nside = mg;  
+
+// Checking for PATH definition card
+  key = "MAPPATH";
+  if (dc.NbParam(key) < 3)  strncpy(mapPath, dc.SParam(key, 0).c_str(), 255);
+  mapPath[255] = '\0';
+  key = "DEBUGLEVEL";
+  debuglev =  dc.IParam(key, 0, 0);
+  return(rc);
+}
+
+static char buff_flnm[1024];   // Mal protege !
+/* Nouvelle-Fonction */
+char* BuildFITSFileName(string const & fname)
+{
+if (mapPath[0] != '\0') sprintf(buff_flnm, "%s/%s", mapPath, fname.c_str());
+else sprintf(buff_flnm, "%s", fname.c_str());
+return(buff_flnm);
+}
+
+/* Nouvelle-Fonction */
 // template <class T> 
 void addComponent(SpectralResponse&  sr, PixelMap<float>& finalMap, 
@@ -144 +235 @@
 }
 
+/* Nouvelle-Fonction */
 void addComponent(SpectralResponse&  sr, PixelMap<double>& finalMap, 
                   PixelMap<double>& mapToAdd, RadSpectra& rs, double K)
@@ -158 +250 @@
     }
 }
+
+
+/* Nouvelle-Fonction */
+void RadSpec2Nt(RadSpectra & rs, POutPersist & so, string name)
+{
+  char *ntn[2] = {"nu","fnu"};
+  NTuple nt(2,ntn);  // Creation NTuple (AVEC new )
+  float xnt[2];
+  double nu;
+  double numin = rs.minFreq();
+  double numax = rs.maxFreq();
+  int nmax = 500;
+  double dnu = (numax-numin)/nmax;
+  for(int k=0; k<nmax; k++) {
+    nu = numin+k*dnu;
+    xnt[0] = nu;
+    xnt[1] = rs.flux(nu);
+    nt.Fill(xnt);
+  }
+  ObjFileIO<NTuple> oiont(nt);
+  oiont.Write(so, name);
+  return;
+}
+
+/* Nouvelle-Fonction */
+void SpectralResponse2Nt(SpectralResponse& sr, POutPersist & so, string name)
+{
+  char *ntn[2] = {"nu","tnu"};
+  NTuple nt(2,ntn);  // Creation NTuple (AVEC new )
+  float xnt[2];
+  double nu;
+  double numin = sr.minFreq();
+  double numax = sr.maxFreq();
+  int nmax = 500;
+  double dnu = (numax-numin)/nmax;
+  for(int k=0; k<nmax; k++) {
+    nu = numin+k*dnu;
+    xnt[0] = nu;
+    xnt[1] = sr.transmission(nu);
+    nt.Fill(xnt);
+  }
+  ObjFileIO<NTuple> oiont(nt);
+  oiont.Write(so, name);
+  return;
+}
+