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