Line data Source code
1 : /**
2 : Bojan Nikolic <b.nikolic@mrao.cam.ac.uk>, <bojan@bnikolic.co.uk>
3 : Initial version January 2010.
4 : Maintained by ESO since 2013.
5 :
6 : This file is part of LibAIR and is licensed under GNU Public
7 : License Version 2
8 :
9 : \file arraydata.cpp
10 : Rename arraydata.cc 2023
11 :
12 : Structure to hold data from all WVRs in an array
13 : */
14 :
15 : #include <cmath>
16 :
17 : #include "arraydata.h"
18 :
19 : namespace LibAIR2 {
20 :
21 :
22 50 : InterpArrayData::InterpArrayData(const std::vector<double> &time,
23 : const std::vector<double> &el,
24 : const std::vector<double> &az,
25 : const std::vector<size_t> &state,
26 : const std::vector<size_t> &field,
27 : const std::vector<size_t> &source,
28 50 : size_t nAnts):
29 50 : time(time),
30 50 : el(el),
31 50 : az(az),
32 50 : state(state),
33 50 : field(field),
34 50 : source(source),
35 50 : wvrdata(time.size(),nAnts,4, 0.),
36 50 : nAnts(nAnts)
37 : {
38 50 : }
39 :
40 25 : void InterpArrayData::offsetTime(double dt)
41 : {
42 949 : for(size_t i=0; i<time.size(); ++i)
43 924 : time[i] += dt;
44 25 : }
45 :
46 0 : void interpBadAnt(InterpArrayData &d,
47 : size_t a,
48 : const AntSet &aset)
49 : {
50 0 : size_t n=aset.size();
51 0 : const InterpArrayData::wvrdata_t &data(d.g_wvrdata());
52 0 : for(size_t i=0; i<d.g_time().size(); ++i)
53 : {
54 0 : for(size_t k=0; k < 4; ++k)
55 : {
56 0 : double p=0;
57 0 : for(AntSet::const_iterator j=aset.begin();
58 0 : j!=aset.end();
59 0 : ++j)
60 : {
61 0 : p+=data(i,*j,k);
62 : }
63 0 : d.set(i, a, k, p/n);
64 : }
65 : }
66 0 : }
67 :
68 0 : void interpBadAntW(InterpArrayData &d,
69 : size_t a,
70 : const AntSetWeight &aset)
71 : {
72 0 : const InterpArrayData::wvrdata_t &data(d.g_wvrdata());
73 0 : for(size_t i=0; i<d.g_time().size(); ++i)
74 : {
75 0 : for(size_t k=0; k < 4; ++k)
76 : {
77 0 : double p=0;
78 0 : for(AntSetWeight::const_iterator j=aset.begin();
79 0 : j!=aset.end();
80 0 : ++j)
81 : {
82 0 : p+=data(i,j->second,k)*j->first;
83 : }
84 0 : d.set(i, a, k, p);
85 : }
86 : }
87 0 : }
88 :
89 25 : InterpArrayData *filterState(InterpArrayData &d,
90 : const std::set<size_t>& states)
91 : {
92 25 : std::vector<double> time;
93 25 : std::vector<double> el;
94 25 : std::vector<double> az;
95 25 : std::vector<size_t> state;
96 25 : std::vector<size_t> field;
97 25 : std::vector<size_t> src;
98 949 : for(size_t i=0; i<d.g_state().size(); ++i)
99 : {
100 924 : if(states.count(d.g_state()[i]))
101 : {
102 924 : time.push_back(d.g_time()[i]);
103 924 : el.push_back(d.g_el()[i]);
104 924 : az.push_back(d.g_az()[i]);
105 924 : state.push_back(d.g_state()[i]);
106 924 : field.push_back(d.g_field()[i]);
107 924 : src.push_back(d.g_source()[i]);
108 : }
109 : }
110 : std::unique_ptr<InterpArrayData> res(new InterpArrayData(time,
111 : el,
112 : az,
113 : state,
114 : field,
115 : src,
116 25 : d.nAnts));
117 :
118 25 : size_t n=0;
119 949 : for(size_t i=0; i<d.g_state().size(); ++i)
120 924 : if(states.count(d.g_state()[i]))
121 : {
122 19746 : for (size_t j=0; j<d.nAnts; ++j)
123 : {
124 94110 : for(size_t k=0; k<4; ++k)
125 : {
126 75288 : res->set(n, j, k,
127 75288 : d.g_wvrdata()(i,j,k));
128 : }
129 : }
130 924 : ++n;
131 : }
132 50 : return res.release();
133 25 : }
134 :
135 0 : void smoothWVR(InterpArrayData &d,
136 : size_t nsample)
137 : {
138 0 : if (nsample<2)
139 0 : return;
140 :
141 0 : const size_t N=d.g_time().size();
142 0 : const size_t delt=(nsample)/2;
143 0 : const bool even=((nsample%2)==0);
144 :
145 0 : size_t starti=0;
146 0 : while (starti < N )
147 : {
148 : // Figure out where the end of this source/state observation is
149 0 : size_t src=d.g_source()[starti];
150 0 : size_t state=d.g_state()[starti];
151 0 : size_t endi=starti;
152 0 : while(endi<N and d.g_source()[endi] == src and d.g_state()[endi] == state)
153 0 : ++endi;
154 :
155 : //go through each wvr channel
156 0 : for(size_t k=0; k<4; ++k)
157 : {
158 : //go through each wvr.
159 0 : for (size_t l=0; l<d.nAnts; ++l)
160 : {
161 : //go through each sample in time
162 0 : for(size_t i=starti; i<endi; ++i)
163 : {
164 : // Check if we are on edge
165 0 : if (i< starti+delt or i>endi - delt-1)
166 0 : continue;
167 :
168 0 : double sum=0;
169 : // Go through elements to average over for new sample i
170 0 : for(size_t j=i-delt; j<i+delt+1; ++j)
171 : {
172 0 : if ( even and (j==i-delt or j==i+delt))
173 : {
174 0 : sum=sum+(d.g_wvrdata()(j,l,k)*0.5);
175 : }
176 : else
177 : {
178 0 : sum=sum+(d.g_wvrdata()(j,l,k));
179 : }
180 : }
181 0 : d.set(i, l, k, sum/(float(nsample)));
182 :
183 : }
184 : }
185 : }
186 0 : starti=endi;
187 : }
188 : }
189 :
190 : }
191 :
192 :
|