Ignore:
Timestamp:
Oct 18, 1999, 4:37:44 PM (26 years ago)
Author:
ansari
Message:

modifs francois : passage en double etc. 18-OCT-99

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/SophyaLib/Samba/spherethetaphi.cc

    r470 r473  
    8787//++
    8888template <class T>
    89 SphereThetaPhi<T>::SphereThetaPhi(int_4 m)
     89SphereThetaPhi<T>::SphereThetaPhi(int m)
    9090
    9191//    Constructeur : m est le nombre de bandes en theta sur un hémisphère
    … …  
    105105  NTheta_= s.NTheta_;
    106106  NPix_  = s.NPix_;
    107   NPhi_  = new int_4[NTheta_];
    108   Theta_ = new r_4[NTheta_+1];
    109   TNphi_ = new int_4[NTheta_+1];
     107  NPhi_  = new int[NTheta_];
     108  Theta_ = new double[NTheta_+1];
     109  TNphi_ = new int[NTheta_+1];
    110110  for(int k = 0; k < NTheta_; k++)
    111111    {
    … …  
    161161//++
    162162template <class T>
    163 void SphereThetaPhi<T>::Resize(int_4 m)
     163void SphereThetaPhi<T>::Resize(int m)
    164164//   re-pixelize the sphere
    165165//--
    … …  
    172172//++
    173173template <class T>
    174 int_4 SphereThetaPhi<T>::NbPixels() const
     174int SphereThetaPhi<T>::NbPixels() const
    175175
    176176//    Retourne le nombre de pixels du découpage
    … …  
    183183//++
    184184template <class T>
    185 T& SphereThetaPhi<T>::PixVal(int_4 k)
     185T& SphereThetaPhi<T>::PixVal(int k)
    186186
    187187//    Retourne la valeur du contenu du pixel d'indice k
    … …  
    190190  if((k < 0) || (k >= NPix_))
    191191    {
    192       //  THROW(out_of_range("SphereThetaPhi::PIxVal Pixel index out of range"));
     192      //THROW(out_of_range("SphereThetaPhi::PIxVal Pixel index out of range"));
     193      cout << " SphereThetaPhi::PIxVal : exceptions a mettre en place" <<endl;
     194      THROW(rangeCheckErr);
     195    }
     196  return pixels_(k);
     197}
     198
     199//++
     200template <class T>
     201T const& SphereThetaPhi<T>::PixVal(int k) const
     202
     203//    Retourne la valeur du contenu du pixel d'indice k
     204//--
     205{
     206  if((k < 0) || (k >= NPix_))
     207    {
    193208      cout << " SphereThetaPhi::PIxVal : exceptions  a mettre en place" <<endl;
    194       THROW(rangeCheckErr);
    195     }
    196   return pixels_(k);
    197 }
    198 
    199 //++
    200 template <class T>
    201 T const& SphereThetaPhi<T>::PixVal(int_4 k) const
    202 
    203 //    Retourne la valeur du contenu du pixel d'indice k
    204 //--
    205 {
    206   if((k < 0) || (k >= NPix_))
    207     {
    208       cout << " SphereThetaPhi::PIxVal : exceptions  a mettre en place" <<endl;
    209       //    THROW(out_of_range("SphereThetaPhi::PIxVal Pixel index out of range"));
     209      //THROW(out_of_range("SphereThetaPhi::PIxVal Pixel index out of range"));
    210210      throw "SphereThetaPhi::PIxVal Pixel index out of range";
    211211    }
    … …  
    216216//++
    217217template <class T>
    218 int_4 SphereThetaPhi<T>::PixIndexSph(r_4 theta, r_4 phi) const
    219 
    220 //    Retourne l'indice du pixel vers lequel pointe une direction définie par
     218int SphereThetaPhi<T>::PixIndexSph(double theta, double phi) const
     219
     220//  Retourne l'indice du pixel vers lequel pointe une direction définie par
    221221//    ses coordonnées sphériques
    222222//--
    223223{
    224   r_4 dphi;
    225   int_4 i,j,k;
     224  double dphi;
     225  int i,j,k;
    226226  bool fgzn = false;
    227227
    228   if( (theta > (r_4)Pi) || (theta < 0. ) ) return(-1);
    229   if( (phi < 0.) || (phi > DeuxPi) ) return(-1);
    230   if( theta > (r_4)Pi*0.5)  { fgzn = true;  theta = Pi-theta; }
     228  if((theta > Pi) || (theta < 0.)) return(-1);
     229  if((phi < 0.) || (phi > DeuxPi)) return(-1);
     230  if(theta > Pi*0.5) {fgzn = true; theta = Pi-theta;}
    231231 
    232232  // La bande d'indice kt est limitée par les valeurs de theta contenues dans
    … …  
    235235    if( theta < Theta_[i] ) break;
    236236 
    237   dphi= (r_4)DeuxPi/(r_4)NPhi_[i-1];
    238  
    239   if (fgzn)  k= NPix_-TNphi_[i]+(int_4)(phi/dphi);
    240   else k= TNphi_[i-1]+(int_4)(phi/dphi);
     237  dphi= DeuxPi/(double)NPhi_[i-1];
     238 
     239  if (fgzn)  k= NPix_-TNphi_[i]+(int)(phi/dphi);
     240  else k= TNphi_[i-1]+(int)(phi/dphi);
    241241  return(k);
    242242}
    … …  
    245245//++
    246246template <class T>
    247 void SphereThetaPhi<T>::PixThetaPhi(int_4 k, r_4& theta, r_4& phi) const
     247void SphereThetaPhi<T>::PixThetaPhi(int k,double& theta,double& phi) const
    248248
    249249//    Retourne les coordonnées (theta,phi) du milieu du pixel d'indice k
    250250//--
    251251{
    252   int_4 i;
     252  int i;
    253253  bool fgzn = false;
    254254 
    255   if( (k < 0) || (k >= NPix_) ) {theta = -99999.; phi = -99999.; return; }
    256   if( k >= NPix_/2)  {fgzn = true;  k = NPix_-1-k; }
     255  if((k < 0) || (k >= NPix_)) {theta = -99999.; phi = -99999.; return; }
     256  if( k >= NPix_/2)  {fgzn = true; k = NPix_-1-k;}
    257257
    258258  // recupère l'indice i de la tranche qui contient le pixel k
    … …  
    266266  // angle phi
    267267  k -= TNphi_[i];
    268   phi= (r_4)DeuxPi/(r_4)NPhi_[i]*(r_4)(k+.5);
     268  phi= DeuxPi/(double)NPhi_[i]*(double)(k+.5);
    269269  if (fgzn) phi= DeuxPi-phi;
    270270}
    … …  
    272272//++
    273273template <class T>
    274 r_8  SphereThetaPhi<T>::PixSolAngle(int_4 dummy) const
     274double SphereThetaPhi<T>::PixSolAngle(int dummy) const
    275275
    276276//    Pixel Solid angle  (steradians)
    … …  
    286286//++
    287287template <class T>
    288 void SphereThetaPhi<T>::Limits(int_4 k, r_4& tetMin, r_4& tetMax, r_4& phiMin, r_4& phiMax)
     288void SphereThetaPhi<T>::Limits(int k,double& tetMin,double& tetMax,double& phiMin,double& phiMax)
    289289
    290290//    Retourne les valeurs de theta et phi limitant le pixel d'indice k
    291291//--
    292292  {
    293   int_4 j;
    294   r_4 dphi;
     293  int j;
     294  double dphi;
    295295  bool fgzn= false;
    296296 
    297   if( (k< 0) || (k >= NPix_) ) {
     297  if((k < 0) || (k >= NPix_)) {
    298298    tetMin= -99999.;
    299299    phiMin= -99999.; 
    … …  
    304304 
    305305  // si on se trouve dans l'hémisphère Sud
    306   if( k >= NPix_/2 ) {
     306  if(k >= NPix_/2) {
    307307    fgzn= true;
    308308    k= NPix_-1-k;
    … …  
    312312  int i;
    313313  for( i=0; i< NTheta_; i++ )
    314     if( k< TNphi_[i+1] ) break;
     314    if(k < TNphi_[i+1]) break;
    315315 
    316316  // valeurs limites de theta dans l'hemisphere Nord
    … …  
    328328  // indice j de discretisation ( phi= j*dphi )
    329329  j= k-TNphi_[i];
    330   dphi= (r_4)DeuxPi/(r_4)NPhi_[i];
     330  dphi= DeuxPi/(double)NPhi_[i];
    331331 
    332332  // valeurs limites de phi
    333   phiMin= dphi*(r_4)(j);
    334   phiMax= dphi*(r_4)(j+1);
     333  phiMin= dphi*(double)(j);
     334  phiMax= dphi*(double)(j+1);
    335335  return;
    336336}
    … …  
    339339//++
    340340template <class T>
    341 int_4 SphereThetaPhi<T>::NbThetaSlices() const
     341int SphereThetaPhi<T>::NbThetaSlices() const
    342342
    343343//    Retourne le nombre de tranches en theta sur la sphere
    344344//--
    345345{
    346   int_4 nbslices;
     346  int nbslices;
    347347  nbslices= 2*NTheta_;
    348348  return(nbslices);
    … …  
    352352//++
    353353template <class T>
    354 int_4 SphereThetaPhi<T>::NPhi(int_4 kt) const
     354int SphereThetaPhi<T>::NPhi(int kt) const
    355355
    356356//    Retourne le nombre de pixels en phi de la tranche kt
    357357//--
    358358{
    359   int_4 nbpix; 
     359  int nbpix; 
    360360  // verification
    361   if( (kt< 0) || (kt>= 2*NTheta_) ) return(-1);
     361  if((kt < 0) || (kt >= 2*NTheta_)) return(-1);
    362362 
    363363  // si on se trouve dans l'hemisphere Sud
    364   if( kt >= NTheta_ ) {
     364  if(kt >= NTheta_) {
    365365    kt= 2*NTheta_-1-kt;
    366366  }
    … …  
    375375//++
    376376template <class T>
    377 void SphereThetaPhi<T>::Theta(int_4 kt,r_4& tetMin,r_4& tetMax)
     377void SphereThetaPhi<T>::Theta(int kt,double& tetMin,double& tetMax)
    378378
    379379//    Retourne les valeurs de theta limitant la tranche kt
    … …  
    407407//++
    408408template <class T>
    409 void SphereThetaPhi<T>::Phi(int_4 kt,int_4 jp,r_4& phiMin,r_4& phiMax)
     409void SphereThetaPhi<T>::Phi(int kt,int jp,double& phiMin,double& phiMax)
    410410
    411411//    Retourne les valeurs de phi limitant le pixel jp de la tranche kt
    … …  
    413413{
    414414  // verification
    415   if( (kt< 0) || (kt>= 2*NTheta_) ) {
     415  if((kt < 0) || (kt >= 2*NTheta_)) {
    416416    phiMin= -99999.;
    417417    phiMax= -99999.;
    … …  
    420420 
    421421  // si on se trouve dans l'hemisphere Sud
    422   if( kt >= NTheta_ ) kt= 2*NTheta_-1-kt;
     422  if(kt >= NTheta_) kt= 2*NTheta_-1-kt;
    423423 
    424424  // verifie si la tranche kt contient au moins jp pixels
    … …  
    429429  }
    430430 
    431   r_4 dphi= (r_4)DeuxPi/(r_4)NPhi_[kt];
    432   phiMin= dphi*(r_4)(jp);
    433   phiMax= dphi*(r_4)(jp+1);
     431  double dphi= DeuxPi/(double)NPhi_[kt];
     432  phiMin= dphi*(double)(jp);
     433  phiMax= dphi*(double)(jp+1);
    434434  return;
    435435}
    … …  
    438438//++
    439439template <class T>
    440 int_4 SphereThetaPhi<T>::Index(int_4 kt,int_4 jp) const
     440int SphereThetaPhi<T>::Index(int kt,int jp) const
    441441
    442442//    Retourne l'indice du pixel d'indice jp dans la tranche kt
    443443//--
    444444{
    445   int_4 k;
     445  int k;
    446446  bool fgzn= false;
    447447 
    448448  // si on se trouve dans l'hemisphere Sud
    449   if( kt >= NTheta_ ) {
     449  if(kt >= NTheta_) {
    450450    fgzn= true;
    451451    kt= 2*NTheta_-1-kt;
    … …  
    471471//++
    472472template <class T>
    473 void SphereThetaPhi<T>::ThetaPhiIndex(int_4 k,int_4& kt,int_4& jp)
     473void SphereThetaPhi<T>::ThetaPhiIndex(int k,int& kt,int& jp)
    474474
    475475//    Retourne les indices kt et jp du pixel d'indice k
    … …  
    478478  bool fgzn= false; 
    479479  // si on se trouve dans l'hemisphere Sud
    480   if( k >= NPix_/2 ) {
    481     fgzn= true;
    482     k= NPix_-1-k;
    483   }
     480  if(k >= NPix_/2)
     481    {
     482      fgzn= true;
     483      k= NPix_-1-k;
     484    }
    484485 
    485486  // on recupere l'indice kt de la tranche qui contient le pixel k
    486487  int i;
    487   for( i=0; i< NTheta_; i++ )
    488     if( k< TNphi_[i+1] ) break;
     488  for(i = 0; i < NTheta_; i++)
     489    if(k < TNphi_[i+1]) break;
    489490 
    490491  // indice  kt de tranche
    … …  
    494495  // indice jp de pixel
    495496  if (fgzn) jp= TNphi_[i+1]-k-1;
    496   else jp= k-TNphi_[i];
    497  
    498 }
    499 //++
    500 template <class T>
    501 void  SphereThetaPhi<T>::Pixelize( int_4 m)
     497  else jp= k-TNphi_[i]; 
     498}
     499//++
     500template <class T>
     501void  SphereThetaPhi<T>::Pixelize(int m)
    502502
    503503//    effectue le découpage en pixels (m et pet ont la même signification
    … …  
    511511//--
    512512{
    513   int_4 ntotpix,i,j;
     513  int ntotpix,i,j;
    514514 
    515515  // Decodage et controle des arguments d'appel :
    … …  
    524524  // Creation des vecteurs contenant :
    525525  // Les valeurs limites de theta (une valeur de plus que le nombre de bandes...)
    526   Theta_= new r_4[m+1];
     526  Theta_= new double[m+1];
    527527 
    528528  // Le nombre de pixels en phi de chacune des bandes en theta
    529   NPhi_ = new int_4[m];
     529  NPhi_ = new int[m];
    530530 
    531531  // Le nombre/Deuxpi total des pixels contenus dans les bandes de z superieur a une
    532532  // bande donnee (mTPphi[m] contient le nombre de pixels total de l'hemisphere)
    533   TNphi_= new int_4[m+1];
     533  TNphi_= new int[m+1];
    534534 
    535535  // Calcul du nombre total de pixels dans chaque bande pour optimiser
    … …  
    544544    {
    545545      TNphi_[j]= TNphi_[j-1]+NPhi_[j-1];
    546       NPhi_[j] = (int_4)(.5+4.*(double)(m-.5)*sin(Pi*(double)j/(double)(2.*m-1.))) ;
     546      NPhi_[j] = (int)(.5+4.*(double)(m-.5)*sin(Pi*(double)j/(double)(2.*m-1.))) ;
    547547    }
    548548 
    … …  
    576576//++
    577577template <class T>
    578 void  SphereThetaPhi<T>::GetThetaSlice(int_4 index, r_4& theta, TVector<float>& phi, TVector<T>& value) const
     578void SphereThetaPhi<T>::GetThetaSlice(int index,double& theta, TVector<double>& phi, TVector<T>& value) const
    579579
    580580//    Retourne, pour la tranche en theta d'indice 'index' le theta
    … …  
    594594    }
    595595
    596   int_4 iring= Index(index,0);
    597   int_4 bid  = this->NPhi(index);
     596  int iring= Index(index,0);
     597  int bid  = this->NPhi(index);
    598598  int lring  = bid;
    599599  cout << " iring= " << iring << " lring= " << lring << endl;
    … …  
    601601  phi.ReSize(lring);
    602602  value.ReSize(lring);
    603   float Te= 0.;
    604   float Fi= 0.;
     603  double Te= 0.;
     604  double Fi= 0.;
    605605  for(int kk = 0; kk < lring; kk++)
    606606    {
    … …  
    613613
    614614template <class T>
    615 void SphereThetaPhi<T>::setmNPhi(int_4* array, int_4 m)
    616   //remplit  le tableau contenant le nombre de pixels en phi de chacune des bandes en theta
     615void SphereThetaPhi<T>::setmNPhi(int* array, int m)
     616  //remplit le tableau contenant le nombre de pixels en phi de chacune des bandes en theta
    617617  //--
    618618{
    619   NPhi_= new int_4[m];
     619  NPhi_= new int[m];
    620620  for(int k = 0; k < m; k++) NPhi_[k]= array[k];
    621621}
    622622
    623623template <class T>
    624 void SphereThetaPhi<T>::setmTNphi(int_4* array, int_4 m)
     624void SphereThetaPhi<T>::setmTNphi(int* array, int m)
    625625  //remplit  le tableau contenant le nombre/Deuxpi total des pixels contenus dans les bandes
    626626  //--
    627627{
    628   TNphi_= new int_4[m];
     628  TNphi_= new int[m];
    629629  for(int k = 0; k < m; k++) TNphi_[k]= array[k];
    630630}
    631631
    632632template <class T>
    633 void SphereThetaPhi<T>::setmTheta(r_4* array, int_4 m)
     633void SphereThetaPhi<T>::setmTheta(double* array, int m)
    634634  //remplit  le tableau contenant les valeurs limites de theta
    635635  //--
    636636{
    637   Theta_= new r_4[m];
     637  Theta_= new double[m];
    638638  for(int k = 0; k < m; k++) Theta_[k]= array[k];
    639639}
    640640
    641641template <class T>
    642 void SphereThetaPhi<T>::setDataBlock(T* data, int_4 m)
     642void SphereThetaPhi<T>::setDataBlock(T* data, int m)
    643643  // remplit  le vecteur des contenus des pixels
    644644{
    … …  
    755755    }
    756756
    757   int_4 mNTheta;
     757  int mNTheta;
    758758  is.GetI4(mNTheta); 
    759759  dobj->setSizeIndex(mNTheta);
    760760
    761   int_4 mNPix;
     761  int mNPix;
    762762  is.GetI4(mNPix);
    763763  dobj->setNbPixels(mNPix);
    764764
    765   r_8 mOmeg;
     765  double mOmeg;
    766766  is.GetR8(mOmeg);
    767767  dobj->setPixSolAngle(mOmeg);
    768768
    769   int_4* mNphi= new int_4[mNTheta];
     769  int* mNphi= new int[mNTheta];
    770770  is.GetI4s(mNphi, mNTheta);
    771771  dobj->setmNPhi(mNphi, mNTheta);
    772772  delete [] mNphi;
    773773
    774   int_4* mTNphi= new int_4[mNTheta+1];
     774  int* mTNphi= new int[mNTheta+1];
    775775  is.GetI4s(mTNphi, mNTheta+1);
    776776  dobj->setmTNphi(mTNphi, mNTheta+1);
    777777  delete [] mTNphi;
    778778
    779   r_4* mTheta= new r_4[mNTheta+1];
    780   is.GetR4s(mTheta, mNTheta+1);
     779  double* mTheta= new double[mNTheta+1];
     780  is.GetR8s(mTheta, mNTheta+1);
    781781  dobj->setmTheta(mTheta, mNTheta+1);
    782782  delete [] mTheta;
    783783
    784784  T* mPixels= new T[mNPix];
    785   //is.Get(mPixels, mNPix);
    786785  PIOSReadArray(is, mPixels, mNPix);
    787786  dobj->setDataBlock(mPixels, mNPix);
    … …  
    801800
    802801  char strg[256];
    803   int_4 mNTheta= dobj->SizeIndex();
    804   int_4 mNPix  = dobj->NbPixels();
     802  int mNTheta= dobj->SizeIndex();
     803  int mNPix  = dobj->NbPixels();
    805804 
    806805  if(dobj->ptrInfo())
    … …  
    821820  os.PutI4s(dobj->getmNPhi() , mNTheta);
    822821  os.PutI4s(dobj->getmTNphi(), mNTheta+1);
    823   os.PutR4s(dobj->getmTheta(), mNTheta+1);
     822  os.PutR8s(dobj->getmTheta(), mNTheta+1);
    824823  //os.Put((dobj->getDataBlock())->Data(), mNPix);
    825824  PIOSWriteArray(os,(dobj->getDataBlock())->Data(), mNPix);
Note: See TracChangeset for help on using the changeset viewer.