Commits

Pam Ford authored 3f05bc82a51 Merge
Merge branch 'master' into feature/CAS-10284

code/synthesis/ImagerObjects/test/tLabelandFindRegions.cc

Added
1 +/*
2 + * tSynthesisUtils.cc
3 + *demo of SynthesisUtil functionality
4 +//# Copyright (C) 2013-2014
5 +//# Associated Universities, Inc. Washington DC, USA.
6 +//#
7 +//# This program is free software; you can redistribute it and/or modify it
8 +//# under the terms of the GNU General Public License as published by the Free
9 +//# Software Foundation; either version 2 of the License, or (at your option)
10 +//# any later version.
11 +//#
12 +//# This program is distributed in the hope that it will be useful, but WITHOUT
13 +//# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14 +//# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
15 +//# more details.
16 +//#
17 +//# You should have received a copy of the GNU General Public License along
18 +//# with this program; if not, write to the Free Software Foundation, Inc.,
19 +//# 675 Massachusetts Ave, Cambridge, MA 02139, USA.
20 +//#
21 +//# Correspondence concerning AIPS++ should be addressed as follows:
22 +//# Internet email: aips2-request@nrao.edu.
23 +//# Postal address: AIPS++ Project Office
24 +//# National Radio Astronomy Observatory
25 +//# 520 Edgemont Road
26 +//# Charlottesville, VA 22903-2475 USA
27 + *
28 + * test labelRegions and findBlobSize
29 + */
30 +
31 +
32 +
33 +
34 +
35 +#include <casacore/images/Images/ImageUtilities.h>
36 +#include <casa/iostream.h>
37 +#include <casa/aips.h>
38 +#include <casa/Exceptions/Error.h>
39 +#include <casa/BasicSL/String.h>
40 +#include <casa/Containers/Block.h>
41 +#include <measures/Measures/MRadialVelocity.h>
42 +#include <coordinates/Coordinates/CoordinateSystem.h>
43 +#include <casa/Logging/LogIO.h>
44 +#include <synthesis/ImagerObjects/SynthesisImager.h>
45 +#include <synthesis/ImagerObjects/SIImageStore.h>
46 +#include <imageanalysis/Utilities/SpectralImageUtil.h>
47 +#include <lattices/Lattices/LatticeConcat.h>
48 +#include <images/Images/PagedImage.h>
49 +#include <images/Images/ImageConcat.h>
50 +#include <images/Images/SubImage.h>
51 +#include <casa/namespace.h>
52 +#include <images/Images/TempImage.h>
53 +#include <coordinates/Coordinates/CoordinateUtil.h>
54 +#include <ms/MSSel/MSSourceIndex.h>
55 +#include <synthesis/ImagerObjects/SDMaskHandler.h>
56 +
57 +int main(int argc, char **argv)
58 +{
59 + using namespace std;
60 + using namespace casacore;
61 + using namespace casa;
62 + try{
63 +
64 +
65 + String imagename;
66 + if (argc > 1) {
67 + imagename=argv[1];
68 + }
69 + else {
70 + imagename="tLabelandFindRegions.temp.mask";
71 + Int imsize = 100;
72 + Int nchan = 2;
73 + IPosition blc(4,45,50,0,0);
74 + IPosition trc(4,55,55,0,0);
75 + IPosition shape(4, imsize, imsize, 1, nchan);
76 + CoordinateSystem csys=CoordinateUtil::defaultCoords4D();
77 + PagedImage<Float> maskImage(TiledShape(shape), csys, imagename);
78 + maskImage.setUnits(Unit("Jy/pixel"));
79 + maskImage.set(0.0);
80 + // sanity check
81 + if (blc(0) <= trc(0) && blc(1) <= trc(1) && blc(3) < nchan) {
82 + Int dx = trc(0) - blc(0);
83 + Int dy = trc(1) - blc(1);
84 + for (uInt i=0+blc(3); i < trc(3)+1; i++) {
85 + for (uInt j=0; j < (uInt)dx+1; j++) {
86 + for (uInt k=0; k < (uInt)dy+1; k++) {
87 + IPosition loc(4,blc(0)+j, blc(1)+k, 0, Int(i));
88 + Float val = 1.0;
89 + maskImage.putAt(val, loc);
90 + }
91 + }
92 + }
93 + }
94 + }
95 +
96 + CountedPtr<ImageInterface<Float> > im;
97 + im = ImageUtilities::openImage<Float>(imagename);
98 + IPosition shp=im->shape();
99 + cerr<<"shp="<<shp<<endl;
100 + IPosition start(4, 0, 0, 0, 0);
101 + IPosition length(4, shp(0),shp(1),shp(2),1);
102 + Slicer sl(start, length);
103 + AxesSpecifier aspec(False);
104 + SubImage<Float>* subIm = new SubImage<Float>(*im, sl, aspec, True);
105 + cerr<<"subIm shape="<<subIm->shape()<<" slice="<<sl<<endl;
106 + SDMaskHandler smh;
107 +
108 + TempImage<Float>* tempIm = new TempImage<Float> (TiledShape(subIm->shape(), subIm->niceCursorShape()), subIm->coordinates() );
109 + tempIm->copyData(LatticeExpr<Float> (*subIm));
110 + cerr<<"tempIm shape="<<tempIm->shape()<<endl;
111 +
112 + cerr<<"setup a blobmap ..."<<endl;
113 + TempImage<Float>* blobMap = new TempImage<Float> (TiledShape(subIm->shape(), subIm->niceCursorShape()), subIm->coordinates() );
114 + blobMap->set(0);
115 + cerr<<"blobMap shape="<<blobMap->shape()<<endl;
116 +
117 + cerr<<"Calling labelRegions()..."<<endl;
118 + //Array<Float> subimarr;
119 + //Array<Float> blobarr;
120 + //smh.getarr(*subIm, False);
121 + //smh.getarr(*blobMap, True);
122 + smh.labelRegions(*subIm, *blobMap);
123 + cerr<<"Calling findBlobSize()..."<<endl;
124 + Vector<Float> bloblist=smh.findBlobSize(*blobMap);
125 + for (uInt iblob=0; iblob<bloblist.nelements(); ++iblob) {
126 + cerr<<iblob<<" : "<<bloblist[iblob]<<endl;
127 + }
128 +
129 + }catch( AipsError e ){
130 + cout << "Exception ocurred." << endl;
131 + cout << e.getMesg() << endl;
132 + }
133 + return 0;
134 +};

Everything looks good. We'll let you know here if there's anything you should know about.

Add shortcut