Line data Source code
1 : //# SIMapper.cc: Implementation of SIMapper.h
2 : //# Copyright (C) 1997-2008
3 : //# Associated Universities, Inc. Washington DC, USA.
4 : //#
5 : //# This program is free software; you can redistribute it and/or modify it
6 : //# under the terms of the GNU General Public License as published by the Free
7 : //# Software Foundation; either version 2 of the License, or (at your option)
8 : //# any later version.
9 : //#
10 : //# This program 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 General Public License for
13 : //# more details.
14 : //#
15 : //# You should have received a copy of the GNU General Public License along
16 : //# with this program; if not, write to the Free Software Foundation, Inc.,
17 : //# 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 <casacore/casa/Exceptions/Error.h>
29 : #include <iostream>
30 : #include <sstream>
31 :
32 : #include <casacore/casa/Arrays/Matrix.h>
33 : #include <casacore/casa/Arrays/ArrayMath.h>
34 : #include <casacore/casa/Arrays/ArrayLogical.h>
35 :
36 : #include <casacore/casa/Logging.h>
37 : #include <casacore/casa/Logging/LogIO.h>
38 : #include <casacore/casa/Logging/LogMessage.h>
39 : #include <casacore/casa/Logging/LogSink.h>
40 : #include <casacore/casa/Logging/LogMessage.h>
41 :
42 : #include <casacore/casa/OS/DirectoryIterator.h>
43 : #include <casacore/casa/OS/File.h>
44 : #include <casacore/casa/OS/Path.h>
45 :
46 : #include <casacore/casa/OS/HostInfo.h>
47 : #include <casacore/lattices/Lattices/LatticeLocker.h>
48 : #include <casacore/ms/MeasurementSets/MSHistoryHandler.h>
49 : #include <casacore/ms/MeasurementSets/MeasurementSet.h>
50 :
51 : #include <synthesis/TransformMachines/BeamSkyJones.h>
52 : #include <synthesis/TransformMachines/SkyJones.h>
53 : #include <synthesis/TransformMachines/SimpleComponentFTMachine.h>
54 : #include <synthesis/TransformMachines2/SimpleComponentFTMachine.h>
55 : #include <synthesis/TransformMachines/SimpCompGridMachine.h>
56 : //#include <synthesis/ImagerObjects/SIMapperBase.h>
57 : #include <synthesis/ImagerObjects/SIMapper.h>
58 :
59 :
60 : #include <sys/types.h>
61 : #include <unistd.h>
62 : using namespace std;
63 :
64 : using namespace casacore;
65 : namespace casa { //# NAMESPACE CASA - BEGIN
66 :
67 0 : SIMapper::SIMapper( CountedPtr<SIImageStore>& imagestore,
68 : CountedPtr<FTMachine>& ftm,
69 0 : CountedPtr<FTMachine>& iftm) : useViVb2_p(false)
70 : {
71 : // LogIO os( LogOrigin("SIMapper","Construct a mapper",WHERE) );
72 0 : ft_p=ftm;
73 0 : ift_p=iftm;
74 0 : ft2_p=NULL;
75 0 : ift2_p=NULL;
76 0 : cft_p=NULL;
77 0 : itsImages=imagestore;
78 0 : }
79 1707 : SIMapper::SIMapper( CountedPtr<SIImageStore>& imagestore,
80 : CountedPtr<refim::FTMachine>& ftm,
81 1707 : CountedPtr<refim::FTMachine>& iftm) : useViVb2_p(true)
82 : {
83 : //LogIO os( LogOrigin("SIMapper","Construct a mapper",WHERE) );
84 1707 : ft2_p=ftm;
85 1707 : ift2_p=iftm;
86 1707 : cft_p=NULL;
87 1707 : itsImages=imagestore;
88 1707 : }
89 :
90 0 : SIMapper::SIMapper(const ComponentList& cl,
91 0 : String& whichMachine)
92 : {
93 0 : ft_p=NULL;
94 0 : ift_p=NULL;
95 0 : ft2_p=NULL;
96 0 : ift2_p=NULL;
97 :
98 0 : itsImages=NULL;
99 :
100 0 : if(whichMachine=="SimpleComponentFTMachine"){
101 0 : cft_p=new SimpleComponentFTMachine();
102 0 : cft2_p=new refim::SimpleComponentFTMachine();
103 : }
104 : else{
105 : //SD style component gridding
106 0 : cft_p=new SimpleComponentGridMachine();
107 : /////NEED to be done here
108 : //cft2_p=
109 : }
110 0 : cl_p=cl;
111 :
112 0 : }
113 :
114 3414 : SIMapper::~SIMapper()
115 : {
116 3414 : }
117 :
118 : // #############################################
119 : // #############################################
120 : // ####### Gridding / De-gridding functions ###########
121 : // #############################################
122 : // #############################################
123 2325 : void SIMapper::initializeGrid(vi::VisBuffer2& vb, Bool dopsf, Bool /*firstaccess*/)
124 : {
125 :
126 : //LogIO os( LogOrigin("SIMapper","initializeGrid",WHERE) );
127 2325 : if(!useViVb2_p)
128 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
129 : //Componentlist FTM has nothing to do
130 :
131 2325 : if(ift2_p.null())
132 0 : return;
133 :
134 2325 : ift2_p->initializeToSkyNew( dopsf, vb, itsImages);
135 :
136 : /////////////DEBUG
137 : // CoordinateSystem csys = itsImages->getCSys();
138 : // cout << "SIMapper : im spectral axis : " << csys.spectralCoordinate().referenceValue() << " at " << csys.spectralCoordinate().referencePixel() << " with increment " << csys.spectralCoordinate().increment() << endl;
139 :
140 : }
141 :
142 0 : void SIMapper::initializeGrid(VisBuffer& vb, Bool dopsf, Bool /*firstaccess*/)
143 : {
144 :
145 : //LogIO os( LogOrigin("SIMapper","initializeGrid",WHERE) );
146 0 : if(useViVb2_p)
147 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
148 : //Componentlist FTM has nothing to do
149 0 : if(ift_p.null())
150 0 : return;
151 :
152 0 : ift_p->initializeToSkyNew( dopsf, vb, itsImages);
153 :
154 : /////////////DEBUG
155 : // CoordinateSystem csys = itsImages->getCSys();
156 : // cout << "SIMapper : im spectral axis : " << csys.spectralCoordinate().referenceValue() << " at " << csys.spectralCoordinate().referencePixel() << " with increment " << csys.spectralCoordinate().increment() << endl;
157 :
158 : }
159 :
160 :
161 681176 : void SIMapper::grid(vi::VisBuffer2& vb, Bool dopsf, refim::FTMachine::Type col,
162 : const Int whichFTM)
163 : {
164 : //LogIO os( LogOrigin("SIMapper","grid",WHERE) );
165 681176 : if(!useViVb2_p)
166 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
167 : //Componentlist FTM has no gridding to do
168 : (void)whichFTM;
169 :
170 681176 : if(ift2_p.null())
171 0 : return;
172 :
173 681176 : ift2_p->put(vb, -1, dopsf, col);
174 :
175 : }
176 :
177 : /////////////////OLD vi/vb version
178 0 : void SIMapper::grid(VisBuffer& vb, Bool dopsf, FTMachine::Type col,
179 : const Int whichFTM)
180 : {
181 : //LogIO os( LogOrigin("SIMapper","grid",WHERE) );
182 0 : if(useViVb2_p)
183 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
184 : //Componentlist FTM has no gridding to do
185 : (void)whichFTM;
186 :
187 0 : if(ift_p.null())
188 0 : return;
189 :
190 0 : ift_p->put(vb, -1, dopsf, col);
191 :
192 : }
193 :
194 2317 : void SIMapper::finalizeGrid(vi::VisBuffer2& vb, Bool dopsf)
195 : {
196 : //LogIO os( LogOrigin("SIMapper","finalizeGrid",WHERE) );
197 2317 : if(!useViVb2_p)
198 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
199 2317 : if(ift2_p.null())
200 0 : return;
201 :
202 2317 : ift2_p->finalizeToSkyNew( dopsf, vb, itsImages );
203 :
204 : }
205 : //////////////OLD VI/VB version
206 0 : void SIMapper::finalizeGrid(VisBuffer& vb, Bool dopsf)
207 : {
208 : //LogIO os( LogOrigin("SIMapper","finalizeGrid",WHERE) );
209 0 : if(useViVb2_p)
210 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
211 0 : if(ift_p.null())
212 0 : return;
213 :
214 0 : ift_p->finalizeToSkyNew( dopsf, vb, itsImages );
215 :
216 : }
217 :
218 813 : void SIMapper::initializeDegrid(vi::VisBuffer2& vb, const Int /*row*/)
219 : {
220 : //LogIO os( LogOrigin("SIMapper", "initializeDegrid",WHERE) );
221 813 : if(!useViVb2_p)
222 0 : throw(AipsError("Programmer Error: using vi2 mode with vii constructor"));
223 813 : if(ft2_p.null() && cft2_p.null())
224 0 : return;
225 :
226 813 : if(ft2_p)
227 813 : ft2_p->initializeToVisNew(vb, itsImages);
228 :
229 : }
230 : //////////////////OLD vi/vb version
231 0 : void SIMapper::initializeDegrid(VisBuffer& vb, const Int /*row*/)
232 : {
233 : //LogIO os( LogOrigin("SIMapper", "initializeDegrid",WHERE) );
234 0 : if(useViVb2_p)
235 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
236 :
237 0 : if(ft_p.null() && cft_p.null())
238 0 : return;
239 0 : if(ft_p)
240 0 : ft_p->initializeToVisNew(vb, itsImages);
241 :
242 : }
243 :
244 162534 : void SIMapper::degrid(vi::VisBuffer2& vb)
245 : {
246 : //LogIO os( LogOrigin("SIMapper","degrid",WHERE) );
247 162534 : if(!useViVb2_p)
248 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
249 :
250 :
251 : ///This should not be called even but heck let's ignore
252 162534 : if(ft2_p.null() and cft2_p.null())
253 0 : return;
254 :
255 325068 : Cube<Complex> origCube;
256 162534 : origCube.assign(vb.visCubeModel());
257 :
258 :
259 :
260 325068 : Cube<Complex> mod(vb.nCorrelations(), vb.nChannels(), vb.nRows(), Complex(0.0));
261 162534 : vb.setVisCubeModel( mod);
262 162534 : if( ! ft2_p.null() ) { ft2_p->get(vb); }
263 162534 : if( ! cft2_p.null() ) { cft2_p->get(vb, cl_p); }
264 :
265 162534 : origCube+=vb.visCubeModel();
266 162534 : vb.setVisCubeModel(origCube);
267 :
268 : }
269 :
270 648 : void SIMapper::addPB(vi::VisBuffer2& vb, PBMath& pbMath, const MDirection& altDir, const Bool useAltDir)
271 : {
272 1296 : CoordinateSystem imageCoord=itsImages->pb()->coordinates();
273 :
274 1296 : IPosition imShape=itsImages->pb()->shape();
275 :
276 :
277 :
278 1296 : MDirection wcenter;
279 648 : if(useAltDir)
280 2 : wcenter=altDir;
281 : else{
282 646 : MSColumns mscol(vb.ms());
283 646 : wcenter=mscol.field().phaseDirMeas(vb.fieldId()(0), vb.time()(0));
284 : }
285 1296 : TempImage<Float> pbTemp(imShape, imageCoord);
286 1296 : TempImage<Complex> ctemp(imShape, imageCoord);
287 648 : ctemp.set(1.0);
288 648 : pbMath.applyPB(ctemp, ctemp, wcenter, Quantity(0.0, "deg"), BeamSquint::NONE);
289 648 : StokesImageUtil::To(pbTemp, ctemp);
290 648 : LatticeLocker lock1(*(itsImages->pb()), FileLocker::Write);
291 648 : itsImages->pb()->copyData( (LatticeExpr<Float>)((*(itsImages->pb()))+pbTemp) );
292 :
293 648 : }//addPB
294 :
295 :
296 :
297 : ////////////////////Old vi/Vb version
298 :
299 0 : void SIMapper::degrid(VisBuffer& vb)
300 : {
301 : //LogIO os( LogOrigin("SIMapper","degrid",WHERE) );
302 :
303 0 : if(useViVb2_p)
304 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
305 : ///This should not be called even but heck let's ignore
306 0 : if(ft_p.null() and cft_p.null())
307 0 : return;
308 :
309 0 : Cube<Complex> origCube;
310 0 : origCube.assign(vb.modelVisCube());
311 :
312 0 : vb.setModelVisCube( Complex(0.0,0.0) );
313 :
314 0 : if( ! ft_p.null() ) { ft_p->get(vb); }
315 0 : if( ! cft_p.null() ) { cft_p->get(vb, cl_p); }
316 :
317 0 : vb.modelVisCube()+=origCube;
318 :
319 : }
320 :
321 :
322 813 : void SIMapper::finalizeDegrid()
323 : {
324 : //LogIO os( LogOrigin("SIMapper","finalizeDegrid",WHERE) );
325 813 : if(!useViVb2_p)
326 0 : ft_p->finalizeToVis();
327 : else
328 813 : ft2_p->finalizeToVis();
329 813 : }
330 :
331 :
332 17 : Bool SIMapper::getFTMRecord(Record& rec, const String diskimage)
333 : {
334 : //LogIO os( LogOrigin("SIMapper","getFTMRecord",WHERE) );
335 17 : if(ft_p.null() && !useViVb2_p)
336 0 : return false;
337 17 : if(ft2_p.null() && useViVb2_p)
338 0 : return false;
339 34 : String err;
340 17 : if(!useViVb2_p)
341 0 : return ft_p->toRecord(err, rec, true, diskimage);
342 17 : return ft2_p->toRecord(err, rec, true, diskimage);
343 : // rec = itsFTM->toRecord();
344 :
345 : }
346 17 : Bool SIMapper::getCLRecord(Record& rec)
347 : {
348 17 : if(cft_p.null() && cft2_p.null())
349 17 : return false;
350 0 : String err;
351 0 : return cl_p.toRecord(err, rec);
352 : }
353 :
354 :
355 616 : void SIMapper::initPB()
356 : {
357 616 : LatticeLocker lock1(*(itsImages->pb()), FileLocker::Write);
358 616 : itsImages->pb()->set(0.0);
359 616 : }
360 :
361 0 : void SIMapper::addPB(VisBuffer& vb, PBMath& pbMath)
362 : {
363 0 : CoordinateSystem imageCoord=itsImages->pb()->coordinates();
364 :
365 0 : IPosition imShape=itsImages->pb()->shape();
366 :
367 0 : MDirection wcenter=vb.msColumns().field().phaseDirMeas(vb.fieldId());
368 0 : TempImage<Float> pbTemp(imShape, imageCoord);
369 0 : TempImage<Complex> ctemp(imShape, imageCoord);
370 0 : ctemp.set(1.0);
371 0 : pbMath.applyPB(ctemp, ctemp, wcenter, Quantity(0.0, "deg"), BeamSquint::NONE);
372 0 : StokesImageUtil::To(pbTemp, ctemp);
373 0 : LatticeLocker lock1(*(itsImages->pb()), FileLocker::Write);
374 0 : itsImages->pb()->copyData( (LatticeExpr<Float>)((*(itsImages->pb()))+pbTemp) );
375 :
376 0 : }//addPB
377 :
378 :
379 : } //# NAMESPACE CASA - END
380 :
381 :
|