Index: trunk/SophyaExt/FitsIOServer/fitsioserver.cc
===================================================================
--- trunk/SophyaExt/FitsIOServer/fitsioserver.cc	(revision 471)
+++ trunk/SophyaExt/FitsIOServer/fitsioserver.cc	(revision 477)
@@ -1,12 +1,10 @@
 
-#include "nbmath.h"
-#include <math.h>
 #include <iostream.h>
+#include <list>
+
 #include "fitsioserver.h" 
-#include "dvlist.h"
 #include "strutil.h"
 
-  typedef map<string, MuTyV, less<string> >  ValList; 
-void FitsIoServer::load(TMatrix<double>& mat, char flnm[]) 
+void FitsIoServer::load(TMatrix<double>& mat,char flnm[]) 
 {
   int nbrows=0;
@@ -29,7 +27,6 @@
   nbcols=n2;
   if (naxis == 1) nbcols=1;
-  // cout << " lecture  du dvlist par load " << endl;
-  dvl.Print();
-  ValList::const_iterator it;
+  //dvl.Print();
+  DVList::ValList::const_iterator it;
   for(it = dvl.Begin(); it != dvl.End(); it++)  
     {
@@ -52,5 +49,5 @@
 	  break;
 	default :
-	  cout << " probleme dans type mot cle optionnel" << endl;
+	  cout << " FitsIOServer : probleme dans type mot cle optionnel" << endl;
 	  break;
 	}
@@ -59,8 +56,8 @@
   if (naxis > 2) 
     {
-      cout << " le fichier fits n'est pas une matrice, naxis= " << naxis << endl;
-    }
+      cout<<"  FitsIOServer : le fichier fits n'est pas une matrice, naxis= " <<naxis<< endl;
+    }
+
   // number of components
-
   if (mat.NRows() != nbrows || mat.NCols() != nbcols )
     {
@@ -75,5 +72,4 @@
   for (int j=0; j< nbcols; j++) 
     for (int i = 0; i < nbrows; i++)  mat(i,j) = dtab_[ij++];
-  // cout << " chargement termine " << endl;
 }
 
@@ -84,4 +80,5 @@
   FITS_tab_typ_ = TDOUBLE;
   read_sphere(flnm, npixels, nside);
+
   // number of pixels in the sphere
   if (sph.NbPixels() != npixels)
@@ -91,5 +88,5 @@
       if (nside <= 0 ) 
 	{
-	  cout << " FITSIOSERVER : no resolution parameter on fits file " << endl;
+	  cout<<" FITSIOSERVER: no resolution parameter on fits file "<<endl;
 	  exit(0);
 	}
@@ -97,5 +94,5 @@
 	{ 
 	  sph.Resize(nside);
-	  cout << " resize the sphere to nside=  " << nside  << endl;
+	  cout << " resizing the sphere to nside=  " << nside  << endl;
 	}
       else
@@ -106,5 +103,4 @@
     }
   for (int j = 0; j < sph.NbPixels(); j++) sph(j)= dtab_[j];
-  // cout << " chargement termine " << endl;
 }
 void FitsIoServer::load(SphericalMap<float>& sph, char flnm[]) 
@@ -121,5 +117,5 @@
       if (nside <= 0 ) 
 	{
-	  cout << " FITSIOSERVER : no resolution parameter on fits file " << endl;
+	  cout<<" FITSIOSERVER: no resolution parameter on fits file "<<endl;
 	  exit(0);
 	}
@@ -127,5 +123,5 @@
 	{ 
 	  sph.Resize(nside);
-	  cout << " resize the sphere to nside=  " << nside  << endl;
+	  cout << " resizing the sphere to nside=  " << nside  << endl;
 	}
       else
@@ -136,5 +132,4 @@
     }
   for (int j = 0; j < sph.NbPixels(); j++) sph(j)= ftab_[j];
-  //  cout << " chargement termine " << endl;
 }
 
@@ -162,9 +157,8 @@
   if (naxis != 2) 
     {
-      cout << " le fichier fits n'est pas une localmap, naxis= " << naxis << endl;
-    }
-  // cout << " lecture  du dvlist par load " << endl;
-  dvl.Print();
-  ValList::const_iterator it;
+      cout<<" FitsIOServer : le fichier fits n'est pas une localmap, naxis= "<<naxis<<endl;
+    }
+  //dvl.Print();
+  DVList::ValList::const_iterator it;
   for(it = dvl.Begin(); it != dvl.End(); it++)  
     {
@@ -187,15 +181,14 @@
 	  break;
 	default :
-	  cout << " probleme dans type mot cle optionnel" << endl;
+	  cout << " FitsIOServer : probleme dans type mot cle optionnel" << endl;
 	  break;
 	}
       if (!strcmp(keyname, "NSIDE") )
 	{
-	  // cout << " j'ai trouve NSIDE " << endl;
 	  nside = ival;
 	} 
       //if (!strcmp(keyname, "ORDERING") )
       //	{
-      //	  cout << " j'ai trouve ORDERING " << endl;
+      //  cout << " j'ai trouve ORDERING " << endl;
       //	} 
     }
@@ -209,11 +202,4 @@
   float angley=dvl.GetD("ANGLEY");
   float angle = dvl.GetD("ANGLE");
-  //cout << " theta0= " << theta0 << endl;
-  //cout << " phi0= " << phi0 << endl;
-  //cout << " x0= " << x0 << endl;
-  //cout << " y0= " << y0 << endl;
-  //cout << " angle= " << angle << endl;
-  //cout << " anglex= " << anglex << endl;
-  //cout << " angley= " << angley << endl;
 
   // number of components
@@ -226,5 +212,5 @@
       cout << "  expected " <<   lcm.Size_y()  << endl;
       lcm.ReSize(nbrows,nbcols);
-      cout << " resize the map to x-size=  " << nbrows  << " y-size= " << nbcols << endl;
+      cout << " resizing the map to x-size=  " << nbrows  << " y-size= " << nbcols << endl;
     }
   lcm.SetSize(anglex, angley);
@@ -233,19 +219,7 @@
   for (int j=0; j< nbcols; j++) 
     for (int i = 0; i < nbrows; i++)  lcm(i,j) = dtab_[ij++];
-  //cout << " chargement termine " << endl;
-}
-
-void FitsIoServer::printerror(int status) const
-{
-  //Print out cfitsio error messages and exit program 
-  if( status )
-    {
-      //print error report  
-      fits_report_error(stderr, status);
-      //terminate the program, returning error status 
-      exit( status ); 
-    }
-  return;
-}
+}
+
+
 void FitsIoServer::save( TMatrix<double>& mat, char filename[])   
 {
@@ -253,5 +227,5 @@
   int nbcols = mat.NCols();
   long naxis = nbcols > 1 ? 2 : 1; 
-  //cout << " nombre de lignes : " << nbrows << " colonnes " << nbcols << endl;
+  cout << " nombre de lignes : " << nbrows << " colonnes " << nbcols << endl;
   FITS_tab_typ_ = TDOUBLE;
   if (dtab_ != NULL ) delete[] dtab_;
@@ -278,7 +252,5 @@
   DVList dvl;
 
-  //  cout << " save  : debut ecriture " << endl;
   planck_write_img(fptr, naxis, nbrows, nbcols, 0, dvl);
-  //  cout << " save  : fin ecriture " << endl;
 
 
@@ -326,9 +298,4 @@
   if( status ) printerror( status );
 
-  char comment[] = "Generated by synfast";
-  char svalue [] = "SYNFAST ";
-  strcat(svalue, version);
-  fits_write_key(fptr, TSTRING, "ORIG-MAP", svalue, comment, &status);
-  if( status )  printerror( status );
 
   strcpy(comment,"HEALPIX Pixelisation");
@@ -343,8 +310,4 @@
   if( status )  printerror( status );
 
-  strcpy(comment,"Maximum multipole l used in map synthesis");
-  int lmax= sph.nlmax();
-  fits_write_key(fptr, TINT, "MAX-LPOL", &lmax, comment, &status);
-  if( status )  printerror( status );
 
   strcpy(comment,"Random generator seed");
@@ -353,8 +316,4 @@
   if( status )  printerror( status );
 
-  strcpy(comment,"FWHM in DEGREES");
-  float fwhm_deg= sph.fwhm()/60.0;
-  fits_write_key(fptr, TFLOAT, "FWHM", &fwhm_deg, comment, &status);
-  if( status )  printerror( status );
 
   strcpy(comment,"--------------------------------------------");
@@ -370,28 +329,4 @@
   if( status ) printerror( status );
   
-  //   write information about the power spectrum
-  char card[]="";
-  sph.powfile(card);
-  strcpy(comment,"Input power spectrum in : ");
-  strcat(comment, card);
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
-
-  strcpy(comment,"Quadrupole :");
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
-
-  strcpy(comment,"C(2)= ");
-  sprintf(card,"%g",sph.quadrupole());
-  strcat(comment, card);
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
-
-  strcpy(comment,"keywords relative to input power spectrum ");
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
-  strcpy(comment,"have not yet been defined for PLANCK Surveyor Project");
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
   */
 
@@ -467,13 +402,13 @@
   dvl["NSIDE"] = locm.SizeIndex();
   dvl["ORDERING"]=locm.TypeOfMap();
-  float theta0;
-  float phi0;
+  double theta0;
+  double phi0;
   int x0;
   int y0;
-  float angle;
-  locm.Origin(theta0, phi0,x0, y0, angle);
-  float anglex;
-  float angley;
-  locm.Aperture(anglex, angley);
+  double angle;
+  locm.Origin(theta0,phi0,x0,y0,angle);
+  double anglex;
+  double angley;
+  locm.Aperture(anglex,angley);
   dvl["THETA0"] = theta0;
   dvl["PHI0"] = phi0;
@@ -491,6 +426,4 @@
 
 }
-
-
 void FitsIoServer::planck_write_img( fitsfile *fptr, int naxis,int n1, int n2, int n3, DVList& dvl)
 {
@@ -500,4 +433,5 @@
   if ( FITS_tab_typ_ == TDOUBLE) bitpix = DOUBLE_IMG;
   if ( FITS_tab_typ_ == TFLOAT)  bitpix = FLOAT_IMG;
+  if ( FITS_tab_typ_ == TINT)   bitpix = LONG_IMG;
   long naxes[3];
   naxes[0] = n1;
@@ -524,6 +458,10 @@
       if( status ) printerror( status );
       break;
+    case TINT :
+      fits_write_img(fptr, TINT, 1, nelements, i_4tab_,  &status); 
+      if( status ) printerror( status );
+      break;
     default :
-      cout << " fitsioserver : type de tableau non traite en ecriture " << endl;
+      cout << " FitsIOServer : type de tableau non traite en ecriture " << endl;
       break;
     }
@@ -535,6 +473,6 @@
   // la date servira de repere pour relire ces mots cles.
 
-  dvl.Print();
-  ValList::const_iterator it;
+  //dvl.Print();
+  DVList::ValList::const_iterator it;
   for(it = dvl.Begin(); it != dvl.End(); it++)  
     {
@@ -568,5 +506,5 @@
 	  break;
 	default :
-	  cout << " probleme dans type mot cle optionnel" << endl;
+	  cout << " FitsIOServer : probleme dans type mot cle optionnel" << endl;
 	  break;
 	}
@@ -574,4 +512,5 @@
     }
 }
+
 void FitsIoServer::planck_read_img( fitsfile *fptr, long& naxis,int& n1, int& n2, int& n3, DVList& dvl)
 {
@@ -607,4 +546,9 @@
     {
     case TDOUBLE :
+      if (bitpix !=  DOUBLE_IMG) 
+	{
+	  cout << " FitsIOServer : the data type on fits file is not double, "
+	       << " conversion to double will be achieved by cfitsio lib " << endl; 
+	}
       if (dtab_ != NULL) delete [] dtab_;
       dtab_=new double[nelements];
@@ -614,4 +558,9 @@
       break;
     case TFLOAT :
+      if (bitpix !=  FLOAT_IMG) 
+	{
+	  cout << " FitsIOServer : the data type on fits file is not float, "
+	       << " conversion to float will be achieved by cfitsio lib " << endl; 
+	}
       if (ftab_ != NULL) delete [] ftab_;
       ftab_=new float[nelements];
@@ -620,6 +569,22 @@
       if( status )  printerror( status );
       break;
+
+
+    case TINT :
+      if (bitpix !=  LONG_IMG) 
+	{
+	  cout << " FitsIOServer : the data type on fits file is not long, "
+	       << " conversion to long will be achieved by cfitsio lib " << endl; 
+	}
+      if (i_4tab_ != NULL) delete [] i_4tab_;
+      i_4tab_=new int_4[nelements];
+      fits_read_img(fptr, TINT, 1, nelements,  &fnullval, i_4tab_, 
+		    &anynull, &status);
+      if( status )  printerror( status );
+      break;
+
+
     default :
-      cout << " fitsio::read_img : type non traite: " << FITS_tab_typ_  << endl;
+      cout << " FitsIOServer::read_img : type non traite: " << FITS_tab_typ_  << endl;
       break;
     }
@@ -806,6 +771,6 @@
 
   npixels=n1;
-  dvl.Print();
-  ValList::const_iterator it;
+  //dvl.Print();
+  DVList::ValList::const_iterator it;
   for(it = dvl.Begin(); it != dvl.End(); it++)  
     {
@@ -828,14 +793,13 @@
 	  break;
 	default :
-	  cout << " probleme dans type mot cle optionnel" << endl;
+	  cout << " FitsIOServer : probleme dans type mot cle optionnel" << endl;
 	  break;
 	}
       if (!strcmp(keyname, "NSIDE") )
 	{
-	  //  cout << " j'ai trouve NSIDE " << endl;
 	  nside = ival;
 	} 
-      // if (!strcmp(keyname, "ORDERING") )
-      ///	{
+      //if (!strcmp(keyname, "ORDERING") )
+      //	{
       //	  cout << " j'ai trouve ORDERING " << endl;
       //	} 
@@ -854,5 +818,5 @@
   fitsfile *fptr;  
   int bitpix = FLOAT_IMG;
-  long naxis = 2;   
+  long naxis = 2;
 
   // delete old file if it already exists
@@ -889,3 +853,623 @@
 }
 
-
+void FitsIoServer::save(NTuple& ntpl,char flnm[])
+
+  //****************************************************/
+  //* read the elements of the NTuple ntpl, and create */ 
+  //* a FITS file with a binary table extension        */
+  //****************************************************/
+{
+
+  // delete old file if it already exists
+  remove(flnm);
+
+  // pointer to the FITS file, defined in fitsio.h 
+  fitsfile *fptr;
+  int status = 0;
+  if( fits_create_file(&fptr, flnm, &status) )
+    printerror( status );
+
+  // table will have tfields columns
+  int tfields= ntpl.NVar();
+
+  // table will have nrows rows
+  int nrows= ntpl.NEntry();
+
+  // extension name
+  char extname[] = "NTuple_Binary_tbl";
+
+  // define the name, and the datatype for the tfields columns
+  char **ttype, **tform;
+  ttype= new char*[tfields];
+  tform= new char*[tfields];
+  for(int i = 0; i < tfields; i++)
+    {
+      ttype[i]= new char[FLEN_VALUE];
+      strcpy(ttype[i], ntpl.NomIndex(i));
+      tform[i]= new char[FLEN_VALUE];
+      strcpy(tform[i], "1E");
+    }
+
+  // create a new empty binary table onto the FITS file
+  // physical units if they exist, are defined in the DVList object
+  // so the null pointer is given for the tunit parameters.
+  if ( fits_create_tbl(fptr,BINARY_TBL,nrows,tfields,ttype,tform,
+		       NULL,extname,&status) )
+    printerror( status );
+
+  // first row in table to write
+  int firstrow  = 1; 
+
+  // first element in row  (ignored in ASCII tables)
+  int firstelem = 1;
+
+  for(int i = 0; i < tfields; i++)
+    {
+      float *dens= new float[nrows];
+      for(int j = 0; j < nrows; j++)
+	{
+	  dens[j]= ntpl.GetVal(j,i);
+	}
+      fits_write_col(fptr,TFLOAT,i+1,firstrow,firstelem,nrows,dens,&status);
+      delete []dens;
+    }
+
+  // number of blocks per event
+  int blk= ntpl.BLock();
+  fits_write_key(fptr,TINT,"BLK",&blk,"number of blocks per evt",&status);
+
+  // get names and values from the join DVList object
+  DVList dvl= ntpl.Info();
+  dvl.Print();
+  DVList::ValList::const_iterator it;
+  for(it = dvl.Begin(); it != dvl.End(); it++)
+    {
+      char keytype= (*it).second.typ;
+      char keyname[9]="";
+      strncpy(keyname,(*it).first.substr(0,64).c_str(),8);
+      char comment[FLEN_COMMENT]="";
+
+      switch (keytype)
+	{
+	case 'I' :
+	  {
+	    int ival=(*it).second.mtv.iv;
+	    strcpy(comment,"I entier");
+	    fits_write_key(fptr,TINT,keyname,&ival,comment,&status);
+	    break;
+	  }
+	case 'D' :
+	  {
+	    double dval=(*it).second.mtv.dv;
+	    strcpy(comment,"D double");
+	    fits_write_key(fptr,TDOUBLE,keyname,&dval,comment,&status);
+	    break;
+	  }
+	case 'S' :
+	  {
+	    char strval[]="";
+	    strcpy(strval,(*it).second.mtv.strv);
+	    strcpy(comment,"S character string");
+	    fits_write_key(fptr,TSTRING,keyname,&strval,comment,&status);
+	    break;
+	  }
+	} 
+    }
+
+  //close the FITS file
+  if ( fits_close_file(fptr, &status) )
+    printerror( status );
+  return;
+}  
+
+void FitsIoServer::load(NTuple& ntpl,char flnm[],int hdunum) 
+
+  //********************************************************/
+  //* move to the HDU which has the specified number hdunum*/
+  //* in the FITS file, perform read operations and write  */
+  //* the elements in an NTuple.                           */
+  //********************************************************/
+{
+
+  // pointer to the FITS file, defined in fitsio.h 
+  fitsfile *fptr;
+  int status = 0;
+  if( fits_open_file(&fptr,flnm,READONLY,&status) ) 
+    printerror( status );
+
+  // move to the HDU
+  int hdutype;
+  if( fits_movabs_hdu(fptr,hdunum,&hdutype,&status) ) 
+    printerror( status );
+
+  // get number of keywords
+  int nkeys,keypos;
+  if( fits_get_hdrpos(fptr,&nkeys,&keypos,&status) )
+    printerror( status );
+
+  /*
+  printf("Move to extensions by name and version number: (ffmnhd)\n");
+  char* binname= "Test-BINTABLE";
+  int extvers = 1;
+  fits_movnam_hdu(fptr,ANY_HDU, binname,extvers,&status);
+  fits_get_hdu_num(fptr,&hdunum);
+  printf(" %s, %ld = hdu %d, %d\n", binname, extvers, hdunum, status);
+  */
+
+  if (hdutype == BINARY_TBL)
+    printf("\nReading binary table in HDU %d:\n",hdunum);
+  else
+    {
+      printf("Error:: this HDU is not a binary table\n");
+      exit ( status );
+    }
+
+  // number of columns
+  int tfields;
+  if( fits_get_num_cols(fptr,&tfields,&status) )
+    printerror( status );
+
+  // to get table size
+  long naxis[2];
+  int nfound;
+   if( fits_read_keys_lng(fptr, "NAXIS", 1, 2, naxis, &nfound, &status) ) 
+     printerror( status );
+   int nrows= naxis[1];
+
+   printf("\nInformation about each column:\n");
+
+   //Information about each column
+   char **ttype, **tform;
+   ttype= new char*[tfields];
+   tform= new char*[tfields];
+   for( int ii = 0; ii < tfields; ii++)
+     {
+       ttype[ii]= new char[FLEN_VALUE];
+       tform[ii]= new char[FLEN_VALUE];
+     }
+
+   // number of TTYPEn,TFORMn and TUNITn keywords found
+   int num= 0;
+
+   // read the column names from the TTYPEn keywords
+   fits_read_keys_str(fptr,"TTYPE",1,tfields,ttype,&nfound,&status);
+   num += nfound;
+   // read the column names from the TFORMn keywords
+   fits_read_keys_str(fptr,"TFORM",1,tfields,tform,&nfound,&status);
+   num += nfound;
+
+   for(int ii = 0; ii < tfields; ii++)
+     printf("\nColumn name & Format %8s %8s", ttype[ii],tform[ii]);
+  
+   // select the columns with float data values (typecode = 42)
+   int typecode;
+   long repeat,width;
+   list<int> column;
+   for( int ii = 0; ii < tfields; ii++)
+     {
+       fits_binary_tform(tform[ii],&typecode, &repeat, &width, &status);
+       //printf("\n%4s %3d %2ld %2ld", tform[ii], typecode, repeat, width);
+       if(typecode == 42) 
+	 column.push_back(ii+1);
+     }
+
+   // get input elements to create the NTuple
+   char **clname;
+   clname= new char*[column.size()];
+   for(int ii = 0; ii < column.size(); ii++)
+     clname[ii]= new char[FLEN_VALUE];
+
+   list<int>::iterator itr;
+   int index= 0;
+   for(itr= column.begin(); itr != column.end(); itr++)
+     strcpy(clname[index++],ttype[*itr-1]);
+
+   // check if the specified keyword BLK exists
+   int blk= 512;
+   if(check_keyword(fptr,nkeys,"BLK"))
+     {
+       if( fits_read_key(fptr,TINT,"BLK",&blk,NULL,&status) )
+	 printerror( status );
+     }
+  
+   // create a NTuple
+   NTuple nt0(column.size(),clname,blk);
+
+   float value[1];
+   long felem = 1;
+   char strnull[10];
+   strcpy(strnull, " ");
+   int anynull= 0;
+   long longnull  = 0;
+
+   // value to represent undefined array elements
+   float floatnull = FLOATNULLVALUE;
+
+   // read elements from columns and fill NTuple
+   for (int k = 0; k < nrows; k++)
+     {
+       int j= 0;
+       float *xnt= new float[column.size()];
+       list<int>::iterator itr;
+       for(itr= column.begin(); itr != column.end(); itr++)
+	 {
+	   fits_read_col(fptr,TFLOAT,*itr,k+1,felem,1,&floatnull,value,
+			 &anynull, &status);
+	   xnt[j++]= value[0];
+	 }
+       nt0.Fill(xnt);
+       delete[] xnt;
+     }
+
+   // the TUNITn keywords are optional, if they exist they are put 
+   // in the DVLIst object
+   char keyname[LEN_KEYWORD]= "";
+   char strval[FLEN_VALUE];
+   for(int ii = 0; ii < tfields; ii++)
+     {
+       fits_make_keyn("TUNIT",ii+1,keyname,&status);
+       if(check_keyword(fptr,nkeys,keyname))
+	 {
+	   num++;
+	   if( fits_read_key_str(fptr,keyname,strval,NULL,&status) )
+	     printerror(status);
+	   strip (keyname, 'B',' '); 
+	   strip(strval, 'B',' ');
+	   nt0.Info()[keyname]= strval;
+	 }
+     }
+
+   // add the number of mandatory keywords of a binary table
+   num += 8;
+
+   // put names and values of other reserved keywords in a DVList object
+   char comment[FLEN_COMMENT];
+   char dtype;
+   for(int j = num+1; j <= nkeys; j++)
+     {
+       fits_read_keyn(fptr,j,keyname,strval,comment,&status);
+       fits_get_keytype(strval,&dtype,&status);
+       strip (keyname, 'B',' '); 
+       strip(strval, 'B',' ');
+       
+       switch( dtype ) 
+	 {
+	 case 'C':
+	   nt0.Info()[keyname]= strval;	  
+	   break;
+	 case 'I':
+	   int ival;
+	   ctoi(strval,&ival);
+	   nt0.Info()[keyname]= ival;
+	   break;
+	 case 'F':
+	   double dval;
+	   ctof(strval,&dval);
+	   nt0.Info()[keyname]= dval;
+	   break;
+	 }
+     }
+
+   // copy in the input NTuple
+   ntpl= nt0;
+
+   if( fits_close_file(fptr, &status) )
+     printerror(status);
+
+   printf("\n");
+   return;
+}
+
+bool FitsIoServer::check_keyword(fitsfile *fptr,int nkeys,char keyword[])
+
+  //*****************************************************/
+  //* check if the specified keyword exits in the CHU   */
+  //*****************************************************/
+{
+
+  bool KEY_EXIST = false;
+  int status = 0;
+  char strbide[FLEN_VALUE];
+  char keybide[LEN_KEYWORD]= "";
+  for(int jj = 1; jj <= nkeys; jj++)
+    {
+      if( fits_read_keyn(fptr,jj,keybide,strbide,NULL,&status) )
+	printerror( status );
+      if( !strcmp(keybide,keyword) )
+	{ 
+	  KEY_EXIST= true;
+	  break;
+	}
+    }
+  return(KEY_EXIST);
+}
+
+void FitsIoServer::readheader ( char filename[] )
+
+  //**********************************************************************/
+  //* Print out all the header keywords in all extensions of a FITS file */
+  //**********************************************************************/
+{
+
+  // standard string lengths defined in fitsioc.h
+  char card[FLEN_CARD];
+
+  // pointer to the FITS file, defined in fitsio.h
+  fitsfile *fptr;
+
+  int status = 0;
+  if ( fits_open_file(&fptr, filename, READONLY, &status) ) 
+    printerror( status );
+  
+  // attempt to move to next HDU, until we get an EOF error
+  int hdutype;
+  for (int ii = 1; !(fits_movabs_hdu(fptr,ii,&hdutype,&status));ii++) 
+    {
+      if (hdutype == ASCII_TBL)
+	printf("\nReading ASCII table in HDU %d:\n", ii);
+      else if (hdutype == BINARY_TBL)
+	printf("\nReading binary table in HDU %d:\n", ii);
+      else if (hdutype == IMAGE_HDU)
+	printf("\nReading FITS image in HDU %d:\n", ii);
+      else
+	{
+          printf("Error: unknown type of this HDU \n");
+          printerror( status );
+	}
+
+      // get the number of keywords
+      int nkeys, keypos;
+      if ( fits_get_hdrpos(fptr, &nkeys, &keypos, &status) )
+	printerror( status );
+      
+      printf("Header listing for HDU #%d:\n", ii);
+      for (int jj = 1; jj <= nkeys; jj++)  
+	{
+	  if ( fits_read_record(fptr, jj, card, &status) )
+	    printerror( status );
+
+	  // print the keyword card
+	  printf("%s\n", card); 
+        }
+      printf("END\n\n");
+    }
+
+  // got the expected EOF error; reset = 0
+  if (status == END_OF_FILE)   
+    status = 0;             
+  else
+    printerror( status );
+
+  if ( fits_close_file(fptr, &status) )
+    printerror( status );
+  
+  return;
+}
+
+void FitsIoServer::printerror(int status) const
+
+  //*****************************************************/
+  //* Print out cfitsio error messages and exit program */
+  //*****************************************************/
+{
+
+  // print out cfitsio error messages and exit program 
+  if( status )
+    {
+      // print error report  
+      fits_report_error(stderr, status);
+      // terminate the program, returning error status 
+      exit( status ); 
+    }
+  return;
+}
+
+void FitsIoServer::save(ImageR4& DpcImg,char flnm[])
+{
+  long naxis=2;
+  int siz_x = DpcImg.XSize();
+  int siz_y = DpcImg.YSize();
+  int nbpixels = siz_x*siz_y;
+  FITS_tab_typ_ = TFLOAT;
+
+  // pointer to the FITS file, defined in fitsio.h
+  fitsfile *fptr;
+
+  // initialize status before calling fitsio routines
+  int status = 0;         
+
+  // delete old file if it already exists
+  remove(flnm);
+
+  // create new FITS file  
+  fits_create_file(&fptr, flnm, &status); 
+  if( status )  printerror( status );
+
+  // get the values of the DpcImage
+  if (ftab_ != NULL) delete [] ftab_;
+  ftab_=DpcImg.ImagePtr();
+
+  // write FITS image
+  DVList dvl;
+
+  dvl["DATATYPE"] = int(DpcImg.PixelType());
+  dvl["ORG_X"] = DpcImg.XOrg();
+  dvl["ORG_Y"] = DpcImg.YOrg();
+  dvl["PIXSZ_X"] = DpcImg.XPxSize();
+  dvl["PIXSZ_Y"] = DpcImg.YPxSize();
+  dvl["IDENT"] = DpcImg.Ident();
+  dvl["NAME"] = DpcImg.Nom();
+  dvl["NBSAT"] = DpcImg.nbSat;
+  dvl["NBNUL"] = DpcImg.nbNul;
+  dvl["MINPIX"] = DpcImg.minPix;
+  dvl["MAXPIX"] = DpcImg.maxPix;
+  dvl["MOYPIX"] = DpcImg.moyPix;
+  dvl["SIGPIX"] = DpcImg.sigPix;
+  dvl["FOND"] = DpcImg.fond;
+  dvl["SIGFON"] = DpcImg.sigmaFond;
+
+
+  planck_write_img(fptr, naxis, siz_x, siz_y, 0, dvl);
+
+  // close the file
+  fits_close_file(fptr, &status);
+  if( status )  printerror( status );
+
+  ftab_=NULL;
+}
+
+void FitsIoServer::load(ImageR4& DpcImg,char flnm[])
+{
+  FITS_tab_typ_ = TFLOAT;
+  long naxis;
+  int siz_x;
+  int siz_y;
+  int n3;
+  DVList dvl;
+
+  // pointer to the FITS file, defined in fitsio.h
+  fitsfile *fptr;
+  // initialize status before calling fitsio routines
+  int status = 0; 
+        
+  fits_open_file(&fptr, flnm, READONLY, &status);
+  if( status )  printerror( status );
+
+  planck_read_img(fptr, naxis, siz_x, siz_y, n3, dvl);
+
+  // close the file
+  fits_close_file(fptr, &status);
+  if( status )  printerror( status );
+  
+  if (naxis != 2) 
+    {
+      cout << " FitsIOServer : naxis= " << naxis << " not equal to 2 ?" << endl;
+    }
+
+  DpcImg.Allocate(siz_x, siz_y, 0., 0);
+  float* pixelPtr=DpcImg.ImagePtr();
+  //for (int i=0; i < siz_x*siz_y; i++) pixelPtr[i]=ftab_[i];
+  // verifications de type
+  PBaseDataTypes  dataT=DataType((r_4)0);
+  int TypeDonnees=dvl.GetI("DATATYPE");
+  if (int(dataT) !=  TypeDonnees)
+    {
+      cout << " FitsIOServer : parameter DATATYPE on file is not float " << endl;
+      cout << "  eventual conversion to float is achieved by cfitsio lib " << endl;
+    }
+
+  memcpy(pixelPtr, ftab_, siz_x*siz_y*DataSize(dataT));
+  const char* nom=dvl.GetS("NAME").c_str();
+  DpcImg.SetNameId(dvl.GetI("IDENT"), nom);
+  DpcImg.SetOrg(dvl.GetI("ORG_X"), dvl.GetI("ORG_Y"));
+  DpcImg.SetPxSize(dvl.GetD("PIXSZ_X"), dvl.GetD("PIXSZ_Y"));    
+  DpcImg.SetAtt(dvl.GetI("NBNUL"), dvl.GetI("NBSAT"),
+       dvl.GetD("MINPIX"), dvl.GetD("MAXPIX"), dvl.GetD("MOYPIX"), 
+       dvl.GetD("SIGPIX"),dvl.GetD("FOND"), dvl.GetD("SIGFON"));
+
+}
+
+void FitsIoServer::save(ImageI4& DpcImg,char flnm[])
+{
+  long naxis=2;
+  int siz_x = DpcImg.XSize();
+  int siz_y = DpcImg.YSize();
+  int nbpixels = siz_x*siz_y;
+  FITS_tab_typ_ = TINT;
+
+  // pointer to the FITS file, defined in fitsio.h
+  fitsfile *fptr;
+
+  // initialize status before calling fitsio routines
+  int status = 0;         
+
+  // delete old file if it already exists
+  remove(flnm);
+
+  // create new FITS file  
+  fits_create_file(&fptr, flnm, &status); 
+  if( status )  printerror( status );
+
+  // get the values of the DpcImage
+  if (i_4tab_ != NULL) delete [] i_4tab_;
+  i_4tab_=DpcImg.ImagePtr();
+
+  // write FITS image
+  DVList dvl;
+
+  dvl["DATATYPE"] = int(DpcImg.PixelType());
+  dvl["ORG_X"] = DpcImg.XOrg();
+  dvl["ORG_Y"] = DpcImg.YOrg();
+  dvl["PIXSZ_X"] = DpcImg.XPxSize();
+  dvl["PIXSZ_Y"] = DpcImg.YPxSize();
+  dvl["IDENT"] = DpcImg.Ident();
+  dvl["NAME"] = DpcImg.Nom();
+  dvl["NBSAT"] = DpcImg.nbSat;
+  dvl["NBNUL"] = DpcImg.nbNul;
+  dvl["MINPIX"] = DpcImg.minPix;
+  dvl["MAXPIX"] = DpcImg.maxPix;
+  dvl["MOYPIX"] = DpcImg.moyPix;
+  dvl["SIGPIX"] = DpcImg.sigPix;
+  dvl["FOND"] = DpcImg.fond;
+  dvl["SIGFON"] = DpcImg.sigmaFond;
+
+
+  planck_write_img(fptr, naxis, siz_x, siz_y, 0, dvl);
+
+  // close the file
+  fits_close_file(fptr, &status);
+  if( status )  printerror( status );
+
+  i_4tab_=NULL;
+}
+
+void FitsIoServer::load(ImageI4& DpcImg,char flnm[])
+{
+  FITS_tab_typ_ = TINT;
+  long naxis;
+  int siz_x;
+  int siz_y;
+  int n3;
+  DVList dvl;
+
+  // pointer to the FITS file, defined in fitsio.h
+  fitsfile *fptr;
+  // initialize status before calling fitsio routines
+  int status = 0; 
+        
+  fits_open_file(&fptr, flnm, READONLY, &status);
+  if( status )  printerror( status );
+  planck_read_img(fptr, naxis, siz_x, siz_y, n3, dvl);
+
+  // close the file
+  fits_close_file(fptr, &status);
+  if( status )  printerror( status );
+  
+  if (naxis != 2) 
+    {
+      cout << " FitsIOServer : naxis= " << naxis << " not equal to 2 ?" << endl;
+    }
+
+  DpcImg.Allocate(siz_x, siz_y, 0., 0);
+  int_4* pixelPtr=DpcImg.ImagePtr();
+
+
+  // verifications de type
+  PBaseDataTypes  dataT=DataType((int_4)0);
+  int TypeDonnees=dvl.GetI("DATATYPE");
+  if (int(dataT) !=  TypeDonnees)
+    {
+      cout << " FitsIOServer : parameter DATATYPE on file is not int_4 " << endl;
+      cout << "  eventual conversion to float is achieved by cfitsio lib " << endl;
+    }
+  memcpy(pixelPtr, i_4tab_, siz_x*siz_y*DataSize(dataT));
+  const char* nom=dvl.GetS("NAME").c_str();
+  DpcImg.SetNameId(dvl.GetI("IDENT"), nom);
+  DpcImg.SetOrg(dvl.GetI("ORG_X"), dvl.GetI("ORG_Y"));
+  DpcImg.SetPxSize(dvl.GetD("PIXSZ_X"), dvl.GetD("PIXSZ_Y"));    
+  DpcImg.SetAtt(dvl.GetI("NBNUL"), dvl.GetI("NBSAT"),
+       dvl.GetD("MINPIX"), dvl.GetD("MAXPIX"), dvl.GetD("MOYPIX"), 
+       dvl.GetD("SIGPIX"),dvl.GetD("FOND"), dvl.GetD("SIGFON"));
+}
+
+
Index: trunk/SophyaExt/FitsIOServer/fitsioserver.h
===================================================================
--- trunk/SophyaExt/FitsIOServer/fitsioserver.h	(revision 471)
+++ trunk/SophyaExt/FitsIOServer/fitsioserver.h	(revision 477)
@@ -6,55 +6,83 @@
 #include "sphericalmap.h"
 #include "localmap.h"
+#include "ntuple.h"
+#include "outilsinit.h"
+#include "cimage.h"
 #include "dvlist.h"
+
+#define LEN_KEYWORD 9
 
 class FitsIoServer 
 {
+
 public :
-FitsIoServer() : FITS_tab_typ_(0), itab_(NULL), ftab_(NULL), dtab_(NULL) {;}
-~FitsIoServer()
-  {
-    if (itab_ != NULL) delete [] itab_;
-    if (ftab_ != NULL) delete [] ftab_;
-    if (dtab_ != NULL) delete [] dtab_;
-  }
-void load(SphericalMap<double>& sph, char flnm[]);
-void load(SphericalMap<float>& sph, char flnm[]);
-void load(TMatrix<double>& mat, char flnm[]);
-void load(LocalMap<double>& sph, char flnm[]);
-void save(SphericalMap<double>& sph, char flnm[]);
-void save(SphericalMap<float>& sph, char flnm[]);
-void save(LocalMap<double>& locm, char flnm[]);
-void save(TMatrix<double>& mat, char flnm[]);
-void sinus_picture_projection(SphericalMap<double>& sph, char flnm[]);
-void sinus_picture_projection(SphericalMap<float>& sph, char flnm[]);
-void picture(LocalMap<double>& lcm, char flnm[]);
-private :
-void read_sphere(char flnm[], int& npixels, int& nside);
-void planck_write_img( fitsfile *fptr, int naxis,int n1, int n2, int n3, DVList& dvl);
-void planck_read_img( fitsfile *fptr, long& naxis,int& n1, int& n2, int& n3, DVList& dvl);
-void  write_picture(long* naxes, float* map,  char* filename) const;
 
-void printerror(int status) const;
-inline int key_type_PL2FITS(char typ) 
-  {
-  switch (typ) 
-    {
-    case 'I' :
-      return TINT;
-    case 'D' : 
-      return TDOUBLE;
-    case 'S' : 
-      return TSTRING;
-    default :
-      cout << " FITSIOSERVER : type de donnee non prevu " << endl;
-      return 0;
-    }
-    
-  }
-int FITS_tab_typ_;
-int* itab_;
-float* ftab_;
-double* dtab_;
+  FitsIoServer() : FITS_tab_typ_(0),i_4tab_(NULL),ftab_(NULL),dtab_(NULL) {;}
+ ~FitsIoServer()
+   {
+     if (i_4tab_ != NULL) delete [] i_4tab_;
+     if (ftab_ != NULL) delete [] ftab_;
+     if (dtab_ != NULL) delete [] dtab_;
+   }
+
+ void load(SphericalMap<double>& sph, char flnm[]);
+ void load(SphericalMap<float>& sph, char flnm[]);
+ void load(TMatrix<double>& mat, char flnm[]);
+ void load(LocalMap<double>& sph, char flnm[]);
+ void save(SphericalMap<double>& sph, char flnm[]);
+ void save(SphericalMap<float>& sph, char flnm[]);
+ void save(LocalMap<double>& locm, char flnm[]);
+ void save(TMatrix<double>& mat, char flnm[]);
+ void sinus_picture_projection(SphericalMap<double>& sph, char flnm[]);
+ void sinus_picture_projection(SphericalMap<float>& sph, char flnm[]);
+ void picture(LocalMap<double>& lcm, char flnm[]);
+
+ void load(NTuple& ntpl,char flnm[],int hdunum);
+ void save(NTuple& ntpl,char flnm[]);
+
+ void load(ImageR4& DpcImg,char flnm[]);
+ void save(ImageR4& DpcImg,char flnm[]);
+ void load(ImageI4& DpcImg,char flnm[]);
+ void save(ImageI4& DpcImg,char flnm[]);
+
+ void readheader(char flnm[]);
+ 
+ private :
+  
+ void read_sphere(char flnm[],int& npixels,int& nside);
+
+ // creer, ecrire une imageformat FITS, a partir des tableaux de donnees 
+ // dtab_, ftab_ prealablement remplis 
+ void planck_write_img(fitsfile *fptr,int naxis,int n1,int n2,int n3,DVList& dvl);
+
+ void planck_read_img(fitsfile *fptr,long& naxis,int& n1,int& n2,int& n3,DVList& dvl);
+ void write_picture(long* naxes,float* map,char* filename) const;
+
+ void printerror(int status) const;
+
+ bool check_keyword(fitsfile*,int,char keyword[]);
+
+ inline int key_type_PL2FITS(char typ) 
+   {
+     switch (typ) 
+       {
+       case 'I' :
+	 return TINT;
+       case 'D' : 
+	 return TDOUBLE;
+       case 'S' : 
+	 return TSTRING;
+       default :
+	 cout << " FITSIOSERVER : type de donnee non prevu " << endl;
+	 return 0;
+       }
+   }
+
+ int FITS_tab_typ_;
+ int* i_4tab_;
+ float* ftab_;
+ double* dtab_;
 };
+
 #endif
 
