Line data Source code
1 : //# FlagAgentShadow.cc: This file contains the implementation of the FlagAgentShadow class.
2 : //#
3 : //# CASA - Common Astronomy Software Applications (http://casa.nrao.edu/)
4 : //# Copyright (C) Associated Universities, Inc. Washington DC, USA 2011, All rights reserved.
5 : //# Copyright (C) European Southern Observatory, 2011, All rights reserved.
6 : //#
7 : //# This library is free software; you can redistribute it and/or
8 : //# modify it under the terms of the GNU Lesser General Public
9 : //# License as published by the Free software Foundation; either
10 : //# version 2.1 of the License, or (at your option) any later version.
11 : //#
12 : //# This library is distributed in the hope that it will be useful,
13 : //# but WITHOUT ANY WARRANTY, without even the implied warranty of
14 : //# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 : //# Lesser General Public License for more details.
16 : //#
17 : //# You should have received a copy of the GNU Lesser General Public
18 : //# License along with this library; if not, write to the Free Software
19 : //# Foundation, Inc., 59 Temple Place, Suite 330, Boston,
20 : //# MA 02111-1307 USA
21 : //# $Id: $
22 :
23 : #include <flagging/Flagging/FlagAgentShadow.h>
24 :
25 : #include <casacore/measures/Measures/MBaseline.h>
26 : #include <casacore/measures/Measures/MCBaseline.h>
27 : #include <casacore/measures/Measures/MCDirection.h>
28 : #include <casacore/measures/Measures/MCuvw.h>
29 : #include <casacore/measures/Measures/MDirection.h>
30 : #include <casacore/measures/Measures/MEpoch.h>
31 : #include <casacore/measures/Measures/Muvw.h>
32 : #include <casacore/measures/TableMeasures/ScalarMeasColumn.h>
33 :
34 : using namespace casacore;
35 : namespace casa { //# NAMESPACE CASA - BEGIN
36 :
37 : // Definition of static members for common pre-processing
38 : vector<Int> FlagAgentShadow::shadowedAntennas_p;
39 : casa::async::Mutex FlagAgentShadow::staticMembersMutex_p;
40 : vector<bool> FlagAgentShadow::startedProcessing_p;
41 : bool FlagAgentShadow::preProcessingDone_p = false;
42 : uShort FlagAgentShadow::nAgents_p = 0;
43 :
44 16 : FlagAgentShadow::FlagAgentShadow(FlagDataHandler *dh, Record config, Bool writePrivateFlagCube, Bool flag):
45 16 : FlagAgentBase(dh,config,ROWS_PREPROCESS_BUFFER,writePrivateFlagCube,flag)
46 : {
47 16 : setAgentParameters(config);
48 :
49 : // Set preProcessingDone_p static member to false
50 16 : preProcessingDone_p = false;
51 :
52 : // Request loading antenna1,antenna2 and uvw
53 16 : flagDataHandler_p->preLoadColumn(VisBufferComponent2::Antenna1);
54 16 : flagDataHandler_p->preLoadColumn(VisBufferComponent2::Antenna2);
55 16 : flagDataHandler_p->preLoadColumn(VisBufferComponent2::Uvw);
56 : /////flagDataHandler_p->preLoadColumn(VisBufferComponent2::Time);
57 16 : flagDataHandler_p->preLoadColumn(VisBufferComponent2::TimeCentroid);
58 16 : flagDataHandler_p->preLoadColumn(VisBufferComponent2::PhaseCenter);
59 : /////flagDataHandler_p->preLoadColumn(vi::Direction1);
60 :
61 : // FlagAgentShadow counters and ids to handle static variables
62 16 : staticMembersMutex_p.acquirelock();
63 16 : agentNumber_p = nAgents_p;
64 16 : nAgents_p += 1;
65 16 : staticMembersMutex_p.unlock();
66 :
67 : // Set timekeeper to zero - this will later detect when the timestep changes.
68 16 : currTime_p=0.0;
69 :
70 : // Append the supplied additional antennas to COPIES of existing base-class lists.
71 :
72 : // Append to existing lists of antenna info.
73 16 : Int nAntsInMS = flagDataHandler_p->antennaNames_p->nelements();
74 16 : Int nNewAnts=0;
75 :
76 :
77 : // antennaNames_p
78 : // antennaDiameters_p
79 : // antennaPositions_p
80 16 : if( additionalAntennas_p.nfields() )
81 : {
82 :
83 : // For debugging...
84 : //ostringstream recprint;
85 : //additionalAntennas_p.print(recprint);
86 : //cout << " Additional Antennas : " << recprint.str() << endl;
87 :
88 : // TODO : Verify input Record. If invalid, print warning and proceed with no extra antennas.
89 3 : Bool validants=true;
90 3 : String errorReason;
91 15 : for(Int anew=0; anew<(Int) additionalAntennas_p.nfields(); anew++)
92 : {
93 : // Extract the record.
94 24 : Record arec = additionalAntennas_p.subRecord(RecordFieldId(String::toString(anew)));
95 :
96 36 : if( ! arec.isDefined("diameter") ||
97 24 : ( arec.type(arec.fieldNumber("diameter")) != casacore::TpFloat &&
98 24 : arec.type(arec.fieldNumber("diameter")) != casacore::TpInt &&
99 24 : arec.type(arec.fieldNumber("diameter")) != casacore::TpDouble ) )
100 : {
101 0 : validants=false;
102 0 : errorReason += String("Input Record [") + String::toString(anew) + ("] needs a field 'diameter' of type <double> \n");
103 : }
104 :
105 36 : if( ! arec.isDefined("position") ||
106 24 : ( arec.type(arec.fieldNumber("position")) != casacore::TpArrayFloat &&
107 24 : arec.type(arec.fieldNumber("position")) != casacore::TpInt &&
108 24 : arec.type(arec.fieldNumber("position")) != casacore::TpArrayDouble ) )
109 : {
110 0 : validants=false;
111 0 : errorReason += String("Input Record [") + String::toString(anew) + ("] needs a field 'position' of type Array<double>\n");
112 : }
113 : else
114 : {
115 12 : Array<Double> tpos;
116 12 : arec.get( RecordFieldId(String("position")) , tpos );
117 12 : if(tpos.shape() != IPosition(1,3))
118 : {
119 0 : validants=false;
120 0 : errorReason += String("'position' for Record [") + String::toString(anew)+ ("] must be a vector of 3 floats or doubles\n");
121 : }
122 12 : }
123 :
124 12 : }// end of valid-ants loop
125 :
126 : // If antenna list is valid, set the number of new antennas to add.
127 3 : if(validants)
128 : {
129 3 : nNewAnts = additionalAntennas_p.nfields();
130 : }
131 : else // warn and continue.
132 : {
133 0 : *logger_p << LogIO::WARN << "NOT using additional antennas for shadow calculations, for the following reason(s) : " << errorReason << LogIO::POST;
134 : }
135 3 : }// if additionalAnts exist.
136 :
137 :
138 : // Make holders for cumulative information
139 16 : shadowAntennaPositions_p.resize(nAntsInMS+nNewAnts);
140 : /// shadowAntennaNames_p.resize(nAntsInMS+nNewAnts);
141 16 : shadowAntennaDiameters_p.resize(nAntsInMS+nNewAnts);
142 :
143 : // Copy existing antennas into these arrays
144 382 : for(Int antid=0;antid<nAntsInMS;antid++)
145 : {
146 366 : shadowAntennaPositions_p[antid] = flagDataHandler_p->antennaPositions_p->operator()(antid);
147 : ///shadowAntennaNames_p[antid] = flagDataHandler_p->antennaNames_p->operator()(antid);
148 366 : shadowAntennaDiameters_p[antid] = flagDataHandler_p->antennaDiameters_p->operator()(antid);
149 : }
150 :
151 : // If any additional antennas are given, and are valid, add them to the lists
152 28 : for(Int antid=0;antid<nNewAnts;antid++)
153 : {
154 : // Extract the record.
155 24 : Record arec = additionalAntennas_p.subRecord(RecordFieldId(String::toString(antid)));
156 :
157 : // Extract and add new positions
158 12 : Array<Double> aposarr;
159 12 : arec.get( RecordFieldId(String("position")) , aposarr );
160 12 : Vector<Double> aposvec(aposarr);
161 12 : MVPosition apos(aposvec(0),aposvec(1),aposvec(2));
162 12 : shadowAntennaPositions_p[nAntsInMS+antid] = MPosition(apos,MPosition::Types(MPosition::ITRF));
163 :
164 : // Extract and add new diameters
165 : Double adia;
166 12 : arec.get( RecordFieldId(String("diameter")) , adia );
167 12 : shadowAntennaDiameters_p[nAntsInMS+antid] = adia;
168 :
169 : // Extract and add new names
170 : ///String aname aname;
171 : ///arec.get( RecordFieldId(String("name")) , aname );
172 : ///shadowAntennaNames_p[nAntsInMS+antid] = aname;
173 :
174 12 : }
175 :
176 :
177 16 : firststep_p=false; // Set to true, to print a debug message (antenna uvw coordinates for the first row in the first visbuffer seen by this code...
178 :
179 16 : }// end of constructor
180 :
181 32 : FlagAgentShadow::~FlagAgentShadow()
182 : {
183 : // Compiler automagically calls FlagAgentBase::~FlagAgentBase()
184 :
185 : // NOTE: The following is necessary because the static variables
186 : // persist even if all the instances of the class were deleted!
187 16 : staticMembersMutex_p.acquirelock();
188 16 : agentNumber_p = nAgents_p;
189 16 : nAgents_p -= 1;
190 16 : staticMembersMutex_p.unlock();
191 32 : }
192 :
193 : void
194 16 : FlagAgentShadow::setAgentParameters(Record config)
195 : {
196 16 : logger_p->origin(LogOrigin(agentName_p,__FUNCTION__,WHERE));
197 : int exists;
198 :
199 : // Amount of shadowing to allow. Float or Double, in units of Meters.
200 16 : exists = config.fieldNumber ("tolerance");
201 16 : if (exists >= 0)
202 : {
203 14 : if( config.type(exists) != casacore::TpDouble && config.type(exists) != casacore::TpFloat && config.type(exists) != casacore::TpInt )
204 : {
205 0 : throw( AipsError ( "Parameter 'tolerance' must be of type 'double'" ) );
206 : }
207 :
208 14 : shadowTolerance_p = config.asDouble("tolerance");
209 : }
210 : else
211 : {
212 2 : shadowTolerance_p = 0.0;
213 : }
214 :
215 16 : *logger_p << logLevel_p << " tolerance is " << shadowTolerance_p << " meters "<< LogIO::POST;
216 :
217 : // A list of antenna parameters, to add to those in the antenna subtable, to calculate shadows.
218 16 : exists = config.fieldNumber ("addantenna");
219 16 : if (exists >= 0)
220 : {
221 3 : if( config.type(exists) != casacore::TpRecord )
222 : {
223 0 : throw( AipsError ( "Parameter 'addantenna' must be of type 'record/dict'" ) );
224 : }
225 :
226 3 : additionalAntennas_p = config.subRecord( RecordFieldId("addantenna") );
227 : }
228 : else
229 : {
230 13 : additionalAntennas_p = Record();
231 : }
232 :
233 16 : ostringstream recprint;
234 16 : additionalAntennas_p.print(recprint);
235 16 : *logger_p << logLevel_p << " addantenna is " << recprint.str() << LogIO::POST;
236 :
237 32 : return;
238 16 : }
239 :
240 : void
241 2014 : FlagAgentShadow::preProcessBuffer(const vi::VisBuffer2 &visBuffer)
242 : {
243 2014 : if (nAgents_p > 1)
244 : {
245 0 : staticMembersMutex_p.acquirelock();
246 :
247 0 : if (!preProcessingDone_p)
248 : {
249 : // Reset processing state variables
250 0 : if (startedProcessing_p.size() != nAgents_p) startedProcessing_p.resize(nAgents_p,false);
251 0 : for (vector<bool>::iterator iter = startedProcessing_p.begin();iter != startedProcessing_p.end();iter++)
252 : {
253 0 : *iter = false;
254 : }
255 :
256 : // Do actual pre-processing
257 0 : preProcessBufferCore(visBuffer);
258 :
259 : // Mark pre-processing as done so that other agents don't redo it
260 0 : preProcessingDone_p = true;
261 : }
262 :
263 0 : staticMembersMutex_p.unlock();
264 : }
265 : else
266 : {
267 2014 : preProcessBufferCore(visBuffer);
268 : }
269 :
270 2014 : return;
271 : }
272 :
273 : void
274 2014 : FlagAgentShadow::preProcessBufferCore(const vi::VisBuffer2 &/*visBuffer*/)
275 : {
276 : // This function is empty, because shadowedAntennas_p needs to be re-calculated for
277 : // every new timestep, and it is done inside computeRowFlags(), whenever the
278 : // timestep changes.
279 2014 : }
280 :
281 : // (1) Go through all listed baselines for the current timestep, use existing uvw values to
282 : // check for shadowing.
283 : // (2) If not ALL baselines exist in the current timestep, or if additional antennas have been
284 : // supplied, calculate u,v,w, for all antennas, and from there, uvw for all remaining baselines
285 : // and check for shadows between them too.
286 : // Note : The calculation of UVW happens per antenna, not baselines. This is an optimization.
287 : // Note : The direction used for UVW re-calculation is the phasecenter, and not the pointing
288 : // direction of each antenna. This was done to prevent a performance hit due to
289 : // accessing vb.direction1() which accesses MS derived columns, which is also thread-unsafe.
290 : // The only situation where phasecenter is inaccurate, is on-the-fly mosaicing, but
291 : // unless one is doing an on-the-fly mosaic of the whole sky, using a single phase-center (!!!)
292 : // this will not adversely affect shadow flags.
293 2014 : void FlagAgentShadow::calculateShadowedAntennas(const vi::VisBuffer2 &visBuffer, Int rownr)
294 : {
295 : // Init the list of antennas.
296 2014 : shadowedAntennas_p.clear();
297 : Double u,v,w, uvDistance;
298 2014 : Int nAnt = shadowAntennaDiameters_p.nelements();
299 :
300 : // Init the list of baselines, to later decide which to read and which to recalculate.
301 2014 : Vector<Bool> listBaselines(nAnt*(nAnt-1)/2);
302 2014 : listBaselines = false;
303 :
304 : // We know the starting row for this timestep. Find the ending row.
305 : // This assumes that all baselines are grouped together.
306 : // This is guaranteed by the sort-order defined for the visIterator.
307 2014 : Int endrownr = rownr;
308 2014 : Double timeval = visBuffer.timeCentroid()(rownr) ;
309 540164 : for (Int row_i=rownr;row_i<visBuffer.nRows();row_i++)
310 : {
311 538150 : if(timeval < visBuffer.timeCentroid()(row_i)) // we have touched the next timestep
312 : {
313 0 : endrownr = row_i-1;
314 0 : break;
315 : }
316 : else
317 : {
318 538150 : endrownr = row_i;
319 : }
320 : }
321 :
322 : // See CAS-12555 for why this loop which was commented out for long time was
323 : // brought back.
324 : // (1) Now, for all rows between 'rownr' and 'endrownr', calculate shadowed Ants.
325 : // This row range represents all listed baselines in the "current" timestep.
326 : // For those, we take the u,v,w from the UVW column of the MS
327 540164 : for (Int row_i=rownr;row_i<=endrownr;row_i++)
328 : {
329 : // Retrieve antenna ids
330 538150 : auto antenna1 = visBuffer.antenna1()(row_i);
331 538150 : auto antenna2 = visBuffer.antenna2()(row_i);
332 :
333 : // Check if this row corresponds to autocorrelation, or radiometer, sqld, etc.
334 : // (Antennas don't shadow themselves)
335 538150 : if (antenna1 == antenna2) continue;
336 : // Record the baseline being processed
337 538150 : listBaselines[baselineIndex(nAnt,antenna1,antenna2)] = true;
338 :
339 : // Compute uv distance
340 538150 : u = visBuffer.uvw()(0,row_i);
341 538150 : v = visBuffer.uvw()(1,row_i);
342 538150 : w = visBuffer.uvw()(2,row_i);
343 538150 : uvDistance = sqrt(u*u + v*v);
344 :
345 538150 : decideBaselineShadow(uvDistance, w, antenna1, antenna2);
346 :
347 : }// end of for 'row'
348 :
349 :
350 : // (2) Now, if there are any untouched baselines, calculate 'uvw' for all antennas,
351 : // using 'computeAntUVW(), and fill in missing baselines.
352 : // This is the part that picks up invisible antennas, whether they come from the antenna_subtable or
353 : // are externally supplied.
354 2014 : if(product(listBaselines)==false) // could use anyEQ(listBaselines, false)
355 : {
356 : // For the current timestep, compute UVWs for all antennas.
357 : // uvwAnt_p will be filled these values.
358 468 : computeAntUVW(visBuffer, rownr);
359 :
360 : // For all untouched baselines, calculate uvw and check for shadows.
361 13212 : for (Int antenna1=0; antenna1<nAnt; antenna1++)
362 : {
363 12744 : Double u1=uvwAnt_p(0,antenna1), v1=uvwAnt_p(1,antenna1), w1=uvwAnt_p(2,antenna1);
364 192780 : for (Int antenna2=antenna1; antenna2<nAnt; antenna2++)
365 : {
366 : // Check if this row corresponds to autocorrelation (Antennas don't shadow themselves)
367 180036 : if (antenna1 == antenna2) continue;
368 :
369 : // Proceed only if we don't already have this.
370 167292 : if(listBaselines[baselineIndex(nAnt,antenna1,antenna2)] == false)
371 : {
372 52704 : Double u2=uvwAnt_p(0,antenna2), v2=uvwAnt_p(1,antenna2), w2=uvwAnt_p(2,antenna2);
373 :
374 52704 : u = u2-u1;
375 52704 : v = v2-v1;
376 52704 : w = w2-w1;
377 52704 : uvDistance = sqrt(u*u + v*v);
378 :
379 52704 : if(firststep_p==true) // this is only a debug message here....
380 : {
381 0 : cout << "Ant1 : " << antenna1 << " : " << u1 << "," << v1 << "," << w1 << " Ant2 : " << antenna2 << " : " << u2 << "," << v2<< "," << w2 << " UVW : " << u << "," << v << "," << w << endl;
382 : }
383 :
384 52704 : decideBaselineShadow(uvDistance, w, antenna1, antenna2);
385 : }
386 : }
387 : }
388 : }
389 :
390 2014 : firststep_p=false;// debug message should happen only once (at most).
391 :
392 2014 : }// end of calculateShadowedAntennas
393 :
394 705442 : uInt FlagAgentShadow::baselineIndex(uInt nAnt, uInt a1, uInt a2)
395 : {
396 705442 : uInt bindex = (nAnt-1)*nAnt/2 - ((nAnt-1)-a1)*((nAnt-1)-a1+1)/2 + a2-a1-1 ;
397 705442 : AlwaysAssert( bindex < nAnt*(nAnt-1)/2 ,AipsError);
398 705442 : return bindex;
399 : }
400 :
401 :
402 590854 : void FlagAgentShadow::decideBaselineShadow(Double uvDistance, Double w, Int antenna1, Int antenna2)
403 : {
404 : Double antennaDiameter1,antennaDiameter2, antennaDistance;
405 :
406 : // Get antenna diameter
407 590854 : antennaDiameter1 = shadowAntennaDiameters_p[antenna1];
408 590854 : antennaDiameter2 = shadowAntennaDiameters_p[antenna2];
409 :
410 : // Compute effective distance for shadowing
411 590854 : antennaDistance = (antennaDiameter1+antennaDiameter2)/2.0;
412 :
413 : // Check if one of the antennas can be shadowed
414 590854 : if (uvDistance < antennaDistance - shadowTolerance_p)
415 : {
416 : ///////////////////////////////////////////////////////
417 : // Conventions.
418 : // (A) For a Right Handed coordinate system, with 'w' pointing towards the source....
419 : // ( as defined here : http://casa.nrao.edu/Memos/229.html#SECTION00041000000000000000 )
420 : // if(w>0) antenna1 is shadowed by antenna2
421 : // if(w<0) antenna2 is shadowed by antenna1
422 : //
423 : // This is implemented in casapy 3.4 and 4.0
424 : //
425 : // (B) For a Left Handed Coordinate system, with 'w' pointing away from the source...
426 : // ( You get B by flipping the sign on all three axes (u,v,w) of A ).
427 : // This is what is present in the data (i.e. filler, simulator, (our use of Measures?)).
428 : // if(w<0) antenna1 is shadowed by antenna2
429 : // if(w>0) antenna2 is shadowed by antenna1
430 : //
431 : // This is implemented in casapy 4.1 ( from 1 Feb 2013 onwards ).
432 : //
433 : ///////////////////////////////////////////////////////
434 :
435 : // if (w>0) ////// as in casapy 3.4 and casapy 4.0
436 1998 : if (w<0) ////// casapy 4.1 onwards.
437 : {
438 1298 : if (std::find (shadowedAntennas_p.begin(), shadowedAntennas_p.end(), antenna1) == shadowedAntennas_p.end())
439 : {
440 956 : shadowedAntennas_p.push_back(antenna1);
441 : }
442 : }
443 : else
444 : {
445 700 : if (std::find (shadowedAntennas_p.begin(), shadowedAntennas_p.end(), antenna2) == shadowedAntennas_p.end())
446 : {
447 436 : shadowedAntennas_p.push_back(antenna2);
448 : }
449 : }
450 : }
451 590854 : }
452 :
453 : /// NOTE : This function is almost a copy of
454 : /// ms/MeasurementSets/NewMSSimulator::calcAntUVW
455 : /// -- TODO : try to re-use that code by moving out all private-member accesses in the simulator.
456 : /// -- TOCHECK : Should we use vb.timeCentroid() ?? This gives closest results so far, for real and simulated data.
457 : /// NOTE : We are using vb.phasecenter() instead of vb.direction() because of a performance hit
458 : /// and thread-safety problems with vb.direction1().
459 468 : Bool FlagAgentShadow::computeAntUVW(const vi::VisBuffer2 &vb, Int rownr)
460 : {
461 : // Get time and timeinterval from the visbuffer.
462 : Double Time;
463 :
464 : // Centroid gives the closest values to uvws in the MS. For simulated data, gives exact values.
465 468 : Time = vb.timeCentroid()(rownr);
466 :
467 : // Make the Time epoch.
468 936 : MEpoch epoch(Quantity((Time), "s"), MEpoch::UT1);
469 :
470 : // Get the MDirection of the feed of antenna 1. Assume all ants point in the same direction.
471 : //MDirection refdir(vb.direction1()(rownr));
472 468 : MDirection refdir(vb.phaseCenter()); // Each visbuf sees only one fieldId
473 :
474 : // read position of first antenna as reference. Does not matter, since uvws are only differences.
475 468 : MPosition obsPos( shadowAntennaPositions_p[0] );
476 :
477 : // Input measure ref. frame
478 468 : MVPosition basePos=obsPos.getValue();
479 468 : MeasFrame measFrame(obsPos);
480 468 : measFrame.set(epoch);
481 468 : measFrame.set(refdir);
482 :
483 : // Convert direction to J2000
484 468 : MDirection::Convert covertionEngine(refdir,MDirection::Ref(MDirection::J2000,measFrame)); // Each visbuf sees only one fieldId
485 468 : MDirection refdirJ2000 = covertionEngine(refdir);
486 :
487 : // Baseline measure
488 468 : MVBaseline mvbl;
489 468 : MBaseline basMeas;
490 468 : MBaseline::Ref basref(MBaseline::ITRF, measFrame);
491 468 : basMeas.set(mvbl, basref);
492 468 : basMeas.getRefPtr()->set(measFrame);
493 :
494 : // going to convert from ITRF vector to J2000 baseline vector I guess !
495 468 : if(refdirJ2000.getRef().getType() != MDirection::J2000)
496 0 : throw(AipsError("Internal FlagAgentShadow restriction : Conversion to J2000 did not work"));
497 :
498 468 : Int nAnt = shadowAntennaDiameters_p.nelements();
499 468 : if(uvwAnt_p.shape() != IPosition(2,3,nAnt))
500 : {
501 5 : uvwAnt_p.resize(3,nAnt);
502 : }
503 :
504 468 : MBaseline::Convert elconv(basMeas, MBaseline::Ref(MBaseline::J2000));
505 936 : Muvw::Convert uvwconv(Muvw(), Muvw::Ref(Muvw::J2000, measFrame));
506 13212 : for(Int k=0; k< nAnt; ++k)
507 : {
508 12744 : MPosition antpos=shadowAntennaPositions_p(k); // msc.antenna().positionMeas()(k);
509 :
510 12744 : MVBaseline mvblA(obsPos.getValue(), antpos.getValue());
511 12744 : basMeas.set(mvblA, basref);
512 12744 : MBaseline bas2000 = elconv(basMeas);
513 12744 : MVuvw uvw2000 (bas2000.getValue(), refdirJ2000.getValue());
514 12744 : const Vector<double>& xyz = uvw2000.getValue();
515 12744 : uvwAnt_p.column(k)=xyz;
516 12744 : }
517 :
518 468 : return true;
519 468 : }
520 :
521 :
522 : bool
523 538150 : FlagAgentShadow::computeRowFlags(const vi::VisBuffer2 &visBuffer, FlagMapper &/*flags*/, uInt row)
524 : {
525 : // If we have advanced to a new timestep, calculate new antenna UVW values and shadowed antennas
526 : // This function resets and fills 'shadowedAntennas_p'.
527 538150 : if( currTime_p != visBuffer.timeCentroid()(row) )
528 : {
529 2014 : currTime_p = visBuffer.timeCentroid()(row) ;
530 2014 : calculateShadowedAntennas(visBuffer, row);
531 : }
532 :
533 538150 : bool flagRow = false;
534 : // Flag row if either antenna1 or antenna2 are in the list of shadowed antennas
535 538150 : Int antenna1 = visBuffer.antenna1()[row];
536 538150 : Int antenna2 = visBuffer.antenna2()[row];
537 1059969 : if ( (std::find (shadowedAntennas_p.begin(), shadowedAntennas_p.end(), antenna1) != shadowedAntennas_p.end()) or
538 1059969 : (std::find (shadowedAntennas_p.begin(), shadowedAntennas_p.end(), antenna2) != shadowedAntennas_p.end()) )
539 : {
540 20848 : flagRow = true;
541 : }
542 :
543 538150 : if ((nAgents_p > 1) and preProcessingDone_p)
544 : {
545 0 : startedProcessing_p[agentNumber_p] = true;
546 0 : if (std::find (startedProcessing_p.begin(), startedProcessing_p.end(), false) == startedProcessing_p.end())
547 : {
548 0 : preProcessingDone_p = false;
549 : }
550 : }
551 :
552 538150 : return flagRow;
553 : }
554 :
555 :
556 : #if 0
557 :
558 : // // Copy of the old version of this function. It has code for "always recalculate UVW, never recalc UVW
559 : // // and 'decide when to calc UVW'. Above, only the 'decide when to calc UVW' part is used.
560 : // void FlagAgentShadow::calculateShadowedAntennas(const VisBuffer &visBuffer, Int rownr)
561 : // {
562 : // // Init the list of antennas.
563 : // shadowedAntennas_p.clear();
564 : // Double u,v,w, uvDistance;
565 : // Int nAnt = shadowAntennaDiameters_p.nelements();
566 :
567 : // // Init the list of baselines, to later decide which to read and which to recalculate.
568 : // Vector<Bool> listBaselines(nAnt*(nAnt-1)/2);
569 : // listBaselines = false;
570 :
571 : // //uInt countread=0;
572 : // //uInt countcalc=0;
573 : // //Double reftime = 4.794e+09;
574 :
575 : // if (decideUVW_p==true)
576 : // {
577 : // // We know the starting row for this timestep. Find the ending row.
578 : // // This assumes that all baselines are grouped together.
579 : // // This is guaranteed by the sort-order defined for the visIterator.
580 : // Int endrownr = rownr;
581 : // Double timeval = visBuffer.timeCentroid()(rownr) ;
582 : // for (Int row_i=rownr;row_i<visBuffer.nRow();row_i++)
583 : // {
584 : // if(timeval < visBuffer.timeCentroid()(row_i)) // we have touched the next timestep
585 : // {
586 : // endrownr = row_i-1;
587 : // break;
588 : // }
589 : // else
590 : // {
591 : // endrownr = row_i;
592 : // }
593 : // }
594 :
595 : // //cout << "For time : " << timeval-4.73423e+09 << " start : " << rownr << " end : " << endrownr << endl;
596 :
597 : // // Now, for all rows between 'rownr' and 'endrownr', calculate shadowed Ants.
598 : // // This row range represents all listed baselines in the "current" timestep.
599 : // Int antenna1, antenna2;
600 : // for (Int row_i=rownr;row_i<=endrownr;row_i++)
601 : // {
602 : // // Retrieve antenna ids
603 : // antenna1 = visBuffer.antenna1()(row_i);
604 : // antenna2 = visBuffer.antenna2()(row_i);
605 :
606 : // // Check if this row corresponds to autocorrelation (Antennas don't shadow themselves)
607 : // if (antenna1 == antenna2) continue;
608 :
609 : // // Record the baseline being processed
610 : // listBaselines[baselineIndex(nAnt,antenna1,antenna2)] = true;
611 :
612 : // // Compute uv distance
613 : // u = visBuffer.uvw()(row_i)(0);
614 : // v = visBuffer.uvw()(row_i)(1);
615 : // w = visBuffer.uvw()(row_i)(2);
616 : // uvDistance = sqrt(u*u + v*v);
617 :
618 : // //if(row_i==0 && rownr==0) cout << " Row : " << row_i << " uvdist : " << uvDistance << " w : " << w << " time-x : " << visBuffer.timeCentroid()(row_i)-reftime << endl;
619 :
620 : // decideBaselineShadow(uvDistance, w, antenna1, antenna2);
621 : // //countread++;
622 :
623 : // }// end of for 'row'
624 :
625 : // // Now, if there are any untouched baselines, calculate 'uvw' for all antennas, and fill in missing baselines.
626 : // if(product(listBaselines)==false)
627 : // {
628 : // // For the current timestep, compute UVWs for all antennas.
629 : // // uvwAnt_p will be filled these values.
630 : // computeAntUVW(visBuffer, rownr);
631 :
632 : // for (Int antenna1=0; antenna1<nAnt; antenna1++)
633 : // {
634 : // Double u1=uvwAnt_p(0,antenna1), v1=uvwAnt_p(1,antenna1), w1=uvwAnt_p(2,antenna1);
635 : // for (Int antenna2=antenna1; antenna2<nAnt; antenna2++)
636 : // {
637 : // // Check if this row corresponds to autocorrelation (Antennas don't shadow themselves)
638 : // if (antenna1 == antenna2) continue;
639 :
640 : // if(listBaselines[baselineIndex(nAnt,antenna1,antenna2)] == false)
641 : // {
642 : // Double u2=uvwAnt_p(0,antenna2), v2=uvwAnt_p(1,antenna2), w2=uvwAnt_p(2,antenna2);
643 :
644 : // u = u2-u1;
645 : // v = v2-v1;
646 : // w = w2-w1;
647 : // uvDistance = sqrt(u*u + v*v);
648 : // //countcalc++;
649 :
650 : // //if(rownr==0 && antenna1==tant1 && antenna2==tant2) cout << " (r)Row : " << rownr << " uvdist : " << uvDistance << " w : " << w << " time-x : " << visBuffer.timeCentroid()(rownr)-reftime << endl;
651 :
652 : // decideBaselineShadow(uvDistance, w, antenna1, antenna2);
653 :
654 : // listBaselines[baselineIndex(nAnt,antenna1,antenna2)] = true;
655 : // }
656 : // }
657 : // }
658 : // }
659 :
660 : // }
661 : // else if(recalculateUVW_p)
662 : // {
663 : // // (1) For the current timestep, compute UVWs for all antennas.
664 : // // uvwAnt_p will be filled these values.
665 : // computeAntUVW(visBuffer, rownr);
666 :
667 : // // debug code.
668 : // // Int tant1 = visBuffer.antenna1()(rownr);
669 : // // Int tant2 = visBuffer.antenna2()(rownr);
670 :
671 : // // (2) For all antenna pairs, calculate UVW of the baselines, and check for shadowing.
672 : // for (Int antenna1=0; antenna1<nAnt; antenna1++)
673 : // {
674 : // Double u1=uvwAnt_p(0,antenna1), v1=uvwAnt_p(1,antenna1), w1=uvwAnt_p(2,antenna1);
675 : // for (Int antenna2=antenna1; antenna2<nAnt; antenna2++)
676 : // {
677 : // // Check if this row corresponds to autocorrelation (Antennas don't shadow themselves)
678 : // if (antenna1 == antenna2) continue;
679 :
680 : // Double u2=uvwAnt_p(0,antenna2), v2=uvwAnt_p(1,antenna2), w2=uvwAnt_p(2,antenna2);
681 :
682 : // u = u2-u1;
683 : // v = v2-v1;
684 : // w = w2-w1;
685 : // uvDistance = sqrt(u*u + v*v);
686 : // //countcalc++;
687 :
688 : // //if(rownr==0 && antenna1==tant1 && antenna2==tant2) cout << " (r)Row : " << rownr << " uvdist : " << uvDistance << " w : " << w << " time-x : " << visBuffer.timeCentroid()(rownr)-reftime << endl;
689 :
690 : // decideBaselineShadow(uvDistance, w, antenna1, antenna2);
691 :
692 : // }// end for antenna2
693 : // }// end for antenna1
694 :
695 : // }// end of recalculateUVW_p==true
696 : // else // recalculateUVW_p = false
697 : // {
698 :
699 : // // We know the starting row for this timestep. Find the ending row.
700 : // // This assumes that all baselines are grouped together.
701 : // // This is guaranteed by the sort-order defined for the visIterator.
702 : // Int endrownr = rownr;
703 : // Double timeval = visBuffer.timeCentroid()(rownr) ;
704 : // for (Int row_i=rownr;row_i<visBuffer.nRow();row_i++)
705 : // {
706 : // if(timeval < visBuffer.timeCentroid()(row_i)) // we have touched the next timestep
707 : // {
708 : // endrownr = row_i-1;
709 : // break;
710 : // }
711 : // else
712 : // {
713 : // endrownr = row_i;
714 : // }
715 : // }
716 :
717 : // //cout << "For time : " << timeval-4.73423e+09 << " start : " << rownr << " end : " << endrownr << endl;
718 :
719 : // // Now, for all rows between 'rownr' and 'endrownr', calculate shadowed Ants.
720 : // // This row range represents all baselines in the "current" timestep.
721 : // Int antenna1, antenna2;
722 : // for (Int row_i=rownr;row_i<=endrownr;row_i++)
723 : // {
724 : // // Retrieve antenna ids
725 : // antenna1 = visBuffer.antenna1()(row_i);
726 : // antenna2 = visBuffer.antenna2()(row_i);
727 :
728 : // // Check if this row corresponds to autocorrelation (Antennas don't shadow themselves)
729 : // if (antenna1 == antenna2) continue;
730 :
731 : // // Compute uv distance
732 : // u = visBuffer.uvw()(row_i)(0);
733 : // v = visBuffer.uvw()(row_i)(1);
734 : // w = visBuffer.uvw()(row_i)(2);
735 : // uvDistance = sqrt(u*u + v*v);
736 : // //countread++;
737 :
738 : // //if(row_i==0 && rownr==0) cout << " Row : " << row_i << " uvdist : " << uvDistance << " w : " << w << " time-x : " << visBuffer.timeCentroid()(row_i)-reftime << endl;
739 :
740 : // decideBaselineShadow(uvDistance, w, antenna1, antenna2);
741 :
742 : // }// end of for 'row'
743 : // }
744 :
745 : // //cout << "Row : " << rownr << " read : " << countread << " calc : " << countcalc << endl;
746 :
747 : // }// end of calculateShadowedAntennas
748 :
749 :
750 : #endif
751 :
752 :
753 : } //# NAMESPACE CASA - END
754 :
755 :
|
