Line data Source code
1 : //# ComponentList.cc: this defines the ComponentList implementation
2 : //# Copyright (C) 1996,1997,1998,1999,2000,2001,2002
3 : //# Associated Universities, Inc. Washington DC, USA.
4 : //#
5 : //# This library is free software; you can redistribute it and/or modify it
6 : //# under the terms of the GNU Library General Public License as published by
7 : //# the Free Software Foundation; either version 2 of the License, or (at your
8 : //# option) any later version.
9 : //#
10 : //# This library is distributed in the hope that it will be useful, but WITHOUT
11 : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 : //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
13 : //# License for more details.
14 : //#
15 : //# You should have received a copy of the GNU Library General Public License
16 : //# along with this library; if not, write to the Free Software Foundation,
17 : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
18 : //#
19 : //# Correspondence concerning AIPS++ should be addressed as follows:
20 : //# Internet email: aips2-request@nrao.edu.
21 : //# Postal address: AIPS++ Project Office
22 : //# National Radio Astronomy Observatory
23 : //# 520 Edgemont Road
24 : //# Charlottesville, VA 22903-2475 USA
25 : //#
26 : //# $Id: ComponentList.cc 21229 2012-04-02 12:00:20Z gervandiepen $
27 :
28 : #include <components/ComponentModels/ComponentList.h>
29 : #include <components/ComponentModels/ComponentType.h>
30 : #include <components/ComponentModels/Flux.h>
31 : #include <casacore/tables/Tables/ScalarColumn.h>
32 : #include <casacore/tables/Tables/ScaColDesc.h>
33 : #include <casacore/tables/Tables/ScaRecordColDesc.h>
34 : #include <components/ComponentModels/ComponentShape.h>
35 : #include <components/ComponentModels/SpectralModel.h>
36 : #include <casacore/measures/TableMeasures/ScalarMeasColumn.h>
37 : #include <casacore/measures/TableMeasures/ArrayQuantColumn.h>
38 : #include <casacore/measures/TableMeasures/ScalarQuantColumn.h>
39 : #include <casacore/measures/TableMeasures/TableMeasDesc.h>
40 : #include <casacore/measures/TableMeasures/TableMeasRefDesc.h>
41 : #include <casacore/measures/TableMeasures/TableMeasValueDesc.h>
42 : #include <casacore/measures/TableMeasures/TableQuantumDesc.h>
43 : #include <casacore/casa/Arrays/Array.h>
44 : #include <casacore/casa/Arrays/ArrayMath.h>
45 : #include <casacore/casa/Arrays/ArrayLogical.h>
46 : #include <casacore/casa/Arrays/Vector.h>
47 : #include <casacore/casa/Arrays/Cube.h>
48 : #include <casacore/casa/Exceptions/Error.h>
49 : #include <casacore/casa/Logging/LogIO.h>
50 : #include <casacore/casa/BasicSL/Constants.h>
51 : #include <casacore/casa/BasicSL/Complex.h>
52 : #include <casacore/casa/BasicMath/Math.h>
53 : #include <casacore/measures/Measures/MDirection.h>
54 : #include <casacore/measures/Measures/MFrequency.h>
55 : #include <casacore/measures/Measures/MCDirection.h>
56 : #include <casacore/measures/Measures/MCFrequency.h>
57 : #include <casacore/casa/Quanta/MVFrequency.h>
58 : #include <casacore/casa/Quanta/MVAngle.h>
59 : #include <casacore/casa/Quanta/MVDirection.h>
60 : #include <casacore/measures/Measures/MeasConvert.h>
61 : #include <casacore/casa/Quanta/Quantum.h>
62 : #include <casacore/casa/Quanta/Unit.h>
63 : #include <casacore/casa/OS/Path.h>
64 : #include <casacore/tables/Tables/ArrColDesc.h>
65 : #include <casacore/tables/Tables/ArrayColumn.h>
66 : #include <casacore/tables/Tables/ColumnDesc.h>
67 : #include <casacore/tables/Tables/ScaColDesc.h>
68 : #include <casacore/tables/Tables/ScalarColumn.h>
69 : #include <casacore/tables/Tables/SetupNewTab.h>
70 : #include <casacore/tables/Tables/TableDesc.h>
71 : #include <casacore/tables/Tables/ColDescSet.h>
72 : #include <casacore/tables/Tables/TableLock.h>
73 : #include <casacore/tables/Tables/TableRecord.h>
74 : #include <casacore/tables/DataMan/TiledCellStMan.h>
75 : #include <casacore/casa/Utilities/Assert.h>
76 : #include <casacore/casa/Utilities/GenSort.h>
77 : #include <casacore/casa/Utilities/Sort.h>
78 : #include <casacore/casa/BasicSL/String.h>
79 : #include <casacore/casa/Containers/Record.h>
80 :
81 : using namespace casacore;
82 : namespace casa { //# NAMESPACE CASA - BEGIN
83 :
84 : const String fluxName = "Flux";
85 : const String fluxUnitName = "Flux_Unit";
86 : const String fluxPolName = "Flux_Polarisation";
87 : const String fluxErrName = "Flux_Error";
88 : const String shapeName = "Shape";
89 : const String refDirName = "Reference_Direction";
90 : const String dirErrName = "Direction_Error";
91 : const String dirErrUnitName = "Direction_Error_Units";
92 : const String shapeParName = "Shape_Parameters";
93 : const String shapeErrName = "Shape_Error";
94 : const String spectrumName = "Spectrum_Shape";
95 : const String refFreqName = "Reference_Frequency";
96 : const String freqErrName = "Frequency_Error";
97 : const String spectralRecordName = "Spectral_Record";
98 : const String freqErrUnitName = "Frequency_Error_Units";
99 : const String spectParName = "Spectral_Parameters";
100 : const String spectErrName = "Spectral_Error";
101 : const String labelName = "Label";
102 : const String optParColName = "Optional_Parameters";
103 :
104 :
105 6538 : ComponentList::ComponentList()
106 : :itsList(),
107 : itsNelements(0),
108 : itsTable(),
109 : itsROFlag(false),
110 : itsSelectedFlags(),
111 : itsOrder(),
112 : itsAddOptCol(false),
113 6538 : itsRewriteTable(True)
114 : {
115 6538 : AlwaysAssert(ok(), AipsError);
116 6538 : }
117 :
118 146 : ComponentList::ComponentList(const Path& fileName, Bool readOnly, Bool rewriteTable)
119 : :itsList(),
120 : itsNelements(0),
121 : itsTable(),
122 : itsROFlag(false),
123 : itsSelectedFlags(),
124 : itsOrder(),
125 : itsAddOptCol(false),
126 146 : itsRewriteTable(rewriteTable)
127 : {
128 146 : readTable(fileName, readOnly);
129 146 : AlwaysAssert(ok(), AipsError);
130 146 : }
131 :
132 165 : ComponentList::ComponentList(const ComponentList& other)
133 165 : :itsList(other.itsList),
134 165 : itsNelements(other.itsNelements),
135 165 : itsTable(other.itsTable),
136 165 : itsROFlag(other.itsROFlag),
137 165 : itsSelectedFlags(other.itsSelectedFlags),
138 165 : itsOrder(other.itsOrder),
139 165 : itsAddOptCol(other.itsAddOptCol),
140 165 : itsRewriteTable(other.itsRewriteTable)
141 : {
142 165 : DebugAssert(ok(), AipsError);
143 165 : }
144 :
145 6827 : ComponentList::~ComponentList() {
146 6827 : if (! itsROFlag && ! itsTable.isNull() && itsRewriteTable) {
147 195 : writeTable();
148 : }
149 6827 : AlwaysAssert(ok(), AipsError);
150 6827 : }
151 :
152 1872 : ComponentList& ComponentList::operator=(const ComponentList& other){
153 1872 : if (this != &other) {
154 1872 : if (! itsROFlag && ! itsTable.isNull() && itsRewriteTable) {
155 0 : writeTable();
156 : }
157 1872 : itsList = other.itsList;
158 1872 : itsNelements = other.itsNelements;
159 1872 : itsTable = other.itsTable;
160 1872 : itsROFlag = other.itsROFlag;
161 1872 : itsSelectedFlags = other.itsSelectedFlags;
162 1872 : itsOrder = other.itsOrder;
163 1872 : itsAddOptCol = other.itsAddOptCol;
164 1872 : itsRewriteTable = other.itsRewriteTable;
165 : }
166 1872 : DebugAssert(ok(), AipsError);
167 1872 : return *this;
168 : }
169 :
170 0 : Bool ComponentList::isPhysical(const Vector<Int>& indices) const {
171 0 : DebugAssert(ok(), AipsError);
172 : // The static_casts are a workaround for an SGI compiler bug
173 0 : DebugAssert(allGE(static_cast<const Vector<Int> &>(indices), 0), AipsError);
174 0 : DebugAssert(allLT(static_cast<const Vector<Int> &>(indices),
175 : static_cast<Int>(nelements())), AipsError);
176 0 : Bool retVal = true;
177 0 : uInt c = indices.nelements();
178 0 : while (retVal && c > 0) {
179 0 : c--;
180 0 : retVal = itsList[indices(c)].isPhysical();
181 : }
182 0 : return retVal;
183 : }
184 :
185 0 : Flux<Double> ComponentList::sample(const MDirection& sampleDir,
186 : const MVAngle& pixelLatSize,
187 : const MVAngle& pixelLongSize,
188 : const MFrequency& centerFreq) const {
189 0 : DebugAssert(ok(), AipsError);
190 0 : const Unit retUnit("Jy");
191 0 : const ComponentType::Polarisation retPol(ComponentType::STOKES);
192 0 : Vector<DComplex> result(4, DComplex(0,0));
193 0 : Flux<Double> compFlux;
194 0 : for (uInt i = 0; i < nelements(); i++) {
195 0 : compFlux = component(i).sample(sampleDir, pixelLatSize, pixelLongSize,
196 0 : centerFreq);
197 0 : compFlux.convertUnit(retUnit);
198 0 : compFlux.convertPol(retPol);
199 0 : result += compFlux.value();
200 : }
201 0 : return Flux<Double>(result, retPol);
202 : }
203 :
204 983 : void ComponentList::sample(Cube<Double>& samples,
205 : const Unit& reqUnit,
206 : const Vector<MVDirection>& directions,
207 : const MeasRef<MDirection>& dirRef,
208 : const MVAngle& pixelLatSize,
209 : const MVAngle& pixelLongSize,
210 : const Vector<MVFrequency>& frequencies,
211 : const MeasRef<MFrequency>& freqRef) const {
212 983 : samples = 0.0;
213 5421 : for (uInt i = 0; i < nelements(); i++) {
214 4438 : component(i).sample(samples, reqUnit, directions, dirRef,
215 4438 : pixelLatSize, pixelLongSize, frequencies, freqRef);
216 : }
217 983 : }
218 :
219 2292 : void ComponentList::add(SkyComponent component) {
220 : // AlwaysAssert(itsROFlag == false, AipsError);
221 2292 : DebugAssert(ok(), AipsError);
222 2292 : uInt blockSize = itsList.nelements();
223 2292 : if (itsNelements == blockSize) {
224 1970 : const uInt newSize
225 1970 : = (blockSize < 50) ? 2 * blockSize + 1 : blockSize + 100;
226 1970 : itsList.resize(newSize);
227 1970 : itsSelectedFlags.resize(newSize);
228 1970 : itsOrder.resize(newSize);
229 : }
230 : // for limb-darkened disk shape, add an optional col
231 2292 : if (component.shape().type()==ComponentType::LDISK) {
232 0 : itsAddOptCol=true;
233 : }
234 2292 : itsList[itsNelements] = component;
235 2292 : itsSelectedFlags[itsNelements] = false;
236 2292 : itsOrder[itsNelements] = itsNelements;
237 2292 : itsNelements++;
238 2292 : }
239 :
240 730 : void ComponentList::addList(const ComponentList& list) {
241 1278 : for (uInt i=0; i<list.nelements(); i++) {
242 548 : add(list.component(i));
243 : }
244 730 : }
245 :
246 :
247 1 : void ComponentList::remove(const uInt& index) {
248 : // AlwaysAssert(itsROFlag == false, AipsError);
249 1 : AlwaysAssert(index < nelements(), AipsError);
250 1 : DebugAssert(ok(), AipsError);
251 1 : uInt realIndex = itsOrder[index];
252 1 : itsSelectedFlags.remove(realIndex, false);
253 1 : itsList.remove(realIndex, false);
254 1 : itsOrder.remove(index, false);
255 1 : itsNelements--;
256 1 : for (uInt i = 0; i < nelements(); i++) {
257 0 : if (itsOrder[i] > realIndex) {
258 0 : itsOrder[i]--;
259 : }
260 : }
261 1 : }
262 :
263 1 : void ComponentList::remove(const Vector<Int>& indices) {
264 2 : Vector<Int> zeroCheck(indices);
265 1 : AlwaysAssert(allGE(zeroCheck, 0), AipsError);
266 1 : uInt c = indices.nelements();
267 2 : Vector<uInt> uIndices(c);
268 1 : convertArray(uIndices, indices);
269 1 : GenSort<uInt>::sort(uIndices);
270 2 : while (c != 0) {
271 1 : c--;
272 1 : remove(uIndices(c));
273 : }
274 1 : }
275 :
276 2305505 : uInt ComponentList::nelements() const {
277 2305505 : return itsNelements;
278 : }
279 :
280 37 : uInt ComponentList::size() const {
281 37 : return itsNelements;
282 : }
283 :
284 0 : void ComponentList::deselect(const Vector<Int>& indexes) {
285 0 : for (uInt i = 0; i < indexes.nelements(); i++) {
286 0 : AlwaysAssert(indexes(i) < Int(nelements()), AipsError);
287 0 : AlwaysAssert(indexes(i) >= 0, AipsError);
288 0 : itsSelectedFlags[itsOrder[indexes(i)]] = false;
289 : }
290 0 : DebugAssert(ok(), AipsError);
291 0 : }
292 :
293 0 : void ComponentList::select(const Vector<Int>& indexes) {
294 0 : for (uInt i = 0; i < indexes.nelements(); i++) {
295 0 : AlwaysAssert(indexes(i) < Int(nelements()), AipsError);
296 0 : AlwaysAssert(indexes(i) >= 0, AipsError);
297 0 : itsSelectedFlags[itsOrder[indexes(i)]] = true;
298 : }
299 0 : DebugAssert(ok(), AipsError);
300 0 : }
301 :
302 0 : Vector<Int> ComponentList::selected() const {
303 0 : DebugAssert(ok(), AipsError);
304 0 : uInt nSelected = 0;
305 0 : for (uInt i = 0; i < nelements(); i++) {
306 0 : if (itsSelectedFlags[i] == true) {
307 0 : nSelected++;
308 : }
309 : }
310 0 : Vector<Int> retVal(nSelected);
311 0 : uInt s = 0;
312 0 : for (uInt j = 0; j < nelements(); j++) {
313 0 : if (itsSelectedFlags[j] == true) {
314 0 : retVal(s) = j;
315 0 : s++;
316 : }
317 : }
318 0 : return retVal;
319 : }
320 :
321 37 : void ComponentList::setLabel(const Vector<Int>& which,
322 : const String& newLabel) {
323 : uInt c;
324 74 : for (uInt i = 0; i < which.nelements(); i++) {
325 37 : AlwaysAssert(which(i) >= 0, AipsError);
326 37 : c = which(i);
327 37 : component(c).label() = newLabel;
328 : }
329 37 : DebugAssert(ok(), AipsError);
330 37 : }
331 :
332 393 : void ComponentList::getFlux(Vector<Quantity>& fluxQuant, const Int& which) const {
333 786 : SkyComponent comp = component(which);
334 : // each element in the returned vector represents a different polarization.
335 : // NumericTraits::Conjugate is just a confusing way of saying Complex if you
336 : // look at how comp.flux().value() is implemented.
337 786 : Vector<Complex> flux(comp.flux().value().nelements());
338 393 : convertArray(flux,comp.flux().value());
339 786 : Unit unit = comp.flux().unit();
340 393 : fluxQuant.resize(flux.nelements());
341 1965 : for (uInt i=0; i<flux.nelements(); ++i) {
342 1572 : fluxQuant[i] = Quantity(real(flux[i]), unit);
343 : }
344 393 : }
345 :
346 0 : void ComponentList::getFlux(Vector<Quantum<Complex> >& fluxQuant, const Int& which) {
347 0 : SkyComponent comp = component(which);
348 0 : Vector<Complex> flux(comp.flux().value().nelements());
349 0 : convertArray(flux,comp.flux().value());
350 0 : Unit unit = comp.flux().unit();
351 0 : fluxQuant.resize(flux.nelements());
352 0 : for (uInt i=0; i<flux.nelements(); ++i) {
353 0 : fluxQuant[i] = Quantum<Complex>(flux[i], unit);
354 : }
355 0 : }
356 :
357 :
358 :
359 37 : void ComponentList::setFlux(const Vector<Int>& which,
360 : const Flux<Double>& newFlux) {
361 : uInt c;
362 74 : for (uInt i = 0; i < which.nelements(); i++) {
363 37 : AlwaysAssert(which(i) >= 0, AipsError);
364 37 : c = which(i);
365 37 : component(c).flux() = newFlux;
366 : }
367 37 : DebugAssert(ok(), AipsError);
368 37 : }
369 :
370 786 : Vector<String> ComponentList::getStokes(const Int& which) const {
371 1572 : SkyComponent comp = component(which);
372 786 : ComponentType::Polarisation stokesType = comp.flux().pol();
373 786 : Vector<String> polarization(4);
374 : // the polarization determination logic needs to be refactored into
375 : // a method in a more appropriate class
376 786 : if (stokesType == ComponentType::STOKES) {
377 786 : polarization[0] = "I";
378 786 : polarization[1] = "Q";
379 786 : polarization[2] = "U";
380 786 : polarization[3] = "V";
381 : }
382 0 : else if (stokesType == ComponentType::LINEAR) {
383 0 : polarization[0] = "XX";
384 0 : polarization[1] = "XY";
385 0 : polarization[2] = "YX";
386 0 : polarization[3] = "YY";
387 : }
388 0 : else if (stokesType == ComponentType::CIRCULAR) {
389 0 : polarization[0] = "RR";
390 0 : polarization[1] = "RL";
391 0 : polarization[2] = "LR";
392 0 : polarization[3] = "LL";
393 : }
394 : else {
395 0 : polarization.set("UNKNOWN");
396 : }
397 1572 : return polarization;
398 : }
399 :
400 0 : void ComponentList::convertFluxUnit(const Vector<Int>& which,
401 : const Unit& unit) {
402 : uInt c;
403 0 : for (uInt i = 0; i < which.nelements(); i++) {
404 0 : AlwaysAssert(which(i) >= 0, AipsError);
405 0 : c = which(i);
406 0 : component(c).flux().convertUnit(unit);
407 : }
408 0 : DebugAssert(ok(), AipsError);
409 0 : }
410 :
411 0 : void ComponentList::convertFluxPol(const Vector<Int>& which,
412 : ComponentType::Polarisation pol) {
413 : uInt c;
414 0 : for (uInt i = 0; i < which.nelements(); i++) {
415 0 : AlwaysAssert(which(i) >= 0, AipsError);
416 0 : c = which(i);
417 0 : component(c).flux().convertPol(pol);
418 : }
419 0 : DebugAssert(ok(), AipsError);
420 0 : }
421 :
422 37 : void ComponentList::setRefDirection(const Vector<Int>& which,
423 : const MVDirection& newDir) {
424 : uInt c;
425 74 : MDirection curDir;
426 74 : for (uInt i = 0; i < which.nelements(); i++) {
427 37 : AlwaysAssert(which(i) >= 0, AipsError);
428 37 : c = which(i);
429 37 : ComponentShape& curShape = component(c).shape();
430 37 : curDir = curShape.refDirection();
431 37 : curDir.set(newDir);
432 37 : curShape.setRefDirection(curDir);
433 : }
434 37 : DebugAssert(ok(), AipsError);
435 37 : }
436 :
437 37 : void ComponentList::setRefDirectionFrame(const Vector<Int>& which,
438 : MDirection::Types newFrame) {
439 : uInt c;
440 74 : MDirection curDir;
441 74 : const MDirection::Ref newRef(newFrame);
442 74 : for (uInt i = 0; i < which.nelements(); i++) {
443 37 : AlwaysAssert(which(i) >= 0, AipsError);
444 37 : c = which(i);
445 37 : ComponentShape& curShape = component(c).shape();
446 37 : curDir = curShape.refDirection();
447 37 : curDir.set(newRef);
448 37 : curShape.setRefDirection(curDir);
449 : }
450 37 : DebugAssert(ok(), AipsError);
451 37 : }
452 :
453 0 : void ComponentList::convertRefDirection(const Vector<Int>& which,
454 : MDirection::Types newFrame) {
455 : uInt c;
456 0 : MDirection::Convert converter;
457 0 : converter.setOut(newFrame);
458 0 : MDirection curDir;
459 0 : for (uInt i = 0; i < which.nelements(); i++) {
460 0 : AlwaysAssert(which(i) >= 0, AipsError);
461 0 : c = which(i);
462 0 : ComponentShape& curShape = component(c).shape();
463 0 : curDir = curShape.refDirection();
464 0 : curShape.setRefDirection(converter(curDir));
465 : }
466 0 : DebugAssert(ok(), AipsError);
467 0 : }
468 :
469 7 : MDirection ComponentList::getRefDirection(Int which) const {
470 7 : const ComponentShape& compShape = component(which).shape();
471 7 : MDirection refDir = compShape.refDirection();
472 7 : return refDir;
473 : }
474 :
475 548 : void ComponentList::setShape(const Vector<Int>& which,
476 : const ComponentShape& newShape) {
477 : uInt c;
478 1096 : for (uInt i = 0; i < which.nelements(); i++) {
479 548 : AlwaysAssert(which(i) >= 0, AipsError);
480 548 : c = which(i);
481 548 : component(c).setShape(newShape);
482 : //for limb-darkened disk shape
483 548 : if (newShape.type()==ComponentType::LDISK) {
484 0 : itsAddOptCol=true;
485 : }
486 : }
487 548 : DebugAssert(ok(), AipsError);
488 548 : }
489 :
490 1497 : const ComponentShape* ComponentList::getShape(Int which) const {
491 1497 : return component(which).shape().getPtr();
492 : }
493 :
494 :
495 37 : void ComponentList::setShapeParms(const Vector<Int>& which,
496 : const ComponentShape& newShape) {
497 : uInt c;
498 74 : MDirection oldDir;
499 74 : for (uInt i = 0; i < which.nelements(); i++) {
500 37 : AlwaysAssert(which(i) >= 0, AipsError);
501 37 : c = which(i);
502 37 : SkyComponent& comp = component(c);
503 37 : oldDir = comp.shape().refDirection();
504 37 : component(c).setShape(newShape);
505 37 : comp.shape().setRefDirection(oldDir);
506 : }
507 37 : DebugAssert(ok(), AipsError);
508 37 : }
509 :
510 0 : void ComponentList::setOptParms(const Vector<Int>& which,
511 : const ComponentShape& newShape) {
512 : uInt c;
513 0 : Vector<Double> optparms;
514 0 : for (uInt i = 0; i < which.nelements(); i++) {
515 0 : AlwaysAssert(which(i) >= 0, AipsError);
516 0 : c = which(i);
517 0 : SkyComponent& comp = component(c);
518 0 : component(c).setShape(newShape);
519 0 : if (comp.shape().type()==ComponentType::LDISK) {
520 0 : optparms = comp.shape().optParameters();
521 : //comp.shape().setOptParameters(optparms);
522 0 : component(c).optionalParameters()=optparms;
523 0 : itsAddOptCol=true;
524 : }
525 : }
526 0 : DebugAssert(ok(), AipsError);
527 0 : }
528 :
529 0 : void ComponentList::setSpectrum(const Vector<Int>& which,
530 : const SpectralModel& newSpectrum) {
531 : uInt c;
532 0 : for (uInt i = 0; i < which.nelements(); i++) {
533 0 : AlwaysAssert(which(i) >= 0, AipsError);
534 0 : c = which(i);
535 0 : component(c).setSpectrum(newSpectrum);
536 : }
537 0 : DebugAssert(ok(), AipsError);
538 0 : }
539 :
540 54 : void ComponentList::setSpectrumParms(const Vector<Int>& which,
541 : const SpectralModel& newSpectrum) {
542 : uInt c;
543 108 : MFrequency oldFreq;
544 108 : for (uInt i = 0; i < which.nelements(); i++) {
545 54 : AlwaysAssert(which(i) >= 0, AipsError);
546 54 : c = which(i);
547 54 : SkyComponent& comp = component(c);
548 54 : oldFreq = comp.spectrum().refFrequency();
549 54 : component(c).setSpectrum(newSpectrum);
550 54 : comp.spectrum().setRefFrequency(oldFreq);
551 : }
552 54 : DebugAssert(ok(), AipsError);
553 54 : }
554 :
555 116 : void ComponentList::setRefFrequency(const Vector<Int>& which,
556 : const MVFrequency& newFreq) {
557 : uInt c;
558 232 : MFrequency curFreq;
559 232 : for (uInt i = 0; i < which.nelements(); i++) {
560 116 : AlwaysAssert(which(i) >= 0, AipsError);
561 116 : c = which(i);
562 116 : SpectralModel& curSpectrum = component(c).spectrum();
563 116 : curFreq = curSpectrum.refFrequency();
564 116 : curFreq.set(newFreq);
565 116 : curSpectrum.setRefFrequency(curFreq);
566 : }
567 116 : DebugAssert(ok(), AipsError);
568 116 : }
569 :
570 37 : void ComponentList::setRefFrequencyFrame(const Vector<Int>& which,
571 : MFrequency::Types newFrame) {
572 : uInt c;
573 74 : MFrequency curFreq;
574 74 : const MFrequency::Ref newRef(newFrame);
575 74 : for (uInt i = 0; i < which.nelements(); i++) {
576 37 : AlwaysAssert(which(i) >= 0, AipsError);
577 37 : c = which(i);
578 37 : SpectralModel& curSpectrum = component(c).spectrum();
579 37 : curFreq = curSpectrum.refFrequency();
580 37 : curFreq.set(newRef);
581 37 : curSpectrum.setRefFrequency(curFreq);
582 : }
583 37 : DebugAssert(ok(), AipsError);
584 37 : }
585 :
586 :
587 37 : void ComponentList::setRefFrequencyUnit(const Vector<Int>& which,
588 : const Unit& unit) {
589 : uInt c;
590 74 : for (uInt i = 0; i < which.nelements(); i++) {
591 37 : AlwaysAssert(which(i) >= 0, AipsError);
592 37 : c = which(i);
593 37 : component(c).spectrum().convertFrequencyUnit(unit);
594 : }
595 37 : DebugAssert(ok(), AipsError);
596 37 : }
597 :
598 2272955 : SkyComponent& ComponentList::component(const uInt& index) {
599 : // AlwaysAssert(itsROFlag == false, AipsError);
600 2272955 : AlwaysAssert(index < nelements(), AipsError);
601 2272955 : DebugAssert(ok(), AipsError);
602 2272955 : return itsList[itsOrder[index]];
603 : }
604 :
605 14867 : const SkyComponent& ComponentList::component(const uInt& index) const {
606 14867 : DebugAssert(ok(), AipsError);
607 14867 : AlwaysAssert(index < nelements(), AipsError);
608 14867 : return itsList[itsOrder[index]];
609 : }
610 :
611 148 : void ComponentList::rename(const Path& fileName,
612 : const Table::TableOption option) {
613 148 : AlwaysAssert(option != Table::Old, AipsError);
614 148 : AlwaysAssert(itsROFlag == false, AipsError);
615 148 : DebugAssert(ok(), AipsError);
616 148 : if (fileName.length() != 0) {
617 : // See if this list is associated with a Table.
618 148 : if (itsTable.isNull()) {
619 : //createTable(fileName, option, addOptCol_p);
620 148 : createTable(fileName, option);
621 : } else {
622 0 : if (!itsTable.isWritable()) itsTable.reopenRW();
623 0 : itsTable.rename(fileName.absoluteName(), option);
624 : }
625 :
626 : // Ensure that the Table::isReadable(fileName) returns true, otherwise the
627 : // ok() function will fail.
628 148 : itsTable.flush();
629 148 : DebugAssert(ok(), AipsError);
630 : } else {
631 0 : if (!itsTable.isNull()) {
632 0 : itsTable.markForDelete();
633 0 : itsTable = Table();
634 0 : itsROFlag = false;
635 : }
636 : }
637 148 : }
638 :
639 38 : ComponentList ComponentList::copy() const {
640 38 : DebugAssert(ok(), AipsError);
641 38 : ComponentList copiedList;
642 76 : SkyComponent currentComp;
643 204 : for (uInt c = 0; c < nelements(); c++) {
644 166 : currentComp = component(c).copy();
645 166 : copiedList.add(currentComp);
646 : }
647 76 : return copiedList;
648 : }
649 :
650 0 : void ComponentList::sort(ComponentList::SortCriteria criteria) {
651 0 : Block<Double> val(nelements());
652 0 : Sort::Order order = Sort::Ascending;
653 0 : Bool doSort = true;
654 0 : switch (criteria) {
655 0 : case ComponentList::FLUX: {
656 0 : for (uInt i = 0; i < nelements(); i++) {
657 0 : itsList[i].flux().convertPol(ComponentType::STOKES);
658 0 : val[i] = abs(itsList[i].flux().value(0u).real());
659 : }
660 0 : order = Sort::Descending;
661 0 : break;
662 : }
663 0 : case ComponentList::POSITION: {
664 0 : MDirection compDir;
665 0 : MVDirection refDir(0.0, 0.0);
666 0 : Vector<Double> position(2);
667 0 : for (uInt i = 0; i < nelements(); i++) {
668 0 : val[i] = refDir.separation(itsList[i].shape().refDirection().getValue());
669 : }
670 0 : order = Sort::Ascending;
671 0 : break;
672 : }
673 0 : case ComponentList::POLARISATION: {
674 0 : Vector<Double> f(4);
675 0 : for (uInt i = 0; i < nelements(); i++) {
676 0 : itsList[i].flux().value(f);
677 0 : if (!nearAbs(f(0), 0.0, DBL_MIN)) {
678 0 : val[i] = sqrt(f(1)*f(1)+f(2)*f(2)+f(3)*f(3))/f(0);
679 : }
680 : else {
681 0 : val[i] = 0.0;
682 : }
683 : }
684 0 : order = Sort::Descending;
685 0 : break;
686 : }
687 0 : case ComponentList::UNSORTED:
688 : case ComponentList::NUMBER_CRITERIA:
689 0 : doSort = false;
690 0 : break;
691 : };
692 : // The genSort function requires a Vector<uInt> and not a Block<uInt> so
693 : // I'll create a temporary Vector here which references the data in the
694 : // 'itsOrder' Block.
695 0 : if (doSort) {
696 0 : Vector<uInt> vecOrder(IPosition(1,nelements()), itsOrder.storage(), SHARE);
697 0 : AlwaysAssert(genSort(vecOrder, val, order) == nelements(),
698 : AipsError);
699 : }
700 0 : }
701 :
702 0 : String ComponentList::name(ComponentList::SortCriteria enumerator) {
703 0 : switch (enumerator) {
704 0 : case ComponentList::FLUX: return "Flux";
705 0 : case ComponentList::POSITION: return "Position";
706 0 : case ComponentList::POLARISATION: return "Polarization";
707 0 : default: return "unknown";
708 : };
709 : }
710 :
711 0 : ComponentList::SortCriteria ComponentList::type(const String& criteria) {
712 0 : String canonicalCase(criteria);
713 0 : canonicalCase.capitalize();
714 0 : for (uInt i = 0; i < ComponentList::NUMBER_CRITERIA; i++) {
715 0 : if (canonicalCase.
716 0 : matches(ComponentList::name((ComponentList::SortCriteria) i))) {
717 0 : return (ComponentList::SortCriteria) i;
718 : }
719 : }
720 0 : return ComponentList::UNSORTED;
721 : }
722 :
723 2306974 : Bool ComponentList::ok() const {
724 : // The LogIO class is only constructed if an Error is detected for
725 : // performance reasons. Both function static and file static variables
726 : // where considered and rejected for this purpose.
727 2306974 : if (itsList.nelements() < itsNelements) {
728 0 : LogIO logErr(LogOrigin("ComponentList", "ok()"));
729 : logErr << LogIO::SEVERE
730 : << "The list size is inconsistant with its cached size"
731 0 : << LogIO::POST;
732 0 : return false;
733 : }
734 2306974 : if (itsROFlag == true && itsTable.isNull() == true) {
735 0 : LogIO logErr(LogOrigin("ComponentList", "ok()"));
736 : logErr << LogIO::SEVERE
737 : << "Only ComponentList's associated with a Table can be readonly"
738 0 : << LogIO::POST;
739 0 : return false;
740 : }
741 2306974 : if (itsTable.isNull() == false) {
742 8267 : String tablename = itsTable.tableName();
743 8267 : if (Table::isReadable(tablename) == false) {
744 0 : LogIO logErr(LogOrigin("ComponentList", "ok()"));
745 : logErr << LogIO::SEVERE
746 : << "Table associated with ComponentList is not readable"
747 0 : << LogIO::POST;
748 0 : return false;
749 : }
750 8267 : if (itsROFlag == false && Table::isWritable(tablename) == false) {
751 0 : LogIO logErr(LogOrigin("ComponentList", "ok()"));
752 : logErr << LogIO::SEVERE
753 : << "Table associated with ComponentList is not writeable"
754 0 : << LogIO::POST;
755 0 : return false;
756 : }
757 : }
758 4631824 : for (uInt i = 0; i < itsNelements; i++) {
759 2324850 : if (itsOrder[i] >= itsNelements) {
760 0 : LogIO logErr(LogOrigin("ComponentList", "ok()"));
761 : logErr << LogIO::SEVERE
762 : << "Cannot index to an element that is outside the list!"
763 0 : << LogIO::POST;
764 0 : return false;
765 : }
766 : }
767 2306974 : return true;
768 : }
769 :
770 148 : void ComponentList::createTable(const Path& fileName,
771 : const Table::TableOption option) {
772 : // const Table::TableOption option,
773 : // const Bool addOptCol) {
774 : // Build a default table description
775 444 : TableDesc td("ComponentListDescription", "4", TableDesc::Scratch);
776 148 : td.comment() = "A description of a component list";
777 : {
778 : {
779 : const ArrayColumnDesc<DComplex>
780 444 : fluxValCol(fluxName, "Flux values", IPosition(1,4),ColumnDesc::Direct);
781 148 : td.addColumn(fluxValCol);
782 : const ScalarColumnDesc<String>
783 296 : fluxUnitCol(fluxUnitName, "Flux units", ColumnDesc::Direct);
784 148 : td.addColumn(fluxUnitCol);
785 : const ScalarColumnDesc<String>
786 : fluxPolCol(fluxPolName, "Flux polarisation representation",
787 296 : ColumnDesc::Direct);
788 148 : td.addColumn(fluxPolCol);
789 : const ArrayColumnDesc<DComplex>
790 148 : fluxErrCol(fluxErrName, "Flux errors", IPosition(1,4),
791 296 : ColumnDesc::Direct);
792 148 : td.addColumn(fluxErrCol);
793 : }
794 : {
795 : const ScalarColumnDesc<String>
796 296 : shapeCol(shapeName, "Shape of the Component", ColumnDesc::Direct);
797 148 : td.addColumn(shapeCol);
798 296 : const String dirValColName = refDirName;
799 : const ArrayColumnDesc<Double>
800 : dirValCol(dirValColName, "Reference direction values",
801 444 : IPosition(1,2), ColumnDesc::Direct);
802 148 : td.addColumn(dirValCol);
803 296 : const String dirRefColName = "Direction_Frame";
804 : const ScalarColumnDesc<String>
805 : dirRefCol(dirRefColName, "The reference direction frame",
806 296 : ColumnDesc::Direct);
807 148 : td.addColumn(dirRefCol);
808 : {
809 296 : const TableMeasRefDesc dirRefTMCol(td, dirRefColName);
810 296 : const TableMeasValueDesc dirValTMCol(td, dirValColName);
811 296 : TableMeasDesc<MDirection> dirTMCol(dirValTMCol, dirRefTMCol);
812 148 : dirTMCol.write(td);
813 : }
814 : const ArrayColumnDesc<Double>
815 : dirErrCol(dirErrName, "Error in the reference direction values",
816 444 : IPosition(1,2), ColumnDesc::Direct);
817 148 : td.addColumn(dirErrCol);
818 : const ArrayColumnDesc<String>
819 : dirErrUnitCol(dirErrUnitName, "Units of the direction error",
820 444 : IPosition(1,2), ColumnDesc::Direct);
821 148 : td.addColumn(dirErrUnitCol);
822 : {
823 296 : TableQuantumDesc dirErrTMCol(td, dirErrName, dirErrUnitName);
824 148 : dirErrTMCol.write(td);
825 : }
826 :
827 : const ArrayColumnDesc<Double>
828 : shapeParmCol(shapeParName,
829 296 : "Parameters specific to the component shape", 1);
830 148 : td.addColumn(shapeParmCol);
831 : const ArrayColumnDesc<Double>
832 : shapeErrCol(shapeErrName,
833 148 : "Error in the shape parameters", 1);
834 148 : td.addColumn(shapeErrCol);
835 : }
836 : {
837 : const ScalarColumnDesc<String>
838 : freqShapeCol(spectrumName, "Shape of the spectrum",
839 296 : ColumnDesc::Direct);
840 148 : td.addColumn (freqShapeCol);
841 296 : const String freqValColName = refFreqName;
842 : const ScalarColumnDesc<Double>
843 : freqValCol(freqValColName, "The reference frequency values",
844 296 : ColumnDesc::Direct);
845 148 : td.addColumn(freqValCol);
846 296 : const String freqRefColName = "Frequency_Frame";
847 : const ScalarColumnDesc<String>
848 : freqRefCol(freqRefColName, "The reference frequency frame",
849 296 : ColumnDesc::Direct);
850 148 : td.addColumn(freqRefCol);
851 : {
852 296 : const TableMeasRefDesc freqRefTMCol(td, freqRefColName);
853 296 : const TableMeasValueDesc freqValTMCol(td, freqValColName);
854 296 : TableMeasDesc<MFrequency> freqTMCol(freqValTMCol, freqRefTMCol);
855 148 : freqTMCol.write(td);
856 : }
857 : const ScalarColumnDesc<Double>
858 : freqErrCol(freqErrName, "Error in the reference frequency",
859 296 : ColumnDesc::Direct);
860 148 : td.addColumn(freqErrCol);
861 : const ScalarColumnDesc<String>
862 : freqErrUnitCol(freqErrUnitName, "Units of the frequency error",
863 296 : ColumnDesc::Direct);
864 148 : td.addColumn(freqErrUnitCol);
865 : {
866 296 : TableQuantumDesc freqErrTMCol(td, freqErrName, freqErrUnitName);
867 148 : freqErrTMCol.write(td);
868 : }
869 : const ArrayColumnDesc<Double>
870 : specParmCol(spectParName,
871 296 : "Parameters specific to the components spectrum", 1);
872 148 : td.addColumn(specParmCol);
873 : const ArrayColumnDesc<Double>
874 : specErrCol(spectErrName,
875 148 : "Errors in the spectral parameters", 1);
876 148 : td.addColumn(specErrCol);
877 : }
878 : {
879 : const ScalarColumnDesc<String>
880 : labelCol(labelName, "An arbitrary label for the user",
881 148 : ColumnDesc::Direct);
882 148 : td.addColumn (labelCol);
883 : }
884 : }
885 148 : if (itsAddOptCol) {
886 : const ArrayColumnDesc<Double>
887 0 : optParCol(optParColName,"optional parameter column",1,ColumnDesc::Undefined);
888 0 : td.addColumn (optParCol);
889 0 : td.defineHypercolumn("TiledOptParms",1,stringToVector(optParColName));
890 : }
891 148 : SetupNewTable newTable(fileName.absoluteName(), td, option);
892 :
893 148 : if (itsAddOptCol) {
894 0 : TiledCellStMan optcolsm("TiledOptParms",IPosition(1,1));
895 0 : newTable.bindColumn(optParColName,optcolsm);
896 : }
897 148 : itsTable = Table(newTable, TableLock::AutoLocking, nelements(), false);
898 : {
899 148 : TableInfo& info(itsTable.tableInfo());
900 148 : info.setType(TableInfo::type(TableInfo::COMPONENTLIST));
901 148 : info.readmeAddLine(String(
902 : "This is a ComponentList Table containing parameterised representations"));
903 148 : info.readmeAddLine(String("of the sky brightness."));
904 : }
905 148 : }
906 :
907 195 : void ComponentList::writeTable() {
908 195 : if (itsTable.isWritable() == false) {
909 0 : itsTable.reopenRW();
910 : }
911 195 : DebugAssert(itsTable.isWritable(), AipsError);
912 : {
913 195 : const casacore::rownr_t nRows = itsTable.nrow();
914 195 : const casacore::rownr_t nelem = nelements();
915 195 : if (nRows < nelem) {
916 0 : itsTable.addRow(nelem - nRows);
917 195 : } else if (nRows > nelem) {
918 0 : Vector<casacore::rownr_t> rows(nRows - nelem);
919 0 : indgen(rows, nelem);
920 0 : itsTable.removeRow(RowNumbers(rows));
921 : }
922 : }
923 390 : ArrayColumn<DComplex> fluxValCol(itsTable, fluxName);
924 390 : ScalarColumn<String> fluxUnitCol(itsTable, fluxUnitName);
925 390 : ScalarColumn<String> fluxPolCol(itsTable, fluxPolName);
926 390 : ScalarColumn<String> shapeCol(itsTable, shapeName);
927 390 : MDirection::ScalarColumn dirCol(itsTable, refDirName);
928 390 : ArrayColumn<Double> shapeParmCol(itsTable, shapeParName);
929 390 : ScalarColumn<String> specShapeCol(itsTable, spectrumName);
930 390 : MFrequency::ScalarColumn freqCol(itsTable, refFreqName);
931 390 : ArrayColumn<Double> specShapeParmCol(itsTable, spectParName);
932 390 : ScalarColumn<String> labelCol(itsTable, labelName);
933 390 : ArrayColumn<DComplex> fluxErrCol;
934 390 : ArrayQuantColumn<Double> dirErrCol;
935 390 : ArrayColumn<Double> shapeErrCol;
936 390 : ScalarQuantColumn<Double> freqErrCol;
937 390 : ArrayColumn<Double> spectErrCol;
938 390 : ScalarColumn<TableRecord> specRecord;
939 390 : ArrayColumn<Double> optParCol;
940 : {
941 195 : const ColumnDescSet& cds=itsTable.tableDesc().columnDescSet();
942 195 : if (!cds.isDefined(spectralRecordName)) {
943 148 : itsTable.addColumn(ScalarRecordColumnDesc(spectralRecordName));
944 : }
945 195 : specRecord.attach(itsTable, spectralRecordName);
946 195 : if (!cds.isDefined(fluxErrName)) {
947 : itsTable.addColumn
948 0 : (ArrayColumnDesc<DComplex>(fluxErrName, "Flux errors", IPosition(1,4),
949 : ColumnDesc::Direct));
950 : }
951 :
952 195 : fluxErrCol.attach(itsTable, fluxErrName);
953 195 : if (!cds.isDefined(dirErrName)) {
954 : itsTable.addColumn
955 0 : (ArrayColumnDesc<Double>(dirErrName,
956 : "Error in the reference direction values",
957 0 : IPosition(1,2), ColumnDesc::Direct));
958 : itsTable.addColumn
959 0 : (ArrayColumnDesc<String>(dirErrUnitName,
960 : "Units of the direction error",
961 0 : IPosition(1,2), ColumnDesc::Direct));
962 : TableQuantumDesc dirErrTMCol(itsTable.tableDesc(), dirErrName,
963 0 : dirErrUnitName);
964 0 : dirErrTMCol.write(itsTable);
965 : }
966 195 : dirErrCol.attach(itsTable, dirErrName);
967 195 : if (!cds.isDefined(shapeErrName)) {
968 : itsTable.addColumn
969 0 : (ArrayColumnDesc<Double>(shapeErrName,
970 : "Error in the shape parameters", 1));
971 : }
972 195 : shapeErrCol.attach(itsTable, shapeErrName);
973 :
974 195 : if (!cds.isDefined(freqErrName)) {
975 : itsTable.addColumn
976 0 : (ScalarColumnDesc<Double>(freqErrName,
977 : "Error in the reference frequency",
978 : ColumnDesc::Direct));
979 : itsTable.addColumn
980 0 : (ScalarColumnDesc<String>(freqErrUnitName,
981 : "Units of the frequency error",
982 : ColumnDesc::Direct));
983 : TableQuantumDesc freqErrTMCol(itsTable.tableDesc(), freqErrName,
984 0 : freqErrUnitName);
985 0 : freqErrTMCol.write(itsTable);
986 : }
987 195 : freqErrCol.attach(itsTable, freqErrName);
988 195 : if (!cds.isDefined(spectErrName)) {
989 : itsTable.addColumn
990 0 : (ArrayColumnDesc<Double>(spectErrName,
991 : "Error in the spectral parameters", 1));
992 : }
993 195 : spectErrCol.attach(itsTable, spectErrName);
994 :
995 195 : if (itsAddOptCol) {
996 0 : if (!cds.isDefined(optParColName)) {
997 : itsTable.addColumn
998 0 : (ArrayColumnDesc<Double>(optParColName,
999 : "Optional parameters", 1));
1000 : //cerr<<"added optional parameter col"<<endl;
1001 : }
1002 0 : optParCol.attach(itsTable, optParColName);
1003 : }
1004 : }
1005 :
1006 390 : Vector<Quantum<Double> > dirErr(2);
1007 580 : for (uInt i = 0; i < nelements(); i++) {
1008 : {
1009 385 : const Flux<Double>& flux = component(i).flux();
1010 385 : fluxValCol.put(i, flux.value());
1011 385 : fluxUnitCol.put(i, flux.unit().getName());
1012 385 : fluxPolCol.put(i, ComponentType::name(flux.pol()));
1013 385 : fluxErrCol.put(i, flux.errors());
1014 : }
1015 : {
1016 385 : const ComponentShape& compShape = component(i).shape();
1017 385 : shapeCol.put(i, compShape.ident());
1018 385 : dirCol.put(i, compShape.refDirection());
1019 385 : if (!dirErrCol.isNull()) {
1020 385 : dirErr(0) = compShape.refDirectionErrorLat();
1021 385 : dirErr(1) = compShape.refDirectionErrorLong();
1022 385 : dirErrCol.put(i, dirErr);
1023 : }
1024 385 : shapeParmCol.put(i, compShape.parameters());
1025 385 : shapeErrCol.put(i, compShape.errors());
1026 : }
1027 : {
1028 385 : const SpectralModel& compSpectrum = component(i).spectrum();
1029 385 : specShapeCol.put(i, compSpectrum.ident());
1030 770 : TableRecord rec;
1031 385 : String err;
1032 385 : if(compSpectrum.toRecord(err, rec))
1033 385 : specRecord.put(i,rec);
1034 385 : freqCol.put(i, compSpectrum.refFrequency());
1035 385 : if (!freqErrCol.isNull()) {
1036 385 : freqErrCol.put(i, compSpectrum.refFrequencyError());
1037 : }
1038 385 : specShapeParmCol.put(i, compSpectrum.parameters());
1039 385 : spectErrCol.put(i, compSpectrum.errors());
1040 : }
1041 : {
1042 385 : labelCol.put(i, component(i).label());
1043 : }
1044 385 : if (itsAddOptCol) {
1045 0 : const ComponentShape& compShape2 = component(i).shape();
1046 0 : if (compShape2.type()==ComponentType::LDISK) {
1047 : //optParCol.put(i,compShape2.optParameters());
1048 0 : optParCol.put(i,component(i).optionalParameters());
1049 : }
1050 : }
1051 : }
1052 195 : }
1053 :
1054 146 : void ComponentList::readTable(const Path& fileName, const Bool readOnly) {
1055 : {
1056 146 : const String& fullName = fileName.absoluteName();
1057 146 : if (readOnly) {
1058 105 : AlwaysAssert(Table::isReadable(fullName), AipsError);
1059 105 : itsTable = Table(fullName, Table::Old);
1060 : }
1061 : else {
1062 41 : AlwaysAssert(Table::isWritable(fullName), AipsError);
1063 82 : itsTable = Table(fullName, TableLock(TableLock::AutoLocking),
1064 41 : Table::Update);
1065 : }
1066 : }
1067 292 : const ArrayColumn<DComplex> fluxValCol(itsTable, fluxName);
1068 292 : const ScalarColumn<String> fluxUnitCol(itsTable, fluxUnitName);
1069 292 : const ScalarColumn<String> fluxPolCol(itsTable, fluxPolName);
1070 292 : const ScalarColumn<String> shapeCol(itsTable, shapeName);
1071 292 : const MDirection::ScalarColumn dirCol(itsTable, refDirName);
1072 292 : const ArrayColumn<Double> shapeParmCol(itsTable, shapeParName);
1073 292 : const ScalarColumn<String> specShapeCol(itsTable, spectrumName);
1074 292 : const MFrequency::ScalarColumn freqCol(itsTable, refFreqName);
1075 292 : const ArrayColumn<Double> spectralParmCol(itsTable,spectParName);
1076 292 : const ScalarColumn<String> labelCol(itsTable, labelName);
1077 :
1078 292 : ArrayColumn<DComplex> fluxErrCol;
1079 292 : ROArrayQuantColumn<Double> dirErrCol;
1080 292 : ArrayColumn<Double> shapeErrCol;
1081 292 : ROScalarQuantColumn<Double> freqErrCol;
1082 292 : ArrayColumn<Double> spectralErrCol;
1083 292 : ScalarColumn<TableRecord> specRecord;
1084 292 : ArrayColumn<Double> optParmCol;
1085 : {// Old componentlist tables may not have the error columns
1086 146 : const ColumnDescSet& cds=itsTable.tableDesc().columnDescSet();
1087 146 : if (cds.isDefined(fluxErrName)) {
1088 146 : fluxErrCol.attach(itsTable, fluxErrName);
1089 : }
1090 146 : if (cds.isDefined(dirErrName)) {
1091 146 : dirErrCol.attach(itsTable, dirErrName);
1092 : }
1093 146 : if (cds.isDefined(shapeErrName)) {
1094 146 : shapeErrCol.attach(itsTable, shapeErrName);
1095 : }
1096 146 : if (cds.isDefined(freqErrName)) {
1097 146 : freqErrCol.attach(itsTable, freqErrName);
1098 : }
1099 146 : if (cds.isDefined(spectErrName)) {
1100 146 : spectralErrCol.attach(itsTable, spectErrName);
1101 : }
1102 146 : if (cds.isDefined(spectralRecordName)) {
1103 146 : specRecord.attach(itsTable, spectralRecordName);
1104 : }
1105 : // new optional parameter column
1106 146 : if (cds.isDefined(optParColName)) {
1107 0 : optParmCol.attach(itsTable,optParColName);
1108 : }
1109 : }
1110 :
1111 292 : SkyComponent currentComp;
1112 146 : const uInt nComp = fluxValCol.nrow();
1113 292 : Vector<DComplex> compFluxValue(4);
1114 292 : Vector<Double> shapeParms, spectralParms, optParms;
1115 292 : String compName, compLabel, compFluxPol, compFluxUnit, compSpectrum;
1116 292 : MDirection compDir;
1117 292 : MFrequency compFreq;
1118 292 : Vector<Quantum<Double> > newDirErr(2);
1119 146 : Quantum<Double> newFreqErr;
1120 454 : for (uInt i = 0; i < nComp; i++) {
1121 308 : shapeCol.get(i, compName);
1122 308 : specShapeCol.get(i, compSpectrum);
1123 : //cout << "CompSpec " << compSpectrum << " shape " <<ComponentType::spectralShape(compSpectrum) << endl;
1124 616 : currentComp = SkyComponent(ComponentType::shape(compName),
1125 616 : ComponentType::spectralShape(compSpectrum));
1126 : {
1127 308 : Flux<Double>& compFlux = currentComp.flux();
1128 308 : fluxValCol.get(i, compFluxValue);
1129 308 : compFlux.setValue(compFluxValue);
1130 308 : fluxUnitCol.get(i, compFluxUnit);
1131 308 : compFlux.setUnit(compFluxUnit);
1132 308 : fluxPolCol.get(i, compFluxPol);
1133 308 : compFlux.setPol(ComponentType::polarisation(compFluxPol));
1134 308 : if (!fluxErrCol.isNull()) {
1135 308 : fluxErrCol.get(i, compFluxValue);
1136 308 : compFlux.setErrors(compFluxValue(0), compFluxValue(1),
1137 308 : compFluxValue(2), compFluxValue(3));
1138 : }
1139 : }
1140 : {
1141 308 : ComponentShape& compShape = currentComp.shape();
1142 308 : dirCol.get(i, compDir);
1143 308 : compShape.setRefDirection(compDir);
1144 308 : if (!dirErrCol.isNull()) {
1145 308 : dirErrCol.get(i, newDirErr);
1146 308 : compShape.setRefDirectionError(newDirErr(0), newDirErr(1));
1147 : }
1148 308 : shapeParmCol.get(i, shapeParms, true);
1149 308 : compShape.setParameters(shapeParms);
1150 308 : if (!shapeErrCol.isNull()) {
1151 308 : shapeErrCol.get(i, shapeParms);
1152 308 : compShape.setErrors(shapeParms);
1153 : }
1154 : }
1155 : {
1156 308 : freqCol.get(i, compFreq);
1157 308 : SpectralModel& compSpectrum = currentComp.spectrum();
1158 :
1159 308 : if (!specRecord.isNull()) {
1160 308 : if(specRecord.isDefined(i)){
1161 616 : TableRecord rec;
1162 308 : specRecord.get(i, rec);
1163 616 : String err;
1164 308 : compSpectrum.fromRecord(err, rec);
1165 : }
1166 : }
1167 308 : compSpectrum.setRefFrequency(compFreq);
1168 308 : if (!freqErrCol.isNull()) {
1169 308 : freqErrCol.get(i, newFreqErr);
1170 308 : compSpectrum.setRefFrequencyError(newFreqErr);
1171 : }
1172 308 : spectralParmCol.get(i, spectralParms, true);
1173 308 : compSpectrum.setParameters(spectralParms);
1174 308 : if (!spectralErrCol.isNull()) {
1175 308 : spectralErrCol.get(i, spectralParms);
1176 308 : compSpectrum.setErrors(spectralParms);
1177 : }
1178 : }
1179 : {
1180 308 : labelCol.get(i, currentComp.label());
1181 : }
1182 : {
1183 308 : if (!optParmCol.isNull()) {
1184 0 : if (optParmCol.isDefined(i)) {
1185 0 : ComponentShape& compShape2 = currentComp.shape();
1186 0 : optParmCol.get(i,optParms,true);
1187 0 : compShape2.setOptParameters(optParms);
1188 0 : currentComp.optionalParameters()=optParms;
1189 : }
1190 : }
1191 : }
1192 :
1193 308 : add(currentComp);
1194 : }
1195 146 : itsROFlag = readOnly;
1196 146 : }
1197 :
1198 948 : Bool ComponentList::toRecord(String& error, RecordInterface& outRec) const {
1199 :
1200 948 : Bool retval=true;
1201 :
1202 948 : outRec.define("nelements", itsNelements);
1203 2100 : for (uInt k=0; k < itsNelements; ++k){
1204 2304 : Record componentContainer;
1205 1152 : retval= (retval && component(k).toRecord(error, componentContainer));
1206 2304 : String componentId=String("component")+String::toString(k);
1207 1152 : outRec.defineRecord(componentId, componentContainer);
1208 :
1209 : }
1210 948 : return retval;
1211 :
1212 : }
1213 :
1214 224 : Bool ComponentList::fromRecord(String& error, const RecordInterface& inRec){
1215 :
1216 224 : Bool retval= true;
1217 224 : if(itsNelements > 0){
1218 0 : LogIO logErr(LogOrigin("ComponentList", "fromRecord()"));
1219 : logErr << LogIO::SEVERE
1220 : << "Trying to overwrite a non-empty componentList from Record"
1221 0 : << LogIO::POST;
1222 0 : return false;
1223 :
1224 : }
1225 :
1226 224 : uInt nelements=0;
1227 224 : if (inRec.isDefined("nelements")) {
1228 224 : inRec.get("nelements", nelements);
1229 224 : if(nelements >0){
1230 516 : for(uInt k=0; k < nelements; ++k){
1231 584 : String componentId=String("component")+String::toString(k);
1232 292 : Record componentRecord=inRec.asRecord(componentId);
1233 292 : SkyComponent tempComponent;
1234 292 : retval=(retval && tempComponent.fromRecord(error, componentRecord));
1235 292 : if(retval){
1236 292 : add(tempComponent);
1237 : }
1238 : else{
1239 0 : return retval;
1240 : }
1241 : }
1242 : }
1243 : }
1244 224 : return retval;
1245 : }
1246 :
1247 0 : String ComponentList::summarize(uInt index) const {
1248 0 : AlwaysAssert(index < nelements(), AipsError);
1249 0 : DebugAssert(ok(), AipsError);
1250 0 : return component(index).summarize();
1251 : }
1252 :
1253 38 : Table& ComponentList::_getTable() {
1254 38 : return itsTable;
1255 : }
1256 :
1257 19 : const Table& ComponentList::getTable() const {
1258 19 : return itsTable;
1259 : }
1260 :
1261 : // Local Variables:
1262 : // compile-command: "gmake ComponentList"
1263 : // End:
1264 :
1265 : } //# NAMESPACE CASA - END
1266 :
|
