Index: trunk/SophyaExt/FitsIOServer/fitsioserver.cc
===================================================================
--- trunk/SophyaExt/FitsIOServer/fitsioserver.cc	(revision 837)
+++ trunk/SophyaExt/FitsIOServer/fitsioserver.cc	(revision 838)
@@ -10,10 +10,10 @@
 //************************************************************************
 
-#include "machdefs.h" 
+#include "machdefs.h"   
 #include <iostream.h>
 
 // #include <sstream>  Ca ne passe pas pour le moment avec g++ - Reza 23/12/99
 
-#include <list>
+//#include <list>
 #include <string>
 
@@ -80,5 +80,5 @@
 	  if(DTYPE == TDOUBLE)
 	    {
-	      SphereGorski<double> *sph= new SphereGorski<double>(nside);
+	      SphereHealpix<double> *sph= new SphereHealpix<double>(nside);
 	      load(*sph,flnm,hdunum);
 	      return sph;
@@ -87,5 +87,5 @@
 	  else if(DTYPE == TFLOAT)
 	    {
-	      SphereGorski<float> *sph= new SphereGorski<float>(nside);
+	      SphereHealpix<float> *sph= new SphereHealpix<float>(nside);
 	      load(*sph,flnm,hdunum);
 	      return sph;
@@ -313,4 +313,5 @@
     }
 
+
   // get the number of columns
   int ncols= 0;
@@ -328,6 +329,5 @@
   long repeat;
   int ntest= 0;
-  int ii;
-  for(ii = 0; ii < ncols; ii++)
+  for(int ii = 0; ii < ncols; ii++)
     {
       int DTYPE;
@@ -354,9 +354,9 @@
   char **clname;
   clname= new char*[tfields.size()];
-  for(ii = 0; ii < tfields.size(); ii++)
+  for(int ii = 0; ii < tfields.size(); ii++)
     {
       //  ostringstream ne passe pas sur tous les compilos .. Reza 23/12/99
       //      ostringstream oss; oss<<"Col_"<<ii+1;
-      //      string ss= oss.str();
+      //      string ss= oss.str(); 
       //      clname[ii]= new char[FLEN_VALUE];
       clname[ii] = new char[16];
@@ -367,5 +367,5 @@
   // create a NTuple
   NTuple nt0(tfields.size(),clname);
-  for(ii = 0; ii < tfields.size(); ii++)
+  for(int ii = 0; ii < tfields.size(); ii++)
     delete [] clname[ii];
   delete [] clname;
@@ -396,4 +396,6 @@
     } 
 
+
+
   // get number of keywords
   int nkeys,keypos;
@@ -411,4 +413,7 @@
   int num= 8;
 
+
+
+
   for(int j = num+1; j <= nkeys; j++)
     {
@@ -427,4 +432,5 @@
 	  strip(strval, 'B',' ');
 	  strip(strval, 'B','\'');
+	  
 	  switch( dtype )  
 	    { 
@@ -435,5 +441,5 @@
 	      int ival;
 	      fits_read_key(fptr,TINT,keyname,&ival,NULL,&status);
-	      nt0.Info()[keyname]= (int_4) ival;		// Portage mac DY
+	        nt0.Info()[keyname]= (int_4) ival;		// Portage mac DY
 	      break;
 	    case 'L':
@@ -449,7 +455,10 @@
 	      break;
 	    }
+	  
 	}
     }
   
+
+
   // copy in the input NTuple
   ntpl= nt0;
@@ -492,5 +501,5 @@
 	{ 
 	  sph.Resize(nside);
-	  cout << "FitsIoServer::load(SphereGorski<double> ...) ReSizing to NSide= " << nside  << endl;
+	  cout << "FitsIoServer::load(SphereHealpix<double> ...) ReSizing to NSide= " << nside  << endl;
 	}
       //      else
@@ -504,9 +513,9 @@
 
 
-void FitsIoServer::load(SphereGorski<float>& sph,char flnm[],int hdunum) 
+void FitsIoServer::load(SphereHealpix<float>& sph,char flnm[],int nth) 
 {
   int npixels= 0;
   int nside  = 0;
-
+  int hdunum=nth+1;
   FITS_tab_typ_= TFLOAT;
   DVList dvl;
@@ -528,10 +537,10 @@
 	{
 	  cout<<" FITSIOSERVER: no resolution parameter on fits file "<<endl;
-	  throw IOExc("FitsIoServer::load(SphereGorski<float>& ," + (string)flnm + ", ) Error No resolution parameter !");
+	  throw IOExc("FitsIoServer::load(SphereHealpix<float>& ," + (string)flnm + ", ) Error No resolution parameter !");
 	}
       if (nside != sph.SizeIndex())
 	{ 
 	  sph.Resize(nside);
-	  cout << " FitsIoServer::load(SphereGorski<float> ...)";
+	  cout << " FitsIoServer::load(SphereHealpix<float> ...)";
 	  cout << "  ReSizing to NSide= " << nside << endl;
 	}
@@ -541,13 +550,12 @@
 }
 
-void FitsIoServer::load(SphereGorski<double>& sph,char flnm[],int hdunum) 
+void FitsIoServer::load(SphereHealpix<double>& sph,char flnm[],int nth) 
 {
   int npixels= 0;
   int nside  = 0;
-
+  int hdunum=nth+1;
   FITS_tab_typ_ = TDOUBLE;
   DVList dvl;
   planck_read_bntbl(flnm,hdunum,npixels,dvl);
-
   nside= dvl.GetI("NSIDE");
 
@@ -562,16 +570,97 @@
     {
       if(nside <= 0) 
-	throw IOExc("FitsIoServer::load(SphereGorski<double>& ," + (string)flnm + ", )  No resol parameter !");
+	throw IOExc("FitsIoServer::load(SphereHealpix<double>& ," + (string)flnm + ", )  No resol parameter !");
 
       if(nside != sph.SizeIndex())
 	{ 
 	  sph.Resize(nside);
-	  cout << " FitsIoServer::load(SphereGorski<double> ...)";
+	  cout << " FitsIoServer::load(SphereHealpix<double> ...)";
 	  cout << "  ReSizing to NSide= " << nside << endl;
 	}
     }
+  //
+
+  for(int j = 0; j < sph.NbPixels(); j++) sph(j)= (double)r_8tab_[j];
+}
+
+void FitsIoServer::load(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char flnm[]) 
+{
+  int npixels= 0;
+  int nside  = 0;
+  int hdunum=2;
+  FITS_tab_typ_ = TDOUBLE;
+  DVList dvl;
+  planck_read_bntbl(flnm,hdunum,npixels,dvl,3);
+  nside= dvl.GetI("NSIDE");
+
+  const char* ordering= dvl.GetS("ORDERING").c_str();
+
+  char* ring = "RING";
+  if(strncmp(ordering,ring,4) != 0) 
+    cout << " numerotation non RING" << endl;
+
+  // number of pixels in the spheres
+  if(sph1.NbPixels() != npixels || sph2.NbPixels() != npixels || sph3.NbPixels() != npixels)
+    {
+      if(nside <= 0) 
+	throw IOExc("FitsIoServer::load(SphereHealpix<double>& ," + (string)flnm + ", )  No resol parameter !");
+
+      if(nside != sph1.SizeIndex() || nside != sph2.SizeIndex() || nside != sph3.SizeIndex())
+	{ 
+	  sph1.Resize(nside);
+	  sph2.Resize(nside);
+	  sph3.Resize(nside);
+	  cout << " FitsIoServer::load(SphereHealpix<double> ...)";
+	  cout << "  ReSizing to NSide= " << nside << endl;
+	}
+    }
   
-  for(int j = 0; j < sph.NbPixels(); j++) sph(j)= (double)r_8tab_[j];
-}
+  for(int j = 0; j < sph1.NbPixels(); j++) sph1(j)= (double)r_8tab_[j];
+  for(int j = 0; j < sph2.NbPixels(); j++) sph2(j)= (double)r_8tab_[j+npixels];
+  for(int j = 0; j < sph3.NbPixels(); j++) sph3(j)= (double)r_8tab_[j+2*npixels];
+}
+
+void FitsIoServer::load(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char flnm[]) 
+{
+  int npixels= 0;
+  int nside  = 0;
+  int hdunum=2;
+  FITS_tab_typ_= TFLOAT;
+  DVList dvl;
+
+  planck_read_bntbl(flnm,hdunum,npixels,dvl,3);
+  nside= dvl.GetI("NSIDE");
+
+  const char* ordering= dvl.GetS("ORDERING").c_str();
+
+  char* ring = "RING";
+  if(strncmp(ordering,ring,4) != 0) 
+    cout << " numerotation non RING" << endl;
+
+
+  // number of pixels in the sphere
+  if(sph1.NbPixels() != npixels || sph2.NbPixels() != npixels || sph3.NbPixels() != npixels)
+    {
+      if (nside <= 0 )  
+	{
+	  cout<<" FITSIOSERVER: no resolution parameter on fits file "<<endl;
+	  throw IOExc("FitsIoServer::load(SphereHealpix<float>& ," + (string)flnm + ", ) Error No resolution parameter !");
+	}
+      if(nside != sph1.SizeIndex() || nside != sph2.SizeIndex() || nside != sph3.SizeIndex())
+	{ 
+	  sph1.Resize(nside);
+	  sph2.Resize(nside);
+	  sph3.Resize(nside);
+	  cout << " FitsIoServer::load(SphereHealpix<float> ...)";
+	  cout << "  ReSizing to NSide= " << nside << endl;
+	}
+    }
+  
+  for(int j = 0; j < sph1.NbPixels(); j++) sph1(j)= (float)r_4tab_[j];
+  for(int j = 0; j < sph2.NbPixels(); j++) sph2(j)= (float)r_4tab_[j+npixels];
+  for(int j = 0; j < sph3.NbPixels(); j++) sph3(j)= (float)r_4tab_[j+2*npixels];
+}
+
+
 
 void FitsIoServer::load(SphericalMap<float>& sph, char flnm[]) 
@@ -954,44 +1043,4 @@
   planck_write_img(filename, 1, npixels, 0, 0, dvl);
 
-  // decider ulterieurement ce qu'on fait de ce qui suit, specifique
-  // pour l'instant, aux spheres gorski.
-
-  /*
-  // write supplementary keywords
-  fits_write_comment(fptr, " ", &status);
-  if( status ) printerror( status );
-
-
-  strcpy(comment,"HEALPIX Pixelisation");
-  strcpy(svalue, "HEALPIX");
-  fits_write_key(fptr, TSTRING, "PIXTYPE", svalue, comment, &status);
-  if( status )  printerror( status );
-
-
-  strcpy(comment,"pixel ordering scheme, either RING or NESTED");
-  strcpy(svalue, "RING");
-  fits_write_key(fptr, TSTRING, "ORDERING", svalue, comment, &status);
-  if( status )  printerror( status );
-
-
-  strcpy(comment,"Random generator seed");
-  int iseed= sph.iseed();
-  fits_write_key(fptr, TINT, "RANDSEED", &iseed, comment, &status);
-  if( status )  printerror( status );
-
-
-  strcpy(comment,"--------------------------------------------");
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
-
-  strcpy(comment," Above keywords are still likely to change !");
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
-
-  strcpy(comment,"--------------------------------------------");
-  fits_write_comment(fptr, comment, &status);
-  if( status ) printerror( status );
-  
-  */
 
  delete[] r_8tab_;  r_8tab_ = NULL;
@@ -1015,5 +1064,5 @@
 
 }
-void FitsIoServer::save(SphereGorski<float>& sph, char filename[])   
+void FitsIoServer::save(SphereHealpix<float>& sph, char filename[], int nth)   
 {
   int npixels = sph.NbPixels();
@@ -1033,9 +1082,9 @@
   char* extname="SIMULATION";
   char* comment1="     Sky Map Pixelisation Specific Keywords";
-  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl);
+  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,1,nth);
   delete[] r_4tab_;  r_4tab_ = NULL;
 
 }
-void FitsIoServer::save(SphereGorski<double>& sph, char filename[])   
+void FitsIoServer::save(SphereHealpix<double>& sph, char filename[], int nth)   
 {
   int npixels = sph.NbPixels();
@@ -1055,10 +1104,72 @@
   char* extname="SIMULATION";
   char* comment1="     Sky Map Pixelisation Specific Keywords";
-  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl);
+  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,1, nth);
   delete[] r_8tab_;  r_8tab_ = NULL;
 
 }
-
-void FitsIoServer::save(LocalMap<double>& locm, char filename[])   
+void FitsIoServer::save(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char filename[])   
+{
+  int npixels = sph1.NbPixels();
+  if (npixels !=  sph2.NbPixels() || npixels !=  sph3.NbPixels())
+    {
+      cout << " trois spheres avec differents nombres de pixels !" << endl;
+    }
+  FITS_tab_typ_ = TDOUBLE;
+  if (r_8tab_ != NULL ) { delete[] r_8tab_;  r_8tab_ = NULL; }
+  r_8tab_=new r_8[3*npixels];
+
+
+  for (int j = 0; j < npixels; j++) r_8tab_[j]= (r_8)sph1(j);
+  for (int j = 0; j < npixels; j++) r_8tab_[j+npixels]= (r_8)sph2(j);
+  for (int j = 0; j < npixels; j++) r_8tab_[j+2*npixels]= (r_8)sph3(j);
+  DVList dvl;
+  dvl.SetS("PIXTYPE","HEALPIX ");
+
+  // ulterieurement: verifier que les 3 spheres ont memes caracterisitiques
+  dvl["ORDERING"]=sph1.TypeOfMap();
+  dvl["NSIDE"] = (int_4)sph1.SizeIndex();
+
+  dvl["FIRSTPIX"]=(int_4)0;
+  dvl["LASTPIX"]=(int_4)(npixels-1);
+  char* typeOfContent="TEMPERATURE+POLAR";
+  char* extname="SIMULATION";
+  char* comment1="     Sky Map Pixelisation Specific Keywords";
+  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,3);
+  delete[] r_8tab_;  r_8tab_ = NULL;
+
+}
+
+void FitsIoServer::save(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char filename[])   
+{
+  int npixels = sph1.NbPixels();
+  if (npixels !=  sph2.NbPixels() || npixels !=  sph3.NbPixels())
+    {
+      cout << " trois spheres avec differents nombres de pixels !" << endl;
+    }
+  FITS_tab_typ_ = TFLOAT;
+  if (r_4tab_ != NULL ) { delete[] r_4tab_;  r_4tab_ = NULL; }
+  r_4tab_=new r_4[3*npixels];
+
+  for (int j = 0; j < npixels; j++) r_4tab_[j]= (r_4)sph1(j);
+  for (int j = 0; j < npixels; j++) r_4tab_[j+npixels]= (r_4)sph2(j);
+  for (int j = 0; j < npixels; j++) r_4tab_[j+2*npixels]= (r_4)sph3(j);
+
+  DVList dvl;
+  dvl.SetS("PIXTYPE","HEALPIX ");
+  // ulterieurement: verifier que les 3 spheres ont memes caracterisitiques
+  dvl["ORDERING"]=sph1.TypeOfMap();
+  dvl["NSIDE"] = (int_4)sph1.SizeIndex();
+  dvl["FIRSTPIX"]=(int_4)0;
+  dvl["LASTPIX"]=(int_4)(npixels-1);
+  char* typeOfContent="TEMPERATURE+POLAR";
+  char* extname="SIMULATION";
+  char* comment1="     Sky Map Pixelisation Specific Keywords";
+  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,3);
+  delete[] r_4tab_;  r_4tab_ = NULL;
+
+}
+
+
+void FitsIoServer::save(const LocalMap<double>& locm, char filename[])   
 {
   int nbrows = locm.Size_x();
@@ -1264,5 +1375,5 @@
       cout << " FitsIOServer : type de tableau non traite en ecriture " << endl;
       break;
-    }
+    }  
 
   // write the current date
@@ -1336,7 +1447,6 @@
 }
 
-void FitsIoServer::planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl)
-{
-
+void FitsIoServer::planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl, int tfields,  int nth)
+{
 
   // pointer to the FITS file, defined in fitsio.h
@@ -1345,34 +1455,61 @@
   // 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);
-  //DBG  cerr << " DBG - Apres fits_create_file status = " << status << endl;
-  if( status )  printerror( status );
-
-  // primary header
-  //  int bitpix=LONG_IMG;
-  //  int naxis=0;
-  // fits_create_img(fptr, bitpix, naxis, NULL, &status);
-  // write the current date
-  //  fits_write_date(fptr, &status);
-  //DBG  cerr << " DBG - Apres write_date status = " << status << endl;
-  // if( status )  printerror( status );
-
-
-  long nrows=npixels/1024;
-  if (nrows==0) nrows=1;
-  if (npixels%1024 !=0) nrows++;
-  int tfields=1;
+  if (nth == 1 ) 
+    {
+      // delete old file if it already exists
+      remove(flnm);
+      // create new FITS file  
+      fits_create_file(&fptr, flnm, &status);
+      //DBG  cerr << " DBG - Apres fits_create_file status = " << status << endl;
+      if( status )  printerror( status );
+    }
+  else
+    {
+      fits_open_file(&fptr, flnm,READWRITE, &status);
+      int hdunum=0;
+      fits_get_num_hdus(fptr, &hdunum, &status); 
+      if( status )  printerror( status, " write_bntbl: probleme dans le nombre de headers");
+      cout << " on a retrouve " << hdunum << " headers" << endl;
+      for (int k=1; k<hdunum; k++) 
+	{
+	  fits_movrel_hdu(fptr, 1, NULL, &status); 
+	  if( status )  printerror( status," write_bntbl: probleme pour avancer au header suivant" );
+	  int hdutype=0;
+	  fits_get_hdu_type(fptr, &hdutype, &status);
+	  if( status )  printerror( status," write_bntbl: probleme pour lire un type de header" );
+	  if (hdutype != BINARY_TBL )
+	    {
+	      cout << " le fichier contient autre chose que des binary table: " << hdutype << endl;
+	    }
+	}
+      if (hdunum != nth) 
+	{
+	  cout <<" il n'y a que " << hdunum-1 << "write_bntbl: spheres stockees, la presente sphere sera en position " << hdunum << endl;
+	}
+    }
+
+
+
+  long nrows;
+  char format[5];
+  if (npixels < 1024) 
+    {
+      nrows = npixels;
+      strcpy(format,"1");
+    }
+  else
+    {
+      nrows = npixels/1024;
+      strcpy(format,"1024");
+      if (npixels%1024 !=0) nrows++;
+    }
+  //  int tfields=1;
   char **ttype, **tform;
   ttype= new char*[tfields];
   tform= new char*[tfields];
-  char* format;
-  if ( FITS_tab_typ_ == TDOUBLE) format="1024D";
-  if ( FITS_tab_typ_ == TFLOAT)  format="1024E";
-  if ( FITS_tab_typ_ == TINT)    format="1024I";
+  char* ff="D";
+  if ( FITS_tab_typ_ == TDOUBLE) strcat(format,"D");
+  if ( FITS_tab_typ_ == TFLOAT)  strcat(format,"E");
+  if ( FITS_tab_typ_ == TINT)    strcat(format,"I");
   for(int i = 0; i < tfields; i++)
     {
@@ -1400,19 +1537,27 @@
     {
     case TDOUBLE :
-      fits_write_col(fptr, TDOUBLE, 1, 1, 1, nelements, r_8tab_, &status);
-      if( status ) 
-        printerror( status, "planck_write_bntbl: Error writing double table" );
+      for (int k=1; k<=tfields; k++)
+	{
+	  fits_write_col(fptr, TDOUBLE, k, 1, 1, nelements, r_8tab_+(k-1)*nelements, &status);
+	  if( status ) 
+	    printerror( status, "planck_write_bntbl: Error writing double table" );
+	}
       break;
 
     case TFLOAT :
  //DBG!!! $CHECK$ Reza      for (int kk=0; kk<10; kk++) cout << r_4tab_[kk] << endl;
-      fits_write_col(fptr, TFLOAT, 1, 1, 1, nelements, r_4tab_, &status);
-      if( status ) 
-        printerror( status, "planck_write_bntbl: Error writing float table" );
+      for (int k=1; k<=tfields; k++)
+	{
+	  fits_write_col(fptr, TFLOAT, k, 1, 1, nelements, r_4tab_+(k-1)*nelements, &status);
+	  if( status ) printerror( status, "planck_write_bntbl: Error writing float table" );
+	}
       break;
     case TINT :
-      fits_write_col(fptr, TINT, 1, 1, 1, nelements, i_4tab_, &status);
-      if( status ) 
-        printerror( status, "planck_write_bntbl: Error writing int table");
+      for (int k=1; k<=tfields; k++)
+	{
+	  fits_write_col(fptr, TINT, k, 1, 1, nelements, i_4tab_+(k-1)*nelements, &status);
+	  if( status ) 
+	    printerror( status, "planck_write_bntbl: Error writing int table");
+	}
       break;
     default :
@@ -1636,10 +1781,10 @@
 }
 
-void FitsIoServer::planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl)
+void FitsIoServer::planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl,int tfields )
 {
   int status=0;
   int nkeys,keypos;
   int hdutype;
-  int tfields;
+  //  int tfields;
   int datype;
   long lastpix;
@@ -1686,10 +1831,9 @@
    //cout << " nombre de mots-cles : " << nkeys << endl;
   fits_get_num_cols(fptr, &tfields, &status);
-  if (tfields != 1)
-    {
-      cout << "FitsIoServer:: il y a plus d'une colonne" << endl;
-      throw IOExc("FitsIoServer::planck_read_bntbl(" + (string)flnm + ") Error >1 column !");
-//      return;
-    }
+  //  if (tfields != 1)
+  //  {
+  //    cout << "FitsIoServer:: il y a plus d'une colonne" << endl;
+  //    throw IOExc("FitsIoServer::planck_read_bntbl(" + (string)flnm + ") Error >1 column !");
+  //  }
   fits_get_num_rows(fptr, &nrows, &status);
   //cout << "nblignes= " << nrows << endl;
@@ -1734,9 +1878,12 @@
 	}
       if (r_8tab_ != NULL) { delete [] r_8tab_;  r_8tab_ = NULL; }
-      r_8tab_=new r_8[ npixels];
-      fits_read_col(fptr, TDOUBLE, 1, 1, 1, npixels,  &dnullval, 
-		    r_8tab_, 
+      r_8tab_=new r_8[ tfields*npixels];
+      for (int k=1; k<=tfields; k++)
+	{
+	  fits_read_col(fptr, TDOUBLE, k, 1, 1, npixels,  &dnullval, 
+		    r_8tab_+(k-1)*npixels, 
 		    &anynull, &status);
       if( status )  printerror( status, "probleme dans lecture du tableau de doubles" );
+	}
       break;
     case TFLOAT :
@@ -1748,8 +1895,11 @@
 	}
       if (r_4tab_ != NULL) { delete [] r_4tab_;  r_4tab_ = NULL; }
-      r_4tab_=new r_4[npixels];
-      fits_read_col(fptr, TFLOAT, 1, 1, 1, npixels,  &fnullval, 
-		    r_4tab_, &anynull, &status);
-      if( status )  printerror( status,"probleme dans lecture du tableau de floats" );
+      r_4tab_=new r_4[ tfields*npixels];
+      for (int k=1; k<=tfields; k++)
+	{
+	  fits_read_col(fptr, TFLOAT, k, 1, 1, npixels,  &fnullval, 
+		    r_4tab_+(k-1)*npixels, &anynull, &status);
+	  if( status )  printerror( status,"probleme dans lecture du tableau de floats" );
+	}
       break;
 
@@ -1762,10 +1912,11 @@
 	}
       if (i_4tab_ != NULL) { delete [] i_4tab_;  i_4tab_ = NULL; }
-      i_4tab_=new int_4[npixels];
-      fits_read_col(fptr, TLONG, 1, 1, 1, npixels,  &inullval, 
-		    i_4tab_, &anynull, &status);
-      //fits_read_img(fptr, TINT, 1, nelements,  &inullval, i_4tab_, 
-      //		    &anynull, &status);
-      if( status )  printerror( status,"probleme dans lecture du tableau de ints" );
+      i_4tab_=new int_4[ tfields*npixels];
+      for (int k=1; k<=tfields; k++)
+	{
+	  fits_read_col(fptr, TLONG, k, 1, 1, npixels,  &inullval, 
+		    i_4tab_+(k-1)*npixels, &anynull, &status);
+	  if( status )  printerror( status,"probleme dans lecture du tableau de ints" );
+	}
       break;
 
@@ -1784,4 +1935,5 @@
   char blank[LEN_KEYWORD]= "";
 
+  // 
    for(int j = 1; j <= nkeys; j++)
      {
Index: trunk/SophyaExt/FitsIOServer/fitsioserver.h
===================================================================
--- trunk/SophyaExt/FitsIOServer/fitsioserver.h	(revision 837)
+++ trunk/SophyaExt/FitsIOServer/fitsioserver.h	(revision 838)
@@ -5,12 +5,12 @@
 #include "tvector.h"
 #include "sphericalmap.h"
-#include "spheregorski.h"
+#include "spherehealpix.h"
 #include "localmap.h"
 #include "ntuple.h"
+//#include "outilsinit.h"
 #include "cimage.h"
 #include "dvlist.h"
 #include "anydataobj.h"
-
-#include "FitsIO/fitsio.h"   
+#include "FitsIO/fitsio.h"    
 
 #define LEN_KEYWORD 9
@@ -39,6 +39,8 @@
  void load(ImageR4& DpcImg,char flnm[]);
  void load(ImageI4& DpcImg,char flnm[]);
- void load(SphereGorski<float>& sph, char flnm[], int hdunum=2); 
- void load(SphereGorski<double>& sph, char flnm[], int hdunum=2); 
+ void load(SphereHealpix<float>& sph, char flnm[], int nth=1); 
+ void load(SphereHealpix<double>& sph, char flnm[], int nth=1); 
+ void load(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char filename[]); 
+ void load(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char filename[]); 
  void save(TMatrix<double>& mat, char flnm[]);
  void save(TMatrix<float>& mat, char flnm[]);
@@ -47,7 +49,15 @@
  void save(SphericalMap<double>& sph, char flnm[]);
  void save(SphericalMap<float>& sph, char flnm[]);
- void save(SphereGorski<float>& sph, char filename[]);   
- void save(SphereGorski<double>& sph, char filename[]);   
- void save(LocalMap<double>& locm, char flnm[]);
+
+ // nth : numero d'ordre de la sphere sur le fichier de sauvegarde
+ // si on veut sauver une seule sphere : nth=1
+ // so veut sauver par exemple 3 spheres, on appelle 3 fois save
+ // avec respectivement nth=1, 2, 3  
+ void save(SphereHealpix<float>& sph, char filename[], int nth=1); 
+ void save(SphereHealpix<double>& sph, char filename[], int nth=1); 
+  
+ void save(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char filename[]);   
+ void save(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char filename[]);   
+ void save(const LocalMap<double>& locm, char flnm[]);
  void save(const ImageR4& DpcImg,char flnm[]);
  void save(const ImageI4& DpcImg,char flnm[]);
@@ -68,6 +78,6 @@
 
  void planck_read_img(char flnm[],int &naxis,int &n1,int &n2,int &n3,DVList &dvl);
- void planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl);
- void planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl);
+ void planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl,int tfields=1);
+ void planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl, int tfields=1, int nth=1);
  void write_picture(long* naxes,float* map,char* filename) const;
 
