Line data Source code
1 : //# Copyright (C) 1998,1999,2000,2001,2003
2 : //# Associated Universities, Inc. Washington DC, USA.
3 : //#
4 : //# This program is free software; you can redistribute it and/or modify it
5 : //# under the terms of the GNU General Public License as published by the Free
6 : //# Software Foundation; either version 2 of the License, or (at your option)
7 : //# any later version.
8 : //#
9 : //# This program is distributed in the hope that it will be useful, but WITHOUT
10 : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11 : //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
12 : //# more details.
13 : //#
14 : //# You should have received a copy of the GNU General Public License along
15 : //# with this program; if not, write to the Free Software Foundation, Inc.,
16 : //# 675 Massachusetts Ave, Cambridge, MA 02139, USA.
17 : //#
18 : //# Correspondence concerning AIPS++ should be addressed as follows:
19 : //# Internet email: aips2-request@nrao.edu.
20 : //# Postal address: AIPS++ Project Office
21 : //# National Radio Astronomy Observatory
22 : //# 520 Edgemont Road
23 : //# Charlottesville, VA 22903-2475 USA
24 : //#
25 : //# $Id: $
26 :
27 : #include <imageanalysis/ImageAnalysis/ImageTask.h>
28 : #include <casacore/coordinates/Coordinates/DirectionCoordinate.h>
29 : #include <casacore/casa/IO/FilebufIO.h>
30 : #include <casacore/casa/OS/Directory.h>
31 : #include <casacore/casa/OS/RegularFile.h>
32 : #include <casacore/casa/OS/SymLink.h>
33 : #include <casacore/images/Images/FITSImage.h>
34 : #include <casacore/images/Images/ImageUtilities.h>
35 : #include <casacore/images/Images/MIRIADImage.h>
36 : #include <casacore/images/Images/PagedImage.h>
37 : #include <casacore/images/Images/TempImage.h>
38 : #include <casacore/tables/Tables/PlainTable.h>
39 :
40 : #include <imageanalysis/ImageAnalysis/ImageHistory.h>
41 : #include <imageanalysis/ImageAnalysis/ImageInputProcessor.h>
42 : #include <imageanalysis/ImageAnalysis/ImageMask.h>
43 : #include <imageanalysis/ImageAnalysis/SubImageFactory.h>
44 : #include <imageanalysis/IO/LogFile.h>
45 : #include <stdcasa/variant.h>
46 :
47 : namespace casa {
48 :
49 0 : template <class T> ImageTask<T>::ImageTask(
50 : const SPCIIT image,
51 : const casacore::String& region, const casacore::Record *const ®ionPtr,
52 : const casacore::String& box, const casacore::String& chanInp,
53 : const casacore::String& stokes, const casacore::String& maskInp,
54 : const casacore::String& outname, casacore::Bool overwrite
55 : ) : _image(image), _regionPtr(regionPtr),_region(region), _box(box),
56 : _chan(chanInp), _stokesString(stokes), _mask(maskInp),
57 : _outname(outname), _overwrite(overwrite), _stretch(false),
58 0 : _logfile() {}
59 :
60 0 : template <class T> ImageTask<T>::ImageTask(
61 : const SPCIIT image, const casacore::Record *const ®ionPtr,
62 : const casacore::String& mask, const casacore::String& outname,
63 : casacore::Bool overwrite
64 : ) : _image(image), _regionPtr(regionPtr),
65 : _region(), _box(), _chan(), _stokesString(), _mask(mask),
66 0 : _outname(outname), _overwrite(overwrite) {}
67 :
68 0 : template <class T> ImageTask<T>::~ImageTask() {}
69 :
70 0 : template <class T> std::vector<OutputDestinationChecker::OutputStruct> ImageTask<T>::_getOutputStruct() {
71 0 : std::vector<OutputDestinationChecker::OutputStruct> outputs;
72 0 : _outname.trim();
73 0 : if (! _outname.empty()) {
74 0 : OutputDestinationChecker::OutputStruct outputImage;
75 0 : outputImage.label = "output image";
76 0 : outputImage.outputFile = &_outname;
77 0 : outputImage.required = true;
78 0 : outputImage.replaceable = _overwrite;
79 0 : outputs.push_back(outputImage);
80 : }
81 0 : return outputs;
82 : }
83 :
84 0 : template <class T> void ImageTask<T>::_construct(casacore::Bool verbose) {
85 0 : ThrowIf(
86 : ! _supportsMultipleBeams() && _image->imageInfo().hasMultipleBeams(),
87 : "This application does not support images with multiple "
88 : "beams. Please convolve your image with a single beam "
89 : "and run this application using that image"
90 : );
91 0 : casacore::String diagnostics;
92 0 : std::vector<OutputDestinationChecker::OutputStruct> outputs = _getOutputStruct();
93 0 : std::vector<OutputDestinationChecker::OutputStruct> *outputPtr = outputs.size() > 0
94 : ? &outputs
95 : : 0;
96 0 : std::vector<casacore::Coordinate::Type> necCoords = _getNecessaryCoordinates();
97 0 : std::vector<casacore::Coordinate::Type> *coordsPtr = necCoords.size() > 0
98 : ? &necCoords
99 : : 0;
100 0 : ThrowIf(
101 : _mustHaveSquareDirectionPixels()
102 : && _image->coordinates().hasDirectionCoordinate()
103 : && ! _image->coordinates().directionCoordinate().hasSquarePixels(),
104 : "This application requires that the input image must have square "
105 : "direction pixels, but the input image does not. Please regrid it "
106 : "so it does and rerun on the regridded image"
107 : );
108 0 : ImageInputProcessor inputProcessor;
109 0 : inputProcessor.process(
110 0 : _regionRecord, diagnostics, outputPtr,
111 0 : _stokesString, _image, _regionPtr,
112 0 : _region, _box, _chan,
113 : _getStokesControl(), _supportsMultipleRegions(),
114 0 : coordsPtr, verbose
115 : );
116 0 : }
117 :
118 : template <class T> void ImageTask<T>::setRegion(const casacore::Record& region) {
119 : ThrowIf(
120 : ! _supportsMultipleRegions() && region.isDefined("regions"),
121 : "This application does not support multiple region selection"
122 : );
123 : _regionRecord = region;
124 : _box = "";
125 : _chan = "";
126 : _stokesString = "";
127 : _region = "";
128 : }
129 :
130 0 : template <class T> void ImageTask<T>::_removeExistingFileIfNecessary(
131 : const casacore::String& filename, casacore::Bool overwrite, casacore::Bool warnOnly
132 : ) const {
133 0 : casacore::File out(filename);
134 0 : if (out.exists()) {
135 0 : if (overwrite) {
136 0 : File f(filename);
137 0 : ThrowIf(
138 : PlainTable::tableCache()(f.path().absoluteName()),
139 : filename + " is currently present in the table cache "
140 : + "and so is being used by another process. Please close "
141 : + "it in the other process first before attempting to "
142 : + "overwrite it"
143 : );
144 0 : if (out.isDirectory()) {
145 0 : casacore::Directory dir(filename);
146 0 : dir.removeRecursive();
147 : }
148 0 : else if (out.isRegular()) {
149 0 : casacore::RegularFile reg(filename);
150 0 : reg.remove();
151 : }
152 0 : else if (out.isSymLink()) {
153 0 : casacore::SymLink link(filename);
154 0 : link.remove();
155 : }
156 : }
157 : else {
158 0 : casacore::String msg = "File " + filename + " exists but overwrite is false "
159 : "so it cannot be overwritten";
160 0 : if (warnOnly) {
161 0 : *_log << casacore::LogIO::WARN << msg << casacore::LogIO::POST;
162 : }
163 : else {
164 0 : ThrowCc(msg);
165 : }
166 : }
167 : }
168 0 : }
169 :
170 : template <class T> void ImageTask<T>::_removeExistingOutfileIfNecessary() const {
171 : _removeExistingFileIfNecessary(_outname, _overwrite);
172 : }
173 :
174 0 : template <class T> casacore::String ImageTask<T>::_summaryHeader() const {
175 0 : casacore::String region = _box.empty() ? _region : "";
176 0 : ostringstream os;
177 0 : os << "Input parameters ---" << endl;
178 0 : os << " --- imagename: " << _image->name() << endl;
179 0 : os << " --- region: " << region << endl;
180 0 : os << " --- box: " << _box << endl;
181 0 : os << " --- channels: " << _chan << endl;
182 0 : os << " --- stokes: " << _stokesString << endl;
183 0 : os << " --- mask: " << _mask << endl;
184 0 : return os.str();
185 : }
186 :
187 : template <class T> void ImageTask<T>::setLogfile(const casacore::String& lf) {
188 : if (lf.empty()) {
189 : return;
190 : }
191 : ThrowIf(
192 : ! _hasLogfileSupport(),
193 : "Logic Error: This task does not support writing of a log file"
194 : );
195 : try {
196 : _logfile.reset(new LogFile(lf));
197 : _logfile->setAppend(_logfileAppend);
198 : }
199 : catch (const casacore::AipsError& x) {}
200 : }
201 :
202 0 : template <class T> const std::shared_ptr<LogFile> ImageTask<T>::_getLogFile() const {
203 0 : ThrowIf(
204 : ! _hasLogfileSupport(),
205 : "Logic Error: This task does not support writing of a log file"
206 : );
207 0 : return _logfile;
208 : }
209 :
210 : template <class T> casacore::Bool ImageTask<T>::_openLogfile() {
211 : if (_logfile.get() == 0) {
212 : return false;
213 : }
214 : ThrowIf(
215 : ! _hasLogfileSupport(),
216 : "Logic Error: This task does not support writing of a log file"
217 : );
218 : return _logfile->open();
219 : }
220 :
221 : template <class T> void ImageTask<T>::_closeLogfile() const {
222 : if (_logfile) {
223 : _logfile->close();
224 : }
225 : }
226 :
227 : template<class T> casacore::Bool ImageTask<T>::_writeLogfile(
228 : const casacore::String& output, const casacore::Bool open, const casacore::Bool close
229 : ) {
230 : ThrowIf(
231 : ! _hasLogfileSupport(),
232 : "Logic Error: This task does not support writing of a log file"
233 : );
234 : if (! _logfile) {
235 : return false;
236 : }
237 : return _logfile->write(output, open, close);
238 : }
239 :
240 : template <class T> void ImageTask<T>::setLogfileAppend(casacore::Bool a) {
241 : ThrowIf(
242 : ! _hasLogfileSupport(),
243 : "Logic Error: This task does not support writing of a log file"
244 : );
245 : _logfileAppend = a;
246 : if (_logfile) {
247 : _logfile->setAppend(a);
248 : }
249 : }
250 :
251 0 : template <class T> void ImageTask<T>::addHistory(
252 : const vector<std::pair<casacore::String, casacore::String> >& msgs
253 : ) const {
254 0 : _newHistory.insert(
255 0 : _newHistory.end(), msgs.begin(), msgs.end()
256 : );
257 0 : }
258 :
259 0 : template <class T> void ImageTask<T>::addHistory(
260 : const casacore::LogOrigin& origin, const casacore::String& msg
261 : ) const {
262 0 : std::pair<casacore::String, casacore::String> x;
263 0 : x.first = origin.fullName();
264 0 : x.second = msg;
265 0 : _newHistory.push_back(x);
266 0 : }
267 :
268 0 : template <class T> void ImageTask<T>::addHistory(
269 : const casacore::LogOrigin& origin, const vector<casacore::String>& msgs
270 : ) const {
271 0 : std::pair<casacore::String, casacore::String> x;
272 0 : x.first = origin.fullName();
273 0 : for( casacore::String m: msgs ) {
274 0 : x.second = m;
275 0 : _newHistory.push_back(x);
276 : }
277 0 : }
278 :
279 : template <class T> void ImageTask<T>::addHistory(
280 : const casacore::LogOrigin& origin, const casacore::String& taskname,
281 : const vector<casacore::String>& paramNames, const vector<casac::variant>& paramValues
282 : ) const {
283 : auto appHistory = ImageHistory<T>::getApplicationHistory(
284 : origin, taskname, paramNames, paramValues, _image->name()
285 : );
286 : _newHistory.insert(_newHistory.end(), appHistory.begin(), appHistory.end());
287 : }
288 :
289 0 : template <class T> void ImageTask<T>::_copyMask(
290 : casacore::Lattice<casacore::Bool>& mask, const casacore::ImageInterface<T>& image
291 : ) {
292 0 : auto cursorShape = image.niceCursorShape(4096*4096);
293 0 : casacore::LatticeStepper stepper(image.shape(), cursorShape, casacore::LatticeStepper::RESIZE);
294 0 : casacore::RO_MaskedLatticeIterator<T> iter(image, stepper);
295 0 : casacore::LatticeIterator<casacore::Bool> miter(mask, stepper);
296 0 : std::unique_ptr<casacore::RO_LatticeIterator<casacore::Bool>> pmiter;
297 0 : if (image.hasPixelMask()) {
298 0 : pmiter.reset(new casacore::RO_LatticeIterator<casacore::Bool>(image.pixelMask(), stepper));
299 : }
300 0 : for (iter.reset(); ! iter.atEnd(); ++iter, ++miter) {
301 0 : auto mymask = iter.getMask();
302 0 : if (pmiter) {
303 0 : mymask = mymask && pmiter->cursor();
304 0 : pmiter->operator++();
305 : }
306 0 : miter.rwCursor() = mymask;
307 : }
308 0 : }
309 :
310 0 : template <class T> void ImageTask<T>::_copyData(
311 : Lattice<T>& data, const Lattice<T>& image
312 : ) {
313 0 : auto cursorShape = image.niceCursorShape(4096*4096);
314 0 : casacore::LatticeStepper stepper(image.shape(), cursorShape, casacore::LatticeStepper::RESIZE);
315 0 : casacore::RO_LatticeIterator<T> iter(image, stepper);
316 0 : casacore::LatticeIterator<T> diter(data, stepper);
317 0 : for (iter.reset(); ! iter.atEnd(); ++iter, ++diter) {
318 0 : diter.rwCursor() = iter.cursor();
319 : }
320 0 : }
321 :
322 0 : template <class T> casacore::Bool ImageTask<T>::_isPVImage() const {
323 0 : const CoordinateSystem& csys = _image->coordinates();
324 0 : if (csys.hasLinearCoordinate() && csys.hasSpectralAxis()) {
325 0 : auto pixelAxes = csys.linearAxesNumbers();
326 0 : auto nPixelAxes = pixelAxes.size();
327 0 : auto names = csys.worldAxisNames();
328 0 : for (uInt j=0; j<nPixelAxes; ++j) {
329 0 : if (pixelAxes[j] >= 0 && names[pixelAxes[j]] == "Offset") {
330 0 : return true;
331 : }
332 : }
333 : }
334 0 : return false;
335 : }
336 :
337 0 : template <class T> SPIIT ImageTask<T>::_prepareOutputImage(
338 : const casacore::ImageInterface<T>& image, const casacore::Array<T> *const values,
339 : const casacore::ArrayLattice<casacore::Bool> *const mask,
340 : const casacore::IPosition *const outShape, const casacore::CoordinateSystem *const coordsys,
341 : const casacore::String *const outname, casacore::Bool overwrite, casacore::Bool dropDegen
342 : ) const {
343 0 : auto oShape = outShape == 0 ? image.shape() : *outShape;
344 0 : casacore::CoordinateSystem csys = coordsys ? *coordsys : image.coordinates();
345 0 : std::shared_ptr<casacore::TempImage<T>> tmpImage(
346 0 : new casacore::TempImage<T>(casacore::TiledShape(oShape), csys)
347 : );
348 0 : if (mask && (! ImageMask::isAllMaskTrue(*mask))) {
349 0 : tmpImage->attachMask(*mask);
350 : }
351 : // because subimages can have two types of masks, a region mask and
352 : // a pixel mask, but most other types of images really just have a
353 : // pixel mask. its very confusing
354 0 : else if (image.hasPixelMask() || image.isMasked()) {
355 : // A paged array is stored on disk and is preferred over an
356 : // ArrayLattice which will exhaust memory for large images.
357 0 : std::unique_ptr<casacore::Lattice<casacore::Bool>> mymask;
358 : const static auto mmasksize = 4096*4096;
359 0 : if (image.size() > mmasksize) {
360 0 : mymask.reset(new casacore::PagedArray<casacore::Bool>(image.shape()));
361 : }
362 : else {
363 0 : mymask.reset(new casacore::ArrayLattice<casacore::Bool>(image.shape()));
364 : }
365 0 : _copyMask(*mymask, image);
366 0 : if (! ImageMask::isAllMaskTrue(image)) {
367 0 : tmpImage->attachMask(*mymask);
368 : }
369 : }
370 0 : auto myOutname = outname ? *outname : _outname;
371 0 : if (! outname) {
372 0 : overwrite = _overwrite;
373 : }
374 0 : SPIIT outImage = tmpImage;
375 0 : if (values) {
376 0 : outImage->put(*values);
377 : }
378 : else {
379 : // FIXME this is a superfluous copy if dropgen || ! myOutname.empty()
380 0 : _copyData(*outImage, image);
381 : }
382 0 : if (dropDegen || ! myOutname.empty()) {
383 0 : if (! myOutname.empty()) {
384 0 : _removeExistingFileIfNecessary(myOutname, overwrite);
385 : }
386 0 : casacore::String emptyMask = "";
387 0 : casacore::Record empty;
388 0 : outImage = SubImageFactory<T>::createImage(
389 0 : *tmpImage, myOutname, empty, emptyMask,
390 : dropDegen, false, true, false
391 : );
392 : }
393 0 : casacore::ImageUtilities::copyMiscellaneous(*outImage, image);
394 0 : _doHistory(outImage);
395 : // CAS-9267 force metadata to be written to disk, in case of PagedImage
396 0 : outImage->flush();
397 0 : return outImage;
398 : }
399 :
400 0 : template <class T> SPIIT ImageTask<T>::_prepareOutputImage(
401 : const casacore::ImageInterface<T>& image, casacore::Bool dropDeg
402 : ) const {
403 0 : if (! _outname.empty()) {
404 0 : _removeExistingFileIfNecessary(_outname, _overwrite);
405 : }
406 0 : static const casacore::Record empty;
407 0 : static const casacore::String emptyString;
408 0 : auto outImage = SubImageFactory<T>::createImage(
409 0 : image, _outname, empty, emptyString, dropDeg,
410 0 : _overwrite, true, false, false
411 : );
412 0 : _doHistory(outImage);
413 0 : return outImage;
414 : }
415 :
416 0 : template <class T> SPIIT ImageTask<T>::_prepareOutputImage(
417 : const casacore::ImageInterface<T>& image, const casacore::String& outname,
418 : casacore::Bool overwrite, casacore::Bool warnOnly
419 : ) const {
420 0 : if (! outname.empty()) {
421 0 : _removeExistingFileIfNecessary(outname, overwrite, warnOnly);
422 : }
423 0 : static const casacore::Record empty;
424 0 : static const casacore::String emptyString;
425 0 : auto outImage = SubImageFactory<T>::createImage(
426 : image, outname, empty, emptyString, false,
427 : overwrite, true, false, false
428 : );
429 0 : _doHistory(outImage);
430 0 : return outImage;
431 : }
432 :
433 :
434 0 : template <class T> SPIIT ImageTask<T>::_prepareOutputImage(
435 : const casacore::ImageInterface<T>& image, const casacore::Lattice<T>& data
436 : ) const {
437 0 : if (! _outname.empty()) {
438 0 : _removeExistingFileIfNecessary(_outname, _overwrite);
439 : }
440 0 : static const casacore::Record empty;
441 0 : static const casacore::String emptyString;
442 0 : auto outImage = SubImageFactory<T>::createImage(
443 0 : image, _outname, empty, emptyString, false,
444 0 : _overwrite, true, false, false, &data
445 : );
446 0 : _doHistory(outImage);
447 0 : return outImage;
448 : }
449 :
450 : template <class T>
451 0 : template <class U> void ImageTask<T>::_doHistory(std::shared_ptr<casacore::ImageInterface<U>>& image) const {
452 0 : if (! _suppressHistory) {
453 0 : ImageHistory<U> history(image);
454 0 : if (history.get(false).empty()) {
455 0 : history.append(_image);
456 : }
457 0 : history.addHistory(_newHistory);
458 : /*
459 : for (const auto& line: _newHistory) {
460 : history.addHistory(line.first, line.second);
461 : }
462 : */
463 : }
464 0 : }
465 :
466 0 : template <class T> void ImageTask<T>::_reportOldNewImageShapes(const IPosition& outShape) const {
467 0 : LogOrigin lor(getClass(), __func__);
468 0 : ostringstream os;
469 0 : os << "Original " << _getImage()->name() << " size => "
470 0 : << _getImage()->shape();
471 0 : addHistory(lor, os.str());
472 0 : *_getLog() << LogIO::NORMAL << os.str() << LogIO::POST;
473 0 : os.str("");
474 0 : os << "New " << _getOutname() << " size => " << outShape;
475 0 : addHistory(lor, os.str());
476 0 : *_getLog() << LogIO::NORMAL << os.str() << LogIO::POST;
477 0 : }
478 :
479 : template <class T> void ImageTask<T>::_reportOldNewImageShapes(const ImageInterface<T>& out) const {
480 : _reportOldNewImageShapes(out.shape());
481 : }
482 :
483 : }
484 :
|