Changeset 838 in Sophya for trunk/SophyaExt


Ignore:
Timestamp:
Apr 7, 2000, 5:43:14 PM (25 years ago)
Author:
ansari
Message:

gorski change en healpix

Location:
trunk/SophyaExt/FitsIOServer
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/SophyaExt/FitsIOServer/fitsioserver.cc

    r833 r838  
    1010//************************************************************************
    1111
    12 #include "machdefs.h"
     12#include "machdefs.h"  
    1313#include <iostream.h>
    1414
    1515// #include <sstream>  Ca ne passe pas pour le moment avec g++ - Reza 23/12/99
    1616
    17 #include <list>
     17//#include <list>
    1818#include <string>
    1919
    … …  
    8080          if(DTYPE == TDOUBLE)
    8181            {
    82               SphereGorski<double> *sph= new SphereGorski<double>(nside);
     82              SphereHealpix<double> *sph= new SphereHealpix<double>(nside);
    8383              load(*sph,flnm,hdunum);
    8484              return sph;
    … …  
    8787          else if(DTYPE == TFLOAT)
    8888            {
    89               SphereGorski<float> *sph= new SphereGorski<float>(nside);
     89              SphereHealpix<float> *sph= new SphereHealpix<float>(nside);
    9090              load(*sph,flnm,hdunum);
    9191              return sph;
    … …  
    313313    }
    314314
     315
    315316  // get the number of columns
    316317  int ncols= 0;
    … …  
    328329  long repeat;
    329330  int ntest= 0;
    330   int ii;
    331   for(ii = 0; ii < ncols; ii++)
     331  for(int ii = 0; ii < ncols; ii++)
    332332    {
    333333      int DTYPE;
    … …  
    354354  char **clname;
    355355  clname= new char*[tfields.size()];
    356   for(ii = 0; ii < tfields.size(); ii++)
     356  for(int ii = 0; ii < tfields.size(); ii++)
    357357    {
    358358      //  ostringstream ne passe pas sur tous les compilos .. Reza 23/12/99
    359359      //      ostringstream oss; oss<<"Col_"<<ii+1;
    360       //      string ss= oss.str();
     360      //      string ss= oss.str(); 
    361361      //      clname[ii]= new char[FLEN_VALUE];
    362362      clname[ii] = new char[16];
    … …  
    367367  // create a NTuple
    368368  NTuple nt0(tfields.size(),clname);
    369   for(ii = 0; ii < tfields.size(); ii++)
     369  for(int ii = 0; ii < tfields.size(); ii++)
    370370    delete [] clname[ii];
    371371  delete [] clname;
    … …  
    396396    }
    397397
     398
     399
    398400  // get number of keywords
    399401  int nkeys,keypos;
    … …  
    411413  int num= 8;
    412414
     415
     416
     417
    413418  for(int j = num+1; j <= nkeys; j++)
    414419    {
    … …  
    427432          strip(strval, 'B',' ');
    428433          strip(strval, 'B','\'');
     434         
    429435          switch( dtype ) 
    430436            {
    … …  
    435441              int ival;
    436442              fits_read_key(fptr,TINT,keyname,&ival,NULL,&status);
    437               nt0.Info()[keyname]= (int_4) ival;                // Portage mac DY
     443                nt0.Info()[keyname]= (int_4) ival;              // Portage mac DY
    438444              break;
    439445            case 'L':
    … …  
    449455              break;
    450456            }
     457         
    451458        }
    452459    }
    453460 
     461
     462
    454463  // copy in the input NTuple
    455464  ntpl= nt0;
    … …  
    492501        {
    493502          sph.Resize(nside);
    494           cout << "FitsIoServer::load(SphereGorski<double> ...) ReSizing to NSide= " << nside  << endl;
     503          cout << "FitsIoServer::load(SphereHealpix<double> ...) ReSizing to NSide= " << nside  << endl;
    495504        }
    496505      //      else
    … …  
    504513
    505514
    506 void FitsIoServer::load(SphereGorski<float>& sph,char flnm[],int hdunum)
     515void FitsIoServer::load(SphereHealpix<float>& sph,char flnm[],int nth)
    507516{
    508517  int npixels= 0;
    509518  int nside  = 0;
    510 
     519  int hdunum=nth+1;
    511520  FITS_tab_typ_= TFLOAT;
    512521  DVList dvl;
    … …  
    528537        {
    529538          cout<<" FITSIOSERVER: no resolution parameter on fits file "<<endl;
    530           throw IOExc("FitsIoServer::load(SphereGorski<float>& ," + (string)flnm + ", ) Error No resolution parameter !");
     539          throw IOExc("FitsIoServer::load(SphereHealpix<float>& ," + (string)flnm + ", ) Error No resolution parameter !");
    531540        }
    532541      if (nside != sph.SizeIndex())
    533542        {
    534543          sph.Resize(nside);
    535           cout << " FitsIoServer::load(SphereGorski<float> ...)";
     544          cout << " FitsIoServer::load(SphereHealpix<float> ...)";
    536545          cout << "  ReSizing to NSide= " << nside << endl;
    537546        }
    … …  
    541550}
    542551
    543 void FitsIoServer::load(SphereGorski<double>& sph,char flnm[],int hdunum)
     552void FitsIoServer::load(SphereHealpix<double>& sph,char flnm[],int nth)
    544553{
    545554  int npixels= 0;
    546555  int nside  = 0;
    547 
     556  int hdunum=nth+1;
    548557  FITS_tab_typ_ = TDOUBLE;
    549558  DVList dvl;
    550559  planck_read_bntbl(flnm,hdunum,npixels,dvl);
    551 
    552560  nside= dvl.GetI("NSIDE");
    553561
    … …  
    562570    {
    563571      if(nside <= 0)
    564         throw IOExc("FitsIoServer::load(SphereGorski<double>& ," + (string)flnm + ", )  No resol parameter !");
     572        throw IOExc("FitsIoServer::load(SphereHealpix<double>& ," + (string)flnm + ", )  No resol parameter !");
    565573
    566574      if(nside != sph.SizeIndex())
    567575        {
    568576          sph.Resize(nside);
    569           cout << " FitsIoServer::load(SphereGorski<double> ...)";
     577          cout << " FitsIoServer::load(SphereHealpix<double> ...)";
    570578          cout << "  ReSizing to NSide= " << nside << endl;
    571579        }
    572580    }
     581  //
     582
     583  for(int j = 0; j < sph.NbPixels(); j++) sph(j)= (double)r_8tab_[j];
     584}
     585
     586void FitsIoServer::load(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char flnm[])
     587{
     588  int npixels= 0;
     589  int nside  = 0;
     590  int hdunum=2;
     591  FITS_tab_typ_ = TDOUBLE;
     592  DVList dvl;
     593  planck_read_bntbl(flnm,hdunum,npixels,dvl,3);
     594  nside= dvl.GetI("NSIDE");
     595
     596  const char* ordering= dvl.GetS("ORDERING").c_str();
     597
     598  char* ring = "RING";
     599  if(strncmp(ordering,ring,4) != 0)
     600    cout << " numerotation non RING" << endl;
     601
     602  // number of pixels in the spheres
     603  if(sph1.NbPixels() != npixels || sph2.NbPixels() != npixels || sph3.NbPixels() != npixels)
     604    {
     605      if(nside <= 0)
     606        throw IOExc("FitsIoServer::load(SphereHealpix<double>& ," + (string)flnm + ", )  No resol parameter !");
     607
     608      if(nside != sph1.SizeIndex() || nside != sph2.SizeIndex() || nside != sph3.SizeIndex())
     609        {
     610          sph1.Resize(nside);
     611          sph2.Resize(nside);
     612          sph3.Resize(nside);
     613          cout << " FitsIoServer::load(SphereHealpix<double> ...)";
     614          cout << "  ReSizing to NSide= " << nside << endl;
     615        }
     616    }
    573617 
    574   for(int j = 0; j < sph.NbPixels(); j++) sph(j)= (double)r_8tab_[j];
    575 }
     618  for(int j = 0; j < sph1.NbPixels(); j++) sph1(j)= (double)r_8tab_[j];
     619  for(int j = 0; j < sph2.NbPixels(); j++) sph2(j)= (double)r_8tab_[j+npixels];
     620  for(int j = 0; j < sph3.NbPixels(); j++) sph3(j)= (double)r_8tab_[j+2*npixels];
     621}
     622
     623void FitsIoServer::load(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char flnm[])
     624{
     625  int npixels= 0;
     626  int nside  = 0;
     627  int hdunum=2;
     628  FITS_tab_typ_= TFLOAT;
     629  DVList dvl;
     630
     631  planck_read_bntbl(flnm,hdunum,npixels,dvl,3);
     632  nside= dvl.GetI("NSIDE");
     633
     634  const char* ordering= dvl.GetS("ORDERING").c_str();
     635
     636  char* ring = "RING";
     637  if(strncmp(ordering,ring,4) != 0)
     638    cout << " numerotation non RING" << endl;
     639
     640
     641  // number of pixels in the sphere
     642  if(sph1.NbPixels() != npixels || sph2.NbPixels() != npixels || sph3.NbPixels() != npixels)
     643    {
     644      if (nside <= 0 ) 
     645        {
     646          cout<<" FITSIOSERVER: no resolution parameter on fits file "<<endl;
     647          throw IOExc("FitsIoServer::load(SphereHealpix<float>& ," + (string)flnm + ", ) Error No resolution parameter !");
     648        }
     649      if(nside != sph1.SizeIndex() || nside != sph2.SizeIndex() || nside != sph3.SizeIndex())
     650        {
     651          sph1.Resize(nside);
     652          sph2.Resize(nside);
     653          sph3.Resize(nside);
     654          cout << " FitsIoServer::load(SphereHealpix<float> ...)";
     655          cout << "  ReSizing to NSide= " << nside << endl;
     656        }
     657    }
     658 
     659  for(int j = 0; j < sph1.NbPixels(); j++) sph1(j)= (float)r_4tab_[j];
     660  for(int j = 0; j < sph2.NbPixels(); j++) sph2(j)= (float)r_4tab_[j+npixels];
     661  for(int j = 0; j < sph3.NbPixels(); j++) sph3(j)= (float)r_4tab_[j+2*npixels];
     662}
     663
     664
    576665
    577666void FitsIoServer::load(SphericalMap<float>& sph, char flnm[])
    … …  
    9541043  planck_write_img(filename, 1, npixels, 0, 0, dvl);
    9551044
    956   // decider ulterieurement ce qu'on fait de ce qui suit, specifique
    957   // pour l'instant, aux spheres gorski.
    958 
    959   /*
    960   // write supplementary keywords
    961   fits_write_comment(fptr, " ", &status);
    962   if( status ) printerror( status );
    963 
    964 
    965   strcpy(comment,"HEALPIX Pixelisation");
    966   strcpy(svalue, "HEALPIX");
    967   fits_write_key(fptr, TSTRING, "PIXTYPE", svalue, comment, &status);
    968   if( status )  printerror( status );
    969 
    970 
    971   strcpy(comment,"pixel ordering scheme, either RING or NESTED");
    972   strcpy(svalue, "RING");
    973   fits_write_key(fptr, TSTRING, "ORDERING", svalue, comment, &status);
    974   if( status )  printerror( status );
    975 
    976 
    977   strcpy(comment,"Random generator seed");
    978   int iseed= sph.iseed();
    979   fits_write_key(fptr, TINT, "RANDSEED", &iseed, comment, &status);
    980   if( status )  printerror( status );
    981 
    982 
    983   strcpy(comment,"--------------------------------------------");
    984   fits_write_comment(fptr, comment, &status);
    985   if( status ) printerror( status );
    986 
    987   strcpy(comment," Above keywords are still likely to change !");
    988   fits_write_comment(fptr, comment, &status);
    989   if( status ) printerror( status );
    990 
    991   strcpy(comment,"--------------------------------------------");
    992   fits_write_comment(fptr, comment, &status);
    993   if( status ) printerror( status );
    994  
    995   */
    9961045
    9971046 delete[] r_8tab_;  r_8tab_ = NULL;
    … …  
    10151064
    10161065}
    1017 void FitsIoServer::save(SphereGorski<float>& sph, char filename[])   
     1066void FitsIoServer::save(SphereHealpix<float>& sph, char filename[], int nth)   
    10181067{
    10191068  int npixels = sph.NbPixels();
    … …  
    10331082  char* extname="SIMULATION";
    10341083  char* comment1="     Sky Map Pixelisation Specific Keywords";
    1035   planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl);
     1084  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,1,nth);
    10361085  delete[] r_4tab_;  r_4tab_ = NULL;
    10371086
    10381087}
    1039 void FitsIoServer::save(SphereGorski<double>& sph, char filename[])   
     1088void FitsIoServer::save(SphereHealpix<double>& sph, char filename[], int nth)   
    10401089{
    10411090  int npixels = sph.NbPixels();
    … …  
    10551104  char* extname="SIMULATION";
    10561105  char* comment1="     Sky Map Pixelisation Specific Keywords";
    1057   planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl);
     1106  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,1, nth);
    10581107  delete[] r_8tab_;  r_8tab_ = NULL;
    10591108
    10601109}
    1061 
    1062 void FitsIoServer::save(LocalMap<double>& locm, char filename[])   
     1110void FitsIoServer::save(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char filename[])   
     1111{
     1112  int npixels = sph1.NbPixels();
     1113  if (npixels !=  sph2.NbPixels() || npixels !=  sph3.NbPixels())
     1114    {
     1115      cout << " trois spheres avec differents nombres de pixels !" << endl;
     1116    }
     1117  FITS_tab_typ_ = TDOUBLE;
     1118  if (r_8tab_ != NULL ) { delete[] r_8tab_;  r_8tab_ = NULL; }
     1119  r_8tab_=new r_8[3*npixels];
     1120
     1121
     1122  for (int j = 0; j < npixels; j++) r_8tab_[j]= (r_8)sph1(j);
     1123  for (int j = 0; j < npixels; j++) r_8tab_[j+npixels]= (r_8)sph2(j);
     1124  for (int j = 0; j < npixels; j++) r_8tab_[j+2*npixels]= (r_8)sph3(j);
     1125  DVList dvl;
     1126  dvl.SetS("PIXTYPE","HEALPIX ");
     1127
     1128  // ulterieurement: verifier que les 3 spheres ont memes caracterisitiques
     1129  dvl["ORDERING"]=sph1.TypeOfMap();
     1130  dvl["NSIDE"] = (int_4)sph1.SizeIndex();
     1131
     1132  dvl["FIRSTPIX"]=(int_4)0;
     1133  dvl["LASTPIX"]=(int_4)(npixels-1);
     1134  char* typeOfContent="TEMPERATURE+POLAR";
     1135  char* extname="SIMULATION";
     1136  char* comment1="     Sky Map Pixelisation Specific Keywords";
     1137  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,3);
     1138  delete[] r_8tab_;  r_8tab_ = NULL;
     1139
     1140}
     1141
     1142void FitsIoServer::save(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char filename[])   
     1143{
     1144  int npixels = sph1.NbPixels();
     1145  if (npixels !=  sph2.NbPixels() || npixels !=  sph3.NbPixels())
     1146    {
     1147      cout << " trois spheres avec differents nombres de pixels !" << endl;
     1148    }
     1149  FITS_tab_typ_ = TFLOAT;
     1150  if (r_4tab_ != NULL ) { delete[] r_4tab_;  r_4tab_ = NULL; }
     1151  r_4tab_=new r_4[3*npixels];
     1152
     1153  for (int j = 0; j < npixels; j++) r_4tab_[j]= (r_4)sph1(j);
     1154  for (int j = 0; j < npixels; j++) r_4tab_[j+npixels]= (r_4)sph2(j);
     1155  for (int j = 0; j < npixels; j++) r_4tab_[j+2*npixels]= (r_4)sph3(j);
     1156
     1157  DVList dvl;
     1158  dvl.SetS("PIXTYPE","HEALPIX ");
     1159  // ulterieurement: verifier que les 3 spheres ont memes caracterisitiques
     1160  dvl["ORDERING"]=sph1.TypeOfMap();
     1161  dvl["NSIDE"] = (int_4)sph1.SizeIndex();
     1162  dvl["FIRSTPIX"]=(int_4)0;
     1163  dvl["LASTPIX"]=(int_4)(npixels-1);
     1164  char* typeOfContent="TEMPERATURE+POLAR";
     1165  char* extname="SIMULATION";
     1166  char* comment1="     Sky Map Pixelisation Specific Keywords";
     1167  planck_write_bntbl(filename, npixels, typeOfContent, extname, comment1, dvl,3);
     1168  delete[] r_4tab_;  r_4tab_ = NULL;
     1169
     1170}
     1171
     1172
     1173void FitsIoServer::save(const LocalMap<double>& locm, char filename[])   
    10631174{
    10641175  int nbrows = locm.Size_x();
    … …  
    12641375      cout << " FitsIOServer : type de tableau non traite en ecriture " << endl;
    12651376      break;
    1266     }
     1377    } 
    12671378
    12681379  // write the current date
    … …  
    13361447}
    13371448
    1338 void FitsIoServer::planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl)
    1339 {
    1340 
     1449void FitsIoServer::planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl, int tfields,  int nth)
     1450{
    13411451
    13421452  // pointer to the FITS file, defined in fitsio.h
    … …  
    13451455  // initialize status before calling fitsio routines
    13461456  int status = 0;         
    1347 
    1348   // delete old file if it already exists
    1349   remove(flnm);
    1350 
    1351   // create new FITS file 
    1352   fits_create_file(&fptr, flnm, &status);
    1353   //DBG  cerr << " DBG - Apres fits_create_file status = " << status << endl;
    1354   if( status )  printerror( status );
    1355 
    1356   // primary header
    1357   //  int bitpix=LONG_IMG;
    1358   //  int naxis=0;
    1359   // fits_create_img(fptr, bitpix, naxis, NULL, &status);
    1360   // write the current date
    1361   //  fits_write_date(fptr, &status);
    1362   //DBG  cerr << " DBG - Apres write_date status = " << status << endl;
    1363   // if( status )  printerror( status );
    1364 
    1365 
    1366   long nrows=npixels/1024;
    1367   if (nrows==0) nrows=1;
    1368   if (npixels%1024 !=0) nrows++;
    1369   int tfields=1;
     1457  if (nth == 1 )
     1458    {
     1459      // delete old file if it already exists
     1460      remove(flnm);
     1461      // create new FITS file 
     1462      fits_create_file(&fptr, flnm, &status);
     1463      //DBG  cerr << " DBG - Apres fits_create_file status = " << status << endl;
     1464      if( status )  printerror( status );
     1465    }
     1466  else
     1467    {
     1468      fits_open_file(&fptr, flnm,READWRITE, &status);
     1469      int hdunum=0;
     1470      fits_get_num_hdus(fptr, &hdunum, &status);
     1471      if( status )  printerror( status, " write_bntbl: probleme dans le nombre de headers");
     1472      cout << " on a retrouve " << hdunum << " headers" << endl;
     1473      for (int k=1; k<hdunum; k++)
     1474        {
     1475          fits_movrel_hdu(fptr, 1, NULL, &status);
     1476          if( status )  printerror( status," write_bntbl: probleme pour avancer au header suivant" );
     1477          int hdutype=0;
     1478          fits_get_hdu_type(fptr, &hdutype, &status);
     1479          if( status )  printerror( status," write_bntbl: probleme pour lire un type de header" );
     1480          if (hdutype != BINARY_TBL )
     1481            {
     1482              cout << " le fichier contient autre chose que des binary table: " << hdutype << endl;
     1483            }
     1484        }
     1485      if (hdunum != nth)
     1486        {
     1487          cout <<" il n'y a que " << hdunum-1 << "write_bntbl: spheres stockees, la presente sphere sera en position " << hdunum << endl;
     1488        }
     1489    }
     1490
     1491
     1492
     1493  long nrows;
     1494  char format[5];
     1495  if (npixels < 1024)
     1496    {
     1497      nrows = npixels;
     1498      strcpy(format,"1");
     1499    }
     1500  else
     1501    {
     1502      nrows = npixels/1024;
     1503      strcpy(format,"1024");
     1504      if (npixels%1024 !=0) nrows++;
     1505    }
     1506  //  int tfields=1;
    13701507  char **ttype, **tform;
    13711508  ttype= new char*[tfields];
    13721509  tform= new char*[tfields];
    1373   char* format;
    1374   if ( FITS_tab_typ_ == TDOUBLE) format="1024D";
    1375   if ( FITS_tab_typ_ == TFLOAT)  format="1024E";
    1376   if ( FITS_tab_typ_ == TINT)    format="1024I";
     1510  char* ff="D";
     1511  if ( FITS_tab_typ_ == TDOUBLE) strcat(format,"D");
     1512  if ( FITS_tab_typ_ == TFLOAT)  strcat(format,"E");
     1513  if ( FITS_tab_typ_ == TINT)    strcat(format,"I");
    13771514  for(int i = 0; i < tfields; i++)
    13781515    {
    … …  
    14001537    {
    14011538    case TDOUBLE :
    1402       fits_write_col(fptr, TDOUBLE, 1, 1, 1, nelements, r_8tab_, &status);
    1403       if( status )
    1404         printerror( status, "planck_write_bntbl: Error writing double table" );
     1539      for (int k=1; k<=tfields; k++)
     1540        {
     1541          fits_write_col(fptr, TDOUBLE, k, 1, 1, nelements, r_8tab_+(k-1)*nelements, &status);
     1542          if( status )
     1543            printerror( status, "planck_write_bntbl: Error writing double table" );
     1544        }
    14051545      break;
    14061546
    14071547    case TFLOAT :
    14081548 //DBG!!! $CHECK$ Reza      for (int kk=0; kk<10; kk++) cout << r_4tab_[kk] << endl;
    1409       fits_write_col(fptr, TFLOAT, 1, 1, 1, nelements, r_4tab_, &status);
    1410       if( status )
    1411         printerror( status, "planck_write_bntbl: Error writing float table" );
     1549      for (int k=1; k<=tfields; k++)
     1550        {
     1551          fits_write_col(fptr, TFLOAT, k, 1, 1, nelements, r_4tab_+(k-1)*nelements, &status);
     1552          if( status ) printerror( status, "planck_write_bntbl: Error writing float table" );
     1553        }
    14121554      break;
    14131555    case TINT :
    1414       fits_write_col(fptr, TINT, 1, 1, 1, nelements, i_4tab_, &status);
    1415       if( status )
    1416         printerror( status, "planck_write_bntbl: Error writing int table");
     1556      for (int k=1; k<=tfields; k++)
     1557        {
     1558          fits_write_col(fptr, TINT, k, 1, 1, nelements, i_4tab_+(k-1)*nelements, &status);
     1559          if( status )
     1560            printerror( status, "planck_write_bntbl: Error writing int table");
     1561        }
    14171562      break;
    14181563    default :
    … …  
    16361781}
    16371782
    1638 void FitsIoServer::planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl)
     1783void FitsIoServer::planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl,int tfields )
    16391784{
    16401785  int status=0;
    16411786  int nkeys,keypos;
    16421787  int hdutype;
    1643   int tfields;
     1788  //  int tfields;
    16441789  int datype;
    16451790  long lastpix;
    … …  
    16861831   //cout << " nombre de mots-cles : " << nkeys << endl;
    16871832  fits_get_num_cols(fptr, &tfields, &status);
    1688   if (tfields != 1)
    1689     {
    1690       cout << "FitsIoServer:: il y a plus d'une colonne" << endl;
    1691       throw IOExc("FitsIoServer::planck_read_bntbl(" + (string)flnm + ") Error >1 column !");
    1692 //      return;
    1693     }
     1833  //  if (tfields != 1)
     1834  //  {
     1835  //    cout << "FitsIoServer:: il y a plus d'une colonne" << endl;
     1836  //    throw IOExc("FitsIoServer::planck_read_bntbl(" + (string)flnm + ") Error >1 column !");
     1837  //  }
    16941838  fits_get_num_rows(fptr, &nrows, &status);
    16951839  //cout << "nblignes= " << nrows << endl;
    … …  
    17341878        }
    17351879      if (r_8tab_ != NULL) { delete [] r_8tab_;  r_8tab_ = NULL; }
    1736       r_8tab_=new r_8[ npixels];
    1737       fits_read_col(fptr, TDOUBLE, 1, 1, 1, npixels,  &dnullval,
    1738                     r_8tab_,
     1880      r_8tab_=new r_8[ tfields*npixels];
     1881      for (int k=1; k<=tfields; k++)
     1882        {
     1883          fits_read_col(fptr, TDOUBLE, k, 1, 1, npixels,  &dnullval,
     1884                    r_8tab_+(k-1)*npixels,
    17391885                    &anynull, &status);
    17401886      if( status )  printerror( status, "probleme dans lecture du tableau de doubles" );
     1887        }
    17411888      break;
    17421889    case TFLOAT :
    … …  
    17481895        }
    17491896      if (r_4tab_ != NULL) { delete [] r_4tab_;  r_4tab_ = NULL; }
    1750       r_4tab_=new r_4[npixels];
    1751       fits_read_col(fptr, TFLOAT, 1, 1, 1, npixels,  &fnullval,
    1752                     r_4tab_, &anynull, &status);
    1753       if( status )  printerror( status,"probleme dans lecture du tableau de floats" );
     1897      r_4tab_=new r_4[ tfields*npixels];
     1898      for (int k=1; k<=tfields; k++)
     1899        {
     1900          fits_read_col(fptr, TFLOAT, k, 1, 1, npixels,  &fnullval,
     1901                    r_4tab_+(k-1)*npixels, &anynull, &status);
     1902          if( status )  printerror( status,"probleme dans lecture du tableau de floats" );
     1903        }
    17541904      break;
    17551905
    … …  
    17621912        }
    17631913      if (i_4tab_ != NULL) { delete [] i_4tab_;  i_4tab_ = NULL; }
    1764       i_4tab_=new int_4[npixels];
    1765       fits_read_col(fptr, TLONG, 1, 1, 1, npixels,  &inullval,
    1766                     i_4tab_, &anynull, &status);
    1767       //fits_read_img(fptr, TINT, 1, nelements,  &inullval, i_4tab_,
    1768       //                    &anynull, &status);
    1769       if( status )  printerror( status,"probleme dans lecture du tableau de ints" );
     1914      i_4tab_=new int_4[ tfields*npixels];
     1915      for (int k=1; k<=tfields; k++)
     1916        {
     1917          fits_read_col(fptr, TLONG, k, 1, 1, npixels,  &inullval,
     1918                    i_4tab_+(k-1)*npixels, &anynull, &status);
     1919          if( status )  printerror( status,"probleme dans lecture du tableau de ints" );
     1920        }
    17701921      break;
    17711922
    … …  
    17841935  char blank[LEN_KEYWORD]= "";
    17851936
     1937  //
    17861938   for(int j = 1; j <= nkeys; j++)
    17871939     {

  • trunk/SophyaExt/FitsIOServer/fitsioserver.h

    r770 r838  
    55#include "tvector.h"
    66#include "sphericalmap.h"
    7 #include "spheregorski.h"
     7#include "spherehealpix.h"
    88#include "localmap.h"
    99#include "ntuple.h"
     10//#include "outilsinit.h"
    1011#include "cimage.h"
    1112#include "dvlist.h"
    1213#include "anydataobj.h"
    13 
    14 #include "FitsIO/fitsio.h"   
     14#include "FitsIO/fitsio.h"   
    1515
    1616#define LEN_KEYWORD 9
    … …  
    3939 void load(ImageR4& DpcImg,char flnm[]);
    4040 void load(ImageI4& DpcImg,char flnm[]);
    41  void load(SphereGorski<float>& sph, char flnm[], int hdunum=2);
    42  void load(SphereGorski<double>& sph, char flnm[], int hdunum=2);
     41 void load(SphereHealpix<float>& sph, char flnm[], int nth=1);
     42 void load(SphereHealpix<double>& sph, char flnm[], int nth=1);
     43 void load(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char filename[]);
     44 void load(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char filename[]);
    4345 void save(TMatrix<double>& mat, char flnm[]);
    4446 void save(TMatrix<float>& mat, char flnm[]);
    … …  
    4749 void save(SphericalMap<double>& sph, char flnm[]);
    4850 void save(SphericalMap<float>& sph, char flnm[]);
    49  void save(SphereGorski<float>& sph, char filename[]);   
    50  void save(SphereGorski<double>& sph, char filename[]);   
    51  void save(LocalMap<double>& locm, char flnm[]);
     51
     52 // nth : numero d'ordre de la sphere sur le fichier de sauvegarde
     53 // si on veut sauver une seule sphere : nth=1
     54 // so veut sauver par exemple 3 spheres, on appelle 3 fois save
     55 // avec respectivement nth=1, 2, 3 
     56 void save(SphereHealpix<float>& sph, char filename[], int nth=1);
     57 void save(SphereHealpix<double>& sph, char filename[], int nth=1);
     58 
     59 void save(SphereHealpix<double>& sph1,SphereHealpix<double>& sph2,SphereHealpix<double>& sph3, char filename[]);   
     60 void save(SphereHealpix<float>& sph1,SphereHealpix<float>& sph2,SphereHealpix<float>& sph3, char filename[]);   
     61 void save(const LocalMap<double>& locm, char flnm[]);
    5262 void save(const ImageR4& DpcImg,char flnm[]);
    5363 void save(const ImageI4& DpcImg,char flnm[]);
    … …  
    6878
    6979 void planck_read_img(char flnm[],int &naxis,int &n1,int &n2,int &n3,DVList &dvl);
    70  void planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl);
    71  void planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl);
     80 void planck_read_bntbl(char flnm[], int hdunum, int& npixels, DVList& dvl,int tfields=1);
     81 void planck_write_bntbl(char flnm[], int npixels, char typeOfContent[], char extname[], char comment1[], DVList& dvl, int tfields=1, int nth=1);
    7282 void write_picture(long* naxes,float* map,char* filename) const;
    7383
Note: See TracChangeset for help on using the changeset viewer.