Line data Source code
1 : //# SynDataSampling.cc: Implementation of SynDataSampling class
2 : //# Copyright (C) 1997,1998,1999,2000,2001,2003
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$
27 :
28 : #include <synthesis/DataSampling/SynDataSampling.h>
29 : #include <casacore/coordinates/Coordinates/CoordinateSystem.h>
30 : #include <casacore/coordinates/Coordinates/DirectionCoordinate.h>
31 : #include <casacore/ms/MeasurementSets/MSColumns.h>
32 : #include <casacore/casa/BasicSL/Constants.h>
33 : #include <synthesis/TransformMachines/SkyJones.h>
34 : #include <msvis/MSVis/VisSet.h>
35 : #include <msvis/MSVis/VisBuffer.h>
36 : #include <msvis/MSVis/VisibilityIterator.h>
37 : #include <casacore/images/Images/ImageInterface.h>
38 : #include <casacore/images/Images/PagedImage.h>
39 : #include <casacore/images/Images/TempImage.h>
40 : #include <casacore/casa/Arrays/ArrayLogical.h>
41 : #include <casacore/casa/Arrays/ArrayMath.h>
42 : #include <casacore/casa/Arrays/MaskedArray.h>
43 : #include <casacore/casa/Arrays/Array.h>
44 : #include <casacore/casa/Arrays/Array.h>
45 : #include <casacore/casa/Arrays/Vector.h>
46 : #include <casacore/casa/Arrays/Matrix.h>
47 : #include <casacore/casa/BasicSL/String.h>
48 : #include <casacore/casa/Utilities/Assert.h>
49 : #include <casacore/casa/Exceptions/Error.h>
50 : #include <casacore/casa/System/ProgressMeter.h>
51 : #include <sstream>
52 :
53 : using namespace casacore;
54 : namespace casa { //# NAMESPACE CASA - BEGIN
55 :
56 0 : SynDataSampling::SynDataSampling(MeasurementSet& ms,
57 : const CoordinateSystem& coords,
58 : const IPosition& shape,
59 0 : const Quantity& sigma)
60 : {
61 :
62 0 : LogIO os(LogOrigin("SynDataSampling", "SynDataSampling"));
63 :
64 0 : DataSampling::IDLScript_p="@app_syn";
65 :
66 : // Now create the VisSet
67 0 : Block<int> sort(4);
68 0 : sort[0] = MS::FIELD_ID;
69 0 : sort[1] = MS::ARRAY_ID;
70 0 : sort[2] = MS::DATA_DESC_ID;
71 0 : sort[3] = MS::TIME;
72 :
73 0 : Matrix<Int> noselection;
74 0 : VisSet vs(ms,sort,noselection);
75 :
76 : // First get the CoordinateSystem for the image and then find
77 : // the DirectionCoordinate
78 0 : Int directionIndex=coords.findCoordinate(Coordinate::DIRECTION);
79 0 : AlwaysAssert(directionIndex>=0, AipsError);
80 0 : DirectionCoordinate directionCoord=coords.directionCoordinate(directionIndex);
81 :
82 0 : dx_p.resize(IPosition(2, 2, 2*ms.nrow())); dx_p=-1.0;
83 0 : data_p.resize(IPosition(1, 2*ms.nrow())); data_p=0.0;
84 0 : sigma_p.resize(IPosition(1, 2*ms.nrow())); sigma_p=-1.0;
85 0 : Int lastRow = 0;
86 :
87 0 : VisIter& vi=vs.iter();
88 0 : VisBuffer vb(vi);
89 0 : vi.originChunks();
90 0 : vi.origin();
91 :
92 0 : Int nx = shape(0);
93 0 : Int ny = shape(1);
94 :
95 0 : Vector<Float> uvScale, uvOffset;
96 0 : uvScale(0)=(Float(nx)*coords.increment()(0));
97 0 : uvScale(1)=(Float(ny)*coords.increment()(1));
98 0 : uvOffset(0)=nx/2;
99 0 : uvOffset(1)=ny/2;
100 :
101 : // Now fill in the data and position columns
102 0 : ProgressMeter pm(1.0, Double(ms.nrow()), "Sampling Data", "", "", "", true);
103 :
104 : // Loop over all visibilities
105 0 : Int cohDone = 0;
106 0 : Float sigmaVal=sigma.getValue();
107 0 : for (vi.originChunks();vi.moreChunks();vi.nextChunk()) {
108 0 : for (vi.origin(); vi.more(); vi++) {
109 :
110 0 : Matrix<Double> uvw(3, vb.uvw().nelements());
111 0 : uvw=0.0;
112 0 : Vector<Double> dphase(vb.uvw().nelements());
113 0 : dphase=0.0;
114 0 : for (Int i=0;i<=vb.nRow();i++) {
115 0 : for (Int idim=0;idim<3;idim++) uvw(idim,i)=vb.uvw()(i)(idim);
116 : }
117 : // rotateUVW(uvw, dphase, vb);
118 :
119 0 : Vector<Double> worldPos(2);
120 0 : Vector<Double> xyPos(2);
121 :
122 0 : for (Int row=0;row<vb.nRow();row++) {
123 :
124 0 : IPosition irow(1, lastRow);
125 0 : IPosition irow1(1, lastRow+1);
126 0 : data_p(irow) = real(vb.correctedVisCube()(0,0,row));
127 0 : data_p(irow1) = imag(vb.correctedVisCube()(0,0,row));
128 :
129 0 : worldPos(0)=uvw(IPosition(2, 0, row));
130 0 : worldPos(1)=uvw(IPosition(2, 1, row));
131 :
132 0 : for (Int i=0;i<2;i++) {
133 0 : xyPos(i) = worldPos(i)*uvScale(i) + uvOffset(i);
134 : }
135 :
136 0 : dx_p(IPosition(2, 0, lastRow)) = Float(xyPos(0));
137 0 : dx_p(IPosition(2, 1, lastRow)) = Float(xyPos(1));
138 0 : dx_p(IPosition(2, 0, lastRow+1)) = Float(xyPos(0));
139 0 : dx_p(IPosition(2, 1, lastRow+1)) = Float(xyPos(1));
140 0 : if(sigmaVal>0.0) {
141 0 : sigma_p(irow) = sigmaVal;
142 0 : sigma_p(irow1) = sigmaVal;
143 : }
144 : else {
145 0 : sigma_p(irow) = vb.sigma()(row);
146 0 : sigma_p(irow1) = vb.sigma()(row);
147 : }
148 0 : lastRow+=2;
149 : }
150 0 : cohDone+=vb.nRow();
151 0 : pm.update(Double(cohDone));
152 : }
153 : }
154 0 : }
155 :
156 : //----------------------------------------------------------------------
157 0 : SynDataSampling& SynDataSampling::operator=(const SynDataSampling& other)
158 : {
159 : if(this!=&other) {
160 : };
161 0 : return *this;
162 : };
163 :
164 : //----------------------------------------------------------------------
165 0 : SynDataSampling::SynDataSampling(const SynDataSampling& other):
166 0 : DataSampling(other)
167 : {
168 0 : operator=(other);
169 0 : }
170 :
171 : //----------------------------------------------------------------------
172 0 : SynDataSampling::~SynDataSampling() {
173 0 : }
174 :
175 0 : void SynDataSampling::ok() {
176 0 : }
177 :
178 : } //# NAMESPACE CASA - END
179 :
|