| [932] | 1 | // This may look like C code, but it is really -*- C++ -*-
 | 
|---|
 | 2 | //--------------------------------------------------------------------------
 | 
|---|
 | 3 | // File and Version Information:
 | 
|---|
 | 4 | //        Simple Operations on Maps 
 | 
|---|
 | 5 | // Description:
 | 
|---|
 | 6 | //
 | 
|---|
 | 7 | // History (add to end):
 | 
|---|
 | 8 | //      Sophie   April, 2000  - creation
 | 
|---|
 | 9 | //        code from Reza
 | 
|---|
 | 10 | //------------------------------------------------------------------------
 | 
|---|
 | 11 | #include "machdefs.h"
 | 
|---|
 | 12 | #include <math.h>
 | 
|---|
 | 13 | #include "mapoperation.h"
 | 
|---|
 | 14 | /*!
 | 
|---|
| [1371] | 15 |  * \class SOPHYA::MapOperations
 | 
|---|
| [932] | 16 |  * This class is for simple operations on maps
 | 
|---|
 | 17 |  */
 | 
|---|
 | 18 | 
 | 
|---|
 | 19 | template <class T>
 | 
|---|
 | 20 | void MapOperations<T>::sphTosph(PixelMap<T>& in, PixelMap<T>& out)
 | 
|---|
 | 21 |   // Cette fonction remplit la sphere out a partir de la sphere in
 | 
|---|
 | 22 |   // Les spheres peuvent etre de type et de pixelisations differentes
 | 
|---|
 | 23 | {
 | 
|---|
 | 24 |   int kin,kout;
 | 
|---|
 | 25 |   double teta,phi;
 | 
|---|
 | 26 |   
 | 
|---|
 | 27 |   int* cnt = new int[out.NbPixels()+1];
 | 
|---|
 | 28 |   for(kout=0; kout<out.NbPixels(); kout++)
 | 
|---|
 | 29 |     { cnt[kout] = 0; out(kout) = 0.; }
 | 
|---|
 | 30 |   
 | 
|---|
 | 31 |   for(kin=0; kin<in.NbPixels(); kin++) {
 | 
|---|
 | 32 |     in.PixThetaPhi(kin, teta, phi);
 | 
|---|
 | 33 |     kout = out.PixIndexSph(teta, phi);
 | 
|---|
| [1950] | 34 |     if (kout < 0) continue;
 | 
|---|
| [932] | 35 |     out(kout) += in(kin);
 | 
|---|
 | 36 |     cnt[kout] ++;
 | 
|---|
 | 37 |   }
 | 
|---|
 | 38 |   
 | 
|---|
 | 39 |   double moy, sig, dcn;
 | 
|---|
 | 40 |   moy = 0.; sig = 0.;
 | 
|---|
 | 41 |   for(kout=0; kout<out.NbPixels(); kout++) {
 | 
|---|
 | 42 |     dcn = cnt[kout];  moy += dcn;        sig += (dcn*dcn);
 | 
|---|
 | 43 |     if (cnt[kout] > 0)  out(kout) /= dcn;
 | 
|---|
 | 44 |     else {
 | 
|---|
 | 45 |       out.PixThetaPhi(kout, teta, phi);
 | 
|---|
 | 46 |       int pixel = in.PixIndexSph(teta,phi);
 | 
|---|
 | 47 |       out(kout) = in.PixVal(pixel);
 | 
|---|
 | 48 |     }
 | 
|---|
 | 49 |   }
 | 
|---|
 | 50 |   
 | 
|---|
 | 51 |   moy /= out.NbPixels();
 | 
|---|
 | 52 |   sig = sig/out.NbPixels() - moy*moy;
 | 
|---|
 | 53 |   if (sig >= 0.)        sig = sqrt(sig);
 | 
|---|
 | 54 |   
 | 
|---|
 | 55 |   delete[] cnt;
 | 
|---|
 | 56 |   
 | 
|---|
 | 57 | }
 | 
|---|
 | 58 | 
 | 
|---|
 | 59 | template <class T>
 | 
|---|
 | 60 | void  MapOperations<T>::meanSig(NDataBlock<T> const & dbl, double& gmoy, double& gsig)
 | 
|---|
 | 61 | {
 | 
|---|
 | 62 |   gmoy=0.;
 | 
|---|
 | 63 |   gsig = 0.; 
 | 
|---|
 | 64 |   double valok;
 | 
|---|
 | 65 |   for(int k=0; k<(int)dbl.Size(); k++) {
 | 
|---|
 | 66 |     valok = dbl(k);
 | 
|---|
 | 67 |     gmoy += valok;  gsig += valok*valok; 
 | 
|---|
 | 68 |   }
 | 
|---|
 | 69 |   gmoy /= (double)dbl.Size();
 | 
|---|
 | 70 |   gsig = gsig/(double)dbl.Size() - gmoy*gmoy;
 | 
|---|
 | 71 |   if (gsig >= 0.) gsig = sqrt(gsig);
 | 
|---|
 | 72 | }
 | 
|---|
 | 73 | 
 | 
|---|
 | 74 | ///////////////////////////////////////////////////////////////
 | 
|---|
 | 75 | #ifdef __CXX_PRAGMA_TEMPLATES__
 | 
|---|
 | 76 | #pragma define_template MapOperations<r_4>
 | 
|---|
 | 77 | #pragma define_template MapOperations<r_8>
 | 
|---|
 | 78 | #endif
 | 
|---|
 | 79 |   
 | 
|---|
 | 80 | #if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
 | 
|---|
 | 81 | template class MapOperations<r_4>;
 | 
|---|
| [940] | 82 | template class MapOperations<r_8>;
 | 
|---|
| [932] | 83 | #endif
 | 
|---|
 | 84 |    
 | 
|---|