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: casa-feedback@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 1760 : SIMapper::SIMapper( CountedPtr<SIImageStore>& imagestore,
80 : CountedPtr<refim::FTMachine>& ftm,
81 1760 : CountedPtr<refim::FTMachine>& iftm) : useViVb2_p(true)
82 : {
83 : //LogIO os( LogOrigin("SIMapper","Construct a mapper",WHERE) );
84 1760 : ft2_p=ftm;
85 1760 : ift2_p=iftm;
86 1760 : cft_p=NULL;
87 1760 : itsImages=imagestore;
88 1760 : }
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 3520 : SIMapper::~SIMapper()
115 : {
116 3520 : }
117 :
118 : // #############################################
119 : // #############################################
120 : // ####### Gridding / De-gridding functions ###########
121 : // #############################################
122 : // #############################################
123 2386 : void SIMapper::initializeGrid(vi::VisBuffer2& vb, Bool dopsf, Bool /*firstaccess*/)
124 : {
125 :
126 : //LogIO os( LogOrigin("SIMapper","initializeGrid",WHERE) );
127 2386 : if(!useViVb2_p)
128 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
129 : //Componentlist FTM has nothing to do
130 :
131 2386 : if(ift2_p.null())
132 0 : return;
133 :
134 2386 : 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 160 : void SIMapper::handleNewMs(const casacore::MeasurementSet &ms)
161 : {
162 160 : if (not ift2_p.null()) {
163 160 : ift2_p->handleNewMs(ms, imageStore());
164 : }
165 160 : }
166 :
167 690122 : void SIMapper::grid(vi::VisBuffer2& vb, Bool dopsf, refim::FTMachine::Type col,
168 : const Int whichFTM)
169 : {
170 : //LogIO os( LogOrigin("SIMapper","grid",WHERE) );
171 690122 : if(!useViVb2_p)
172 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
173 : //Componentlist FTM has no gridding to do
174 : (void)whichFTM;
175 :
176 690122 : if(ift2_p.null())
177 0 : return;
178 :
179 690122 : ift2_p->put(vb, -1, dopsf, col);
180 :
181 : }
182 :
183 : /////////////////OLD vi/vb version
184 0 : void SIMapper::grid(VisBuffer& vb, Bool dopsf, FTMachine::Type col,
185 : const Int whichFTM)
186 : {
187 : //LogIO os( LogOrigin("SIMapper","grid",WHERE) );
188 0 : if(useViVb2_p)
189 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
190 : //Componentlist FTM has no gridding to do
191 : (void)whichFTM;
192 :
193 0 : if(ift_p.null())
194 0 : return;
195 :
196 0 : ift_p->put(vb, -1, dopsf, col);
197 :
198 : }
199 :
200 2377 : void SIMapper::finalizeGrid(vi::VisBuffer2& vb, Bool dopsf)
201 : {
202 : //LogIO os( LogOrigin("SIMapper","finalizeGrid",WHERE) );
203 2377 : if(!useViVb2_p)
204 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
205 2377 : if(ift2_p.null())
206 0 : return;
207 :
208 2377 : ift2_p->finalizeToSkyNew( dopsf, vb, itsImages );
209 :
210 : }
211 : //////////////OLD VI/VB version
212 0 : void SIMapper::finalizeGrid(VisBuffer& vb, Bool dopsf)
213 : {
214 : //LogIO os( LogOrigin("SIMapper","finalizeGrid",WHERE) );
215 0 : if(useViVb2_p)
216 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
217 0 : if(ift_p.null())
218 0 : return;
219 :
220 0 : ift_p->finalizeToSkyNew( dopsf, vb, itsImages );
221 :
222 : }
223 :
224 829 : void SIMapper::initializeDegrid(vi::VisBuffer2& vb, const Int /*row*/)
225 : {
226 : //LogIO os( LogOrigin("SIMapper", "initializeDegrid",WHERE) );
227 829 : if(!useViVb2_p)
228 0 : throw(AipsError("Programmer Error: using vi2 mode with vii constructor"));
229 829 : if(ft2_p.null() && cft2_p.null())
230 0 : return;
231 :
232 829 : if(ft2_p)
233 829 : ft2_p->initializeToVisNew(vb, itsImages);
234 :
235 : }
236 : //////////////////OLD vi/vb version
237 0 : void SIMapper::initializeDegrid(VisBuffer& vb, const Int /*row*/)
238 : {
239 : //LogIO os( LogOrigin("SIMapper", "initializeDegrid",WHERE) );
240 0 : if(useViVb2_p)
241 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
242 :
243 0 : if(ft_p.null() && cft_p.null())
244 0 : return;
245 0 : if(ft_p)
246 0 : ft_p->initializeToVisNew(vb, itsImages);
247 :
248 : }
249 :
250 164368 : void SIMapper::degrid(vi::VisBuffer2& vb)
251 : {
252 : //LogIO os( LogOrigin("SIMapper","degrid",WHERE) );
253 164368 : if(!useViVb2_p)
254 0 : throw(AipsError("Programmer Error: using vi2 mode with vi1 constructor"));
255 :
256 :
257 : ///This should not be called even but heck let's ignore
258 164368 : if(ft2_p.null() and cft2_p.null())
259 0 : return;
260 :
261 164368 : Cube<Complex> origCube;
262 164368 : origCube.assign(vb.visCubeModel());
263 :
264 :
265 :
266 164368 : Cube<Complex> mod(vb.nCorrelations(), vb.nChannels(), vb.nRows(), Complex(0.0));
267 164368 : vb.setVisCubeModel( mod);
268 164368 : if( ! ft2_p.null() ) { ft2_p->get(vb); }
269 164368 : if( ! cft2_p.null() ) { cft2_p->get(vb, cl_p); }
270 :
271 164368 : origCube+=vb.visCubeModel();
272 164368 : vb.setVisCubeModel(origCube);
273 :
274 164368 : }
275 :
276 665 : void SIMapper::addPB(vi::VisBuffer2& vb, PBMath& pbMath, const MDirection& altDir, const Bool useAltDir)
277 : {
278 665 : CoordinateSystem imageCoord=itsImages->pb()->coordinates();
279 :
280 665 : IPosition imShape=itsImages->pb()->shape();
281 :
282 :
283 :
284 665 : MDirection wcenter;
285 665 : if(useAltDir)
286 5 : wcenter=altDir;
287 : else{
288 660 : MSColumns mscol(vb.ms());
289 660 : wcenter=mscol.field().phaseDirMeas(vb.fieldId()(0), vb.time()(0));
290 660 : }
291 665 : TempImage<Float> pbTemp(imShape, imageCoord);
292 665 : TempImage<Complex> ctemp(imShape, imageCoord);
293 665 : ctemp.set(1.0);
294 665 : pbMath.applyPB(ctemp, ctemp, wcenter, Quantity(0.0, "deg"), BeamSquint::NONE);
295 665 : StokesImageUtil::To(pbTemp, ctemp);
296 665 : LatticeLocker lock1(*(itsImages->pb()), FileLocker::Write);
297 665 : itsImages->pb()->copyData( (LatticeExpr<Float>)((*(itsImages->pb()))+pbTemp) );
298 :
299 665 : }//addPB
300 :
301 :
302 :
303 : ////////////////////Old vi/Vb version
304 :
305 0 : void SIMapper::degrid(VisBuffer& vb)
306 : {
307 : //LogIO os( LogOrigin("SIMapper","degrid",WHERE) );
308 :
309 0 : if(useViVb2_p)
310 0 : throw(AipsError("Programmer Error: using vi1 mode with vi2 constructor"));
311 : ///This should not be called even but heck let's ignore
312 0 : if(ft_p.null() and cft_p.null())
313 0 : return;
314 :
315 0 : Cube<Complex> origCube;
316 0 : origCube.assign(vb.modelVisCube());
317 :
318 0 : vb.setModelVisCube( Complex(0.0,0.0) );
319 :
320 0 : if( ! ft_p.null() ) { ft_p->get(vb); }
321 0 : if( ! cft_p.null() ) { cft_p->get(vb, cl_p); }
322 :
323 0 : vb.modelVisCube()+=origCube;
324 :
325 0 : }
326 :
327 :
328 829 : void SIMapper::finalizeDegrid()
329 : {
330 : //LogIO os( LogOrigin("SIMapper","finalizeDegrid",WHERE) );
331 829 : if(!useViVb2_p)
332 0 : ft_p->finalizeToVis();
333 : else
334 829 : ft2_p->finalizeToVis();
335 829 : }
336 :
337 :
338 17 : Bool SIMapper::getFTMRecord(Record& rec, const String diskimage)
339 : {
340 : //LogIO os( LogOrigin("SIMapper","getFTMRecord",WHERE) );
341 17 : if(ft_p.null() && !useViVb2_p)
342 0 : return false;
343 17 : if(ft2_p.null() && useViVb2_p)
344 0 : return false;
345 17 : String err;
346 17 : if(!useViVb2_p)
347 0 : return ft_p->toRecord(err, rec, true, diskimage);
348 17 : return ft2_p->toRecord(err, rec, true, diskimage);
349 : // rec = itsFTM->toRecord();
350 :
351 17 : }
352 17 : Bool SIMapper::getCLRecord(Record& rec)
353 : {
354 17 : if(cft_p.null() && cft2_p.null())
355 17 : return false;
356 0 : String err;
357 0 : return cl_p.toRecord(err, rec);
358 0 : }
359 :
360 :
361 631 : void SIMapper::initPB()
362 : {
363 631 : LatticeLocker lock1(*(itsImages->pb()), FileLocker::Write);
364 631 : itsImages->pb()->set(0.0);
365 631 : }
366 :
367 0 : void SIMapper::addPB(VisBuffer& vb, PBMath& pbMath)
368 : {
369 0 : CoordinateSystem imageCoord=itsImages->pb()->coordinates();
370 :
371 0 : IPosition imShape=itsImages->pb()->shape();
372 :
373 0 : MDirection wcenter=vb.msColumns().field().phaseDirMeas(vb.fieldId());
374 0 : TempImage<Float> pbTemp(imShape, imageCoord);
375 0 : TempImage<Complex> ctemp(imShape, imageCoord);
376 0 : ctemp.set(1.0);
377 0 : pbMath.applyPB(ctemp, ctemp, wcenter, Quantity(0.0, "deg"), BeamSquint::NONE);
378 0 : StokesImageUtil::To(pbTemp, ctemp);
379 0 : LatticeLocker lock1(*(itsImages->pb()), FileLocker::Write);
380 0 : itsImages->pb()->copyData( (LatticeExpr<Float>)((*(itsImages->pb()))+pbTemp) );
381 :
382 0 : }//addPB
383 :
384 :
385 : } //# NAMESPACE CASA - END
386 :
387 :
|