source: branches/Release2.1/src/RowAccumulator.cpp@ 2222

Last change on this file since 2222 was 979, checked in by mar637, 19 years ago

removed debug cout

File size: 3.3 KB
Line 
1//
2// C++ Implementation: RowAccumulator
3//
4// Description:
5//
6//
7// Author: Malte Marquarding <Malte.Marquarding@csiro.au>, (C) 2005
8//
9// Copyright: See COPYING file that comes with this distribution
10//
11//
12#include <casa/iomanip.h>
13
14#include <casa/Arrays/MaskArrMath.h>
15#include <casa/Arrays/ArrayMath.h>
16#include "RowAccumulator.h"
17
18
19using namespace casa;
20using namespace asap;
21
22RowAccumulator::RowAccumulator(WeightType wt) :
23 weightType_(wt),
24 initialized_(False)
25{
26 reset();
27}
28
29RowAccumulator::~RowAccumulator()
30{
31}
32
33
34void RowAccumulator::add( const Vector< Float >& v,
35 const Vector< Bool >& m,
36 const Vector< Float >& tsys,
37 Double interval,
38 Double time )
39{
40 if (!initialized_) {
41 Vector<Float> dummy(v.nelements(), 0.0);
42 spectrum_.setData(dummy, Vector<Bool>(m.nelements(), True));
43 n_.setData(Vector<Float>(v.nelements(), 0.0), m);
44 weightSum_.setData(Vector<Float>(v.nelements(), 0.0), m);
45 tsysSum_.resize(tsys.nelements()); tsysSum_=0.0;
46 }
47 // add spectrum related weights, so far it is variance only.
48 Float totalweight = 1.0;
49 totalweight *= addTsys(tsys);
50 totalweight *= addInterval(interval);
51 addTime(time);
52 addSpectrum(v, m, totalweight);
53 initialized_ = True;
54}
55
56void RowAccumulator::addSpectrum( const Vector< Float >& v,
57 const Vector< Bool >& m,
58 Float weight)
59{
60 Float totalweight = weight;
61 MaskedArray<Float> data(v,m);
62 if ( weightType_ == asap::VAR ) {
63 if (m.nelements() == userMask_.nelements()) {
64 Float fac = 1.0/variance(data(userMask_));
65 totalweight *= fac;
66 } else {
67 Float fac = 1.0/variance(data);
68 totalweight *= fac;
69 }
70 }
71 data *= totalweight;
72 MaskedArray<Float> wadd(Vector<Float>(m.nelements(),totalweight), m);
73 weightSum_ += wadd;
74 spectrum_ += data;
75 const MaskedArray<Float> inc(Vector<Float>(m.nelements(),1.0), m);
76 n_ += inc;
77}
78
79Float RowAccumulator::addTsys( const casa::Vector< casa::Float > & v )
80{
81 // @fixme this assume tsys is the same for all channels
82
83 Float w = 1.0;
84 tsysSum_ += v[0];
85 if ( weightType_ == asap::TSYS || weightType_ == asap::TINTSYS ) {
86 w /= (v[0]*v[0]);
87 }
88 return w;
89}
90
91void asap::RowAccumulator::addTime( casa::Double t )
92{
93 timeSum_ += t;
94}
95
96Float asap::RowAccumulator::addInterval( casa::Double inter )
97{
98 Float w = 1.0;
99 intervalSum_ += inter;
100 if ( weightType_ == asap::TINT || weightType_ == asap::TINTSYS ) {
101 w /= Float(inter);
102 }
103 return w;
104}
105
106void asap::RowAccumulator::reset( )
107{
108 initialized_ = False;
109 intervalSum_ = 0.0;
110 tsysSum_.resize();
111 timeSum_ = 0.0;
112}
113
114casa::Vector< casa::Float > RowAccumulator::getSpectrum( ) const
115{
116 return (spectrum_/weightSum_).getArray();
117}
118
119casa::Double asap::RowAccumulator::getTime( ) const
120{
121 return timeSum_/max(n_);
122}
123
124casa::Double asap::RowAccumulator::getInterval( ) const
125{
126 return intervalSum_;
127}
128
129casa::Vector< casa::Bool > RowAccumulator::getMask( ) const
130{
131 return spectrum_.getMask();
132}
133
134casa::Vector< casa::Float > asap::RowAccumulator::getTsys( ) const
135{
136 // @fixme this assummes tsys.nelements() == 1
137 return tsysSum_/max(n_);
138}
139
140void asap::RowAccumulator::setUserMask( const casa::Vector< casa::Bool > & m )
141{
142 userMask_.resize();
143 userMask_ = m;
144}
Note: See TracBrowser for help on using the repository browser.