Line data Source code
1 : // -*- C++ -*-
2 : //# AWProjectWBFTHPG.cc: Implementation of AWProjectWBFTHPG class
3 : //# Copyright (C) 2021
4 : //# Associated Universities, Inc. Washington DC, USA.
5 : //#
6 : //# This library is free software; you can redistribute it and/or modify it
7 : //# under the terms of the GNU Library General Public License as published by
8 : //# the Free Software Foundation; either version 2 of the License, or (at your
9 : //# option) any later version.
10 : //#
11 : //# This library is distributed in the hope that it will be useful, but WITHOUT
12 : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 : //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
14 : //# License for more details.
15 : //#
16 : //# You should have received a copy of the GNU Library General Public License
17 : //# along with this library; if not, write to the Free Software Foundation,
18 : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
19 : //#
20 : //# Correspondence concerning AIPS++ should be addressed as follows:
21 : //# Internet email: aips2-request@nrao.edu.
22 : //# Postal address: AIPS++ Project Office
23 : //# National Radio Astronomy Observatory
24 : //# 520 Edgemont Road
25 : //# Charlottesville, VA 22903-2475 USA
26 : //#
27 : //# $Id$
28 :
29 : #include <synthesis/TransformMachines2/AWProjectWBFTHPG.h>
30 : #include <casacore/coordinates/Coordinates/CoordinateSystem.h>
31 : #include <synthesis/ImagerObjects/SIImageStore.h>
32 : #include <synthesis/ImagerObjects/SimpleSIImageStore.h>
33 : #include <synthesis/TransformMachines/StokesImageUtil.h>
34 : #include <synthesis/TransformMachines2/AWVisResampler.h>
35 : #include <synthesis/TransformMachines2/SimplePBConvFunc.h>
36 : #include <casacore/casa/Arrays/Array.h>
37 : #include <casacore/casa/Arrays/ArrayMath.h>
38 : #include <casacore/casa/Arrays/Slice.h>
39 : #include <casacore/casa/Arrays/Vector.h>
40 : #include <casacore/casa/OS/HostInfo.h>
41 : #include <casacore/casa/Utilities/CompositeNumber.h>
42 : #include <casacore/images/Images/ImageInterface.h>
43 : #include <casacore/images/Images/PagedImage.h>
44 : #include <msvis/MSVis/VisBuffer2.h>
45 : #include <sstream>
46 :
47 : using namespace casacore;
48 : namespace casa { //# NAMESPACE CASA - BEGIN
49 : namespace refim {
50 : //---------------------------------------------------------------
51 : //
52 0 : ImageInterface<Complex> &AWProjectWBFTHPG::getImage(Matrix<Float> &weights,
53 : Bool normalize) {
54 0 : LogIO log_l(LogOrigin("AWProjectWBFTHPG", "getImage[R&D]"));
55 : //
56 0 : AlwaysAssert(image, AipsError);
57 :
58 0 : weights.resize(sumWeight.shape());
59 0 : convertArray(weights, sumWeight);
60 : //
61 : // If the weights are all zero then we cannot normalize otherwise
62 : // we don't care.
63 : //
64 0 : if (max(weights) == 0.0)
65 0 : log_l << "No useful data in " << name() << ". Weights all zero"
66 0 : << LogIO::POST;
67 : else {
68 : /*log_l << "Sum of weights: " << weights << " ";
69 : if (griddedData2.nelements() > 0) {
70 : log_l << max(griddedData2) << " " << min(griddedData2);
71 : };
72 : log_l << LogIO::POST;*/
73 0 : cerr << "Sum of weights: " << setprecision(20) << weights << endl;
74 : }
75 :
76 : //
77 : // x and y transforms (lattice has the gridded vis. Make the
78 : // dirty images)
79 : //
80 0 : if (useDoubleGrid_p) {
81 : // ArrayLattice<DComplex> darrayLattice(griddedData2);
82 : // {
83 : // griddedData.resize(griddedData2.shape());
84 : // convertArray(griddedData, griddedData2);
85 : // storeArrayAsImage(String("cgrid_"+visResampler_p->name()+".im"),
86 : // image->coordinates(), griddedData);
87 : // }
88 : // LatticeFFT::cfft2d(darrayLattice,false);
89 :
90 0 : griddedData.resize(griddedData2.shape());
91 0 : convertArray(griddedData, griddedData2);
92 0 : SynthesisUtilMethods::getResource("mem peak in getImage");
93 :
94 : // Don't need the double-prec grid anymore...
95 0 : griddedData2.resize();
96 0 : lattice = new ArrayLattice<Complex>(griddedData);
97 : } else {
98 0 : lattice = new ArrayLattice<Complex>(griddedData);
99 : // // cerr << "##### " << griddedData2.shape() << endl;
100 : // lattice=arrayLattice;
101 : // LatticeFFT::cfft2d(*lattice,false);
102 : }
103 :
104 : //
105 : // Now normalize the dirty image.
106 : //
107 : // Since *lattice is not copied to *image till the end of this
108 : // method, normalizeImage also needs to work with Lattices
109 : // (rather than ImageInterface).
110 : //
111 : // //normalizeImage(*lattice,sumWeight,*avgPB_p,fftNormalization);
112 : // normalizeImage(*lattice,sumWeight,*avgPB_p, *avgPBSq_p,
113 : //fftNormalization);
114 :
115 : // nx ny normalization from GridFT...
116 : {
117 0 : Int inx = lattice->shape()(0);
118 0 : Int iny = lattice->shape()(1);
119 0 : Vector<Complex> correction(inx);
120 0 : correction = Complex(1.0, 0.0);
121 0 : Vector<Float> sincConvX(inx);
122 0 : for (Int ix = 0; ix < inx; ix++) {
123 0 : Float x = M_PI * Float(ix - inx / 2) / (Float(nx) * Float(convSampling));
124 0 : if (ix == inx / 2) {
125 0 : sincConvX(ix) = 1.0;
126 : } else {
127 0 : sincConvX(ix) = sin(x) / x;
128 : }
129 : }
130 0 : Vector<Float> sincConvY(iny);
131 0 : for (Int ix = 0; ix < iny; ix++) {
132 0 : Float x = M_PI * Float(ix - iny / 2) / (Float(ny) * Float(convSampling));
133 0 : if (ix == iny / 2) {
134 0 : sincConvY(ix) = 1.0;
135 : } else {
136 0 : sincConvY(ix) = sin(x) / x;
137 : }
138 : }
139 :
140 : //cerr << "NORM " << normalize << " min correction " << min(sincConvX) << " " << min(sincConvY) << endl;
141 :
142 : // Do the Grid-correction
143 0 : IPosition cursorShape(4, inx, 1, 1, 1);
144 0 : IPosition axisPath(4, 0, 1, 2, 3);
145 0 : LatticeStepper lsx(lattice->shape(), cursorShape, axisPath);
146 0 : LatticeIterator<Complex> lix(*lattice, lsx);
147 :
148 0 : for (lix.reset(); !lix.atEnd(); lix++) {
149 0 : Int pol = lix.position()(2);
150 0 : Int chan = lix.position()(3);
151 0 : if (weights(pol, chan) != 0.0) {
152 0 : Int iy = lix.position()(1);
153 0 : for (Int ix = 0; ix < nx; ix++)
154 0 : correction(ix) = 1 / (sincConvX(ix) * sincConvY(iy));
155 : // cerr <<"Min max correction " << min(correction) << " " <<
156 : // max(correction) << endl;
157 0 : lix.rwVectorCursor() *= correction;
158 0 : if (normalize) {
159 0 : Complex rnorm(Float(inx) * Float(iny) / weights(pol, chan));
160 0 : lix.rwCursor() *= rnorm;
161 : }
162 : } else {
163 0 : lix.woCursor() = 0.0;
164 : }
165 : }
166 :
167 : // for(lix.reset();!lix.atEnd();lix++)
168 : // {
169 : // Int pol=lix.position()(2);
170 : // Int chan=lix.position()(3);
171 : // if (normalize)
172 : // {
173 : // if(weights(pol,chan)!=0.0)
174 : // {
175 : // Complex rnorm(Float(inx)*Float(iny)/(sincConv(inx)*sincConv(iny)*
176 : // weights(pol,chan) )); lix.rwCursor()*=rnorm;
177 : // }
178 : // else
179 : // lix.woCursor()=0.0;
180 : // }
181 : // else
182 : // lix.rwCursor() /= sincConv(inx)*sincConv(iny);
183 : // }
184 0 : }
185 0 : if (!isTiled) {
186 : //
187 : // Check the section from the image BEFORE converting to a lattice
188 : //
189 0 : IPosition blc(4, (nx - image->shape()(0) + (nx % 2 == 0)) / 2,
190 0 : (ny - image->shape()(1) + (ny % 2 == 0)) / 2, 0, 0);
191 0 : IPosition stride(4, 1);
192 0 : IPosition trc(blc + image->shape() - stride);
193 : //
194 : // Do the copy
195 : //
196 : //cerr << "blc" << blc << " trc " << trc << " min max " << min(griddedData) << " max " << max(griddedData) << endl;
197 0 : image->put(griddedData(blc, trc));
198 :
199 0 : if (!lattice.null())
200 0 : lattice = 0;
201 0 : griddedData.resize(IPosition(1, 0));
202 0 : }
203 : {
204 : // TempImage<Complex> tt(lattice->shape(), image->coordinates());
205 : // tt.put(lattice->get());
206 : // storeImg(String("uvgrid"+visResampler_p->name()+".im"), *image,true);
207 : }
208 :
209 0 : return *image;
210 0 : }
211 : //////////////////////////////////////////
212 :
213 0 : void AWProjectWBFTHPG::getWeightImage(ImageInterface<Float> &weightImage,
214 : Matrix<Float> &weights) {
215 : /// This is a mess nothing is initialized properly..
216 : // lets make a guess :)
217 : // SB: Review the guess below (not by SB) to see if it is still
218 : // required after the code cleanup.
219 0 : if (avgPB_p && (avgPB_p->shape()).size() == 0)
220 0 : avgPB_p = nullptr;
221 0 : if (!avgPB_p) {
222 0 : (getImage(weights, false));
223 0 : IPosition cursorShape(4, image->shape()[0], image->shape()[1], 1, 1);
224 0 : IPosition axisPath(4, 0, 1, 2, 3);
225 0 : LatticeStepper lsx(image->shape(), cursorShape, axisPath);
226 0 : LatticeIterator<Complex> lix(*image, lsx);
227 0 : for (lix.reset(); !lix.atEnd(); lix++) {
228 0 : Int pol = lix.position()(2);
229 0 : Int chan = lix.position()(3);
230 0 : if (weights(pol, chan) != 0.0) {
231 0 : Complex rnorm(1 / weights(pol, chan));
232 0 : lix.rwCursor() *= rnorm;
233 : } else {
234 0 : lix.rwCursor() = 0;
235 : }
236 : }
237 0 : StokesImageUtil::ToStokesPSF(weightImage, *image);
238 0 : setWeightImage(weightImage);
239 0 : } else {
240 0 : weightImage.resize(avgPB_p->shape());
241 0 : weightImage.copyData(*avgPB_p);
242 0 : weights.resize(sumWeight.shape());
243 0 : convertArray(weights, sumWeight);
244 : }
245 0 : avgPBReady_p = True;
246 0 : }
247 : //
248 : //---------------------------------------------------------------
249 : // Methods to accumulate data on the grid. These trigger gridding
250 : // to make *one* type of image: weights, PSF or residual depending
251 : // on the vbs.ftmType_p. This is required with HPG since it grids
252 : // for only one type of image at a time for efficiency. This is
253 : // different from the pattern on the CPU where weights are gridded
254 : // along with the gridding for the first image (PSF or residual).
255 : // Hence the specialization here.
256 : //
257 0 : void AWProjectWBFTHPG::resampleDataToGrid(Array<Complex> &griddedData_l,
258 : VBStore &vbs, const VisBuffer2 &vb,
259 : Bool &dopsf) {
260 0 : if (vbs.ftmType_p != casa::refim::FTMachine::WEIGHT)
261 0 : AWProjectFT::resampleDataToGrid(griddedData_l, vbs, vb, dopsf);
262 : // if (!avgPBReady_p)
263 : else {
264 : //
265 : // Get a reference to the pixels of griddedWeights (a
266 : // TempImage!)
267 : //
268 : //Array<Complex> gwts;
269 : //Bool removeDegenerateAxis = false;
270 : //griddedWeights.get(gwts, removeDegenerateAxis);
271 : //resampleCFToGrid(gwts, vbs, vb);
272 0 : vbs.ftmType_p=casa::refim::FTMachine::WEIGHT;
273 0 : Int nDataChan = vbs.flagCube_p.shape()[1];
274 0 : vbs.startChan_p = 0; vbs.endChan_p = nDataChan;
275 0 : Bool locdopsf=true;
276 0 : AWProjectFT::resampleDataToGrid(griddedData_l, vbs, vb, locdopsf);
277 :
278 : }
279 0 : };
280 : //
281 : //---------------------------------------------------------------
282 : //
283 0 : void AWProjectWBFTHPG::resampleDataToGrid(Array<DComplex> &griddedData_l,
284 : VBStore &vbs, const VisBuffer2 &vb,
285 : Bool &dopsf) {
286 0 : if (vbs.ftmType_p != casa::refim::FTMachine::WEIGHT)
287 0 : AWProjectFT::resampleDataToGrid(griddedData_l, vbs, vb, dopsf);
288 : // if (!avgPBReady_p)
289 : else {
290 : //
291 : // Get a reference to the pixels of griddedWeights (a
292 : // TempImage!)
293 : //
294 :
295 : //Array<DComplex> gwts;
296 : //Bool removeDegenerateAxis = false;
297 : //griddedWeights_D.get(gwts, removeDegenerateAxis);
298 : //resampleCFToGrid(griddedData_l, vbs, vb);
299 0 : vbs.ftmType_p=casa::refim::FTMachine::WEIGHT;
300 0 : Int nDataChan = vbs.flagCube_p.shape()[1];
301 0 : vbs.startChan_p = 0; vbs.endChan_p = nDataChan;
302 0 : Bool locdopsf=true;
303 0 : AWProjectFT::resampleDataToGrid(griddedData_l, vbs, vb, locdopsf);
304 : }
305 0 : };
306 : //
307 : //---------------------------------------------------------------
308 : //
309 0 : void AWProjectWBFTHPG::initializeToVisNew(const VisBuffer2 &vb,
310 : CountedPtr<SIImageStore> imstore) {
311 :
312 0 : Matrix<Float> tempWts;
313 :
314 0 : if (!(imstore->forwardGrid()).get())
315 0 : throw(AipsError("FTMAchine::InitializeToVisNew error imagestore has no "
316 0 : "valid grid initialized"));
317 : // Convert from Stokes planes to Correlation planes
318 0 : LatticeLocker lock1(*(imstore->model()), FileLocker::Read);
319 : // cerr << "###Max of imstore-> model " << max((imstore->model())->get())
320 : // << endl;
321 0 : stokesToCorrelation(*(imstore->model()), *(imstore->forwardGrid()));
322 :
323 0 : if (vb.polarizationFrame() == MSIter::Linear) {
324 0 : StokesImageUtil::changeCStokesRep(*(imstore->forwardGrid()),
325 : StokesImageUtil::LINEAR);
326 : } else {
327 0 : StokesImageUtil::changeCStokesRep(*(imstore->forwardGrid()),
328 : StokesImageUtil::CIRCULAR);
329 : }
330 0 : setFTMType(refim::FTMachine::RESIDUAL);
331 0 : visResampler_p->setModelImage((imstore->forwardGrid()));
332 0 : }
333 : //-------------------------------------------------------------------------
334 : //
335 0 : void AWProjectWBFTHPG::setupVBStore(VBStore& vbs,
336 : const VisBuffer2& vb,
337 : const Matrix<Float>& imagingweight,
338 : const Cube<Complex>& visData,
339 : const Matrix<Double>& uvw,
340 : const Cube<Int>& flagCube,
341 : const Vector<Double>& dphase,
342 : const Bool& dopsf,
343 : const Vector<Int>& /*gridShape*/)
344 : {
345 0 : vbs.vb_p = &vb;
346 0 : vbs.wbAWP_p=wbAWP_p;
347 0 : vbs.ftmType_p=ftmType_p;
348 0 : vbs.nWPlanes_p = nWPlanes_p;
349 : //cerr << "HPG setupvbstore " << endl;
350 :
351 0 : visResampler_p->setParams(uvScale,uvOffset,dphase);
352 0 : visResampler_p->setMaps(chanMap, polMap);
353 :
354 : //
355 : // Set up VBStore object to point to the relavent info. of the VB.
356 : //
357 0 : vbs.imRefFreq_p = imRefFreq_p;
358 0 : vbs.nRow_p = vb.nRows();
359 0 : vbs.beginRow_p = 0;
360 0 : vbs.endRow_p = vbs.nRow_p;
361 0 : vbs.spwID_p = vb.spectralWindows()(0);
362 0 : vbs.nDataPol_p = flagCube.shape()[0];
363 0 : vbs.nDataChan_p = flagCube.shape()[1];
364 :
365 0 : vbs.antenna1_p.reference(vb.antenna1());
366 0 : vbs.antenna2_p.reference(vb.antenna2());
367 : //vbs.paQuant_p = Quantity(getPA(vb),"rad");
368 :
369 0 : vbs.corrType_p.reference(vb.correlationTypes());
370 :
371 0 : vbs.uvw_p=uvw;
372 0 : vbs.imagingWeight_p.reference(imagingweight);
373 0 : vbs.visCube_p.reference(visData);
374 :
375 0 : vbs.freq_p.reference(vb.getFrequencies(0));
376 :
377 0 : vbs.rowFlag_p.reference(vb.flagRow());
378 0 : if(!usezero_p)
379 0 : for (Int rownr=0; rownr<vbs.nRow_p; rownr++)
380 0 : if(vb.antenna1()(rownr)==vb.antenna2()(rownr)) vbs.rowFlag_p(rownr)=true;
381 :
382 0 : vbs.flagCube_p.resize(flagCube.shape()); vbs.flagCube_p = false; vbs.flagCube_p(flagCube!=0) = true;
383 :
384 0 : vbs.conjBeams_p=conjBeams_p;
385 :
386 :
387 : // The following code is required only for GPU or multi-threaded
388 : //gridder. Currently does not work without the rest of the
389 : //GPU/multi-threaded infrastructure (though, I (SB) thought this
390 : //was designed to be benign for normal gridding).
391 : //
392 : //visResampler_p->initializeDataBuffers(vbs);
393 0 : }
394 0 : void AWProjectWBFTHPG::init(const vi::VisBuffer2& vb)
395 : {
396 0 : LogIO log_l(LogOrigin("AWProjectFT2", "init[R&D]"));
397 :
398 0 : nx = image->shape()(0);
399 0 : ny = image->shape()(1);
400 0 : npol = image->shape()(2);
401 0 : nchan = image->shape()(3);
402 :
403 :
404 0 : sumWeight.resize(npol, nchan);
405 0 : sumCFWeight.resize(npol, nchan);
406 :
407 0 : wConvSize=max(1, nWPlanes_p);
408 :
409 0 : CoordinateSystem cs=image->coordinates();
410 0 : uvScale.resize(3);
411 0 : uvScale=0.0;
412 0 : uvScale(0)=Float(nx)*cs.increment()(0);
413 0 : uvScale(1)=Float(ny)*cs.increment()(1);
414 0 : uvScale(2)=Float(wConvSize)*abs(cs.increment()(0));
415 :
416 0 : Int index= cs.findCoordinate(Coordinate::SPECTRAL);
417 0 : SpectralCoordinate spCS = cs.spectralCoordinate(index);
418 0 : imRefFreq_p = spCS.referenceValue()(0);
419 : double f1, f2;
420 0 : spCS.toWorld(f1, double(-0.5));
421 0 : spCS.toWorld(f2, double(nchan)-0.5);
422 0 : auto frange=std::make_pair(f1, f2);
423 0 : uvOffset.resize(3);
424 0 : uvOffset(0)=nx/2;
425 0 : uvOffset(1)=ny/2;
426 0 : uvOffset(2)=0;
427 :
428 0 : if(gridder) delete gridder;
429 0 : gridder=0;
430 0 : gridder = new ConvolveGridder<Double, Complex>(IPosition(2, nx, ny),
431 0 : uvScale, uvOffset,
432 0 : "SF");
433 0 : makingPSF = false;
434 :
435 :
436 :
437 : ///We'll always use oversampling of 4 for GPU gridder
438 : //convSampling=4;
439 : //TESTOO
440 :
441 0 : convSampling=10;
442 0 : if (min(nx, ny) > 600)
443 0 : convSampling = 4;
444 : ///////
445 :
446 : // we are not doing parallactiv angle here
447 0 : Double painc=2*M_PI;
448 :
449 0 : if(awConvs_p.use_count()==0){
450 0 : String observatory=(vb.subtableColumns().observation()).telescopeName()(0);
451 0 : awConvs_p=std::make_shared<AWConvFuncHolder>((*image).coordinates(), nx, ny,
452 0 : False, painc, observatory, convSampling);
453 0 : vi::VisibilityIterator2 *vi= const_cast<VisibilityIterator2 *>(vb.getVi());
454 :
455 :
456 0 : std::vector<Double> freqs;
457 0 : std::set<Int> fields;
458 0 : std::vector<Double> pAs={0.0};
459 : //int validspw=-1;
460 0 : Double maxW=0.0;
461 0 : for (vi->originChunks(); vi->moreChunks(); vi->nextChunk()) {
462 0 : for (vi->origin(); vi->more(); vi->next()) {
463 0 : std::vector<Double> chunkfreq;
464 0 : fields.insert(vb.fieldId()(0));
465 : //matchChannel(vb);
466 : //cerr << "MAX chanMap" << chanMap << endl;
467 : //if (max(chanMap) > -1)
468 : {
469 :
470 0 : SimplePBConvFunc::findUsefulChannels(chunkfreq, vb, frange);
471 : //cerr << "vbnchan " << vb.nChannels() << "chunkfreq "
472 : // << chunkfreq << endl;
473 0 : if (chunkfreq.size() > 0) {
474 : // validspw=vb.spectralWindows()(0);
475 : // cerr << "SPW " << vb.spectralWindows()(0) << " freqs " <<
476 : // Vector<Double>(chunkfreq) << endl;
477 0 : std::move(chunkfreq.begin(), chunkfreq.end(),
478 : std::back_inserter(freqs));
479 : double maxfreqused =
480 0 : *(std::max_element(chunkfreq.begin(), chunkfreq.end()));
481 0 : if (nWPlanes_p > 1) {
482 : // maxW=max(maxW,
483 : // max(abs(vb.uvw().row(2)*max(vb.getFrequencies(0))))/C::c);
484 0 : maxW = max(maxW,
485 0 : max(abs(vb.uvw().row(2) * maxfreqused)) / C::c);
486 : }
487 : }
488 : }
489 0 : }
490 : }
491 :
492 :
493 : //return vi to origin
494 0 : vi->originChunks(); vi->origin();
495 :
496 0 : std::sort(freqs.begin(), freqs.end());
497 0 : auto last = std::unique(freqs.begin(), freqs.end());
498 0 : freqs.erase(last, freqs.end());
499 : // tell holder it is a single field or not
500 0 : (*awConvs_p).setSingleField((fields.size() == 1));
501 :
502 0 : if (nWPlanes_p == 0)
503 0 : nWPlanes_p = 1;
504 0 : Vector<Double> wVals(nWPlanes_p,0);
505 0 : if(nWPlanes_p >1){
506 0 : Double st=maxW/(Double(nWPlanes_p-1)*Double(nWPlanes_p-1));
507 0 : for (int k=0; k <nWPlanes_p; ++k)
508 0 : wVals[k]=Double(k*k)*st;
509 : }
510 : //cerr << "XXXXXINIT wVals " << wVals << endl;
511 : //cerr << "XXXXXfreqs " << Vector<Double>(freqs) << endl;
512 0 : (*awConvs_p).addConvFunc(Vector<Double>(freqs), wVals, 0.0);
513 : /////TESTOO
514 : /*{
515 : Vector<Double>pixW(wVals.nelements());
516 : indgen(pixW);
517 : CoordinateSystem fiveAxis=image->coordinates();
518 : Vector<Int> stoks(4);
519 : stoks(0) = Stokes::RR;
520 : stoks(1) = Stokes::RL;
521 : stoks(2) = Stokes::LR;
522 : stoks(3) = Stokes::LL;
523 : StokesCoordinate stokesCoords(stoks);
524 : fiveAxis.replaceCoordinate(stokesCoords, 1);
525 :
526 : TabularCoordinate tab(pixW, wVals, "m", "W");
527 : fiveAxis.addCoordinate(tab);
528 : PagedImage<Complex> noo((awConvs_p->getConvFunc()).shape(), fiveAxis, "AWConvVals_"+String::toString(validspw));
529 : noo.put((awConvs_p->getConvFunc()));
530 :
531 : }*/
532 :
533 :
534 : ///////////////////////////
535 0 : visResampler_p->setConvFunc(awConvs_p);
536 0 : }
537 0 : }
538 0 : void AWProjectWBFTHPG::initializeToSky(ImageInterface<Complex>& iimage,
539 : Matrix<Float>& weight,
540 : const VisBuffer2& vb)
541 : {
542 0 : LogIO log_l(LogOrigin("AWProjectWBFT2","initializeToSky[R&D]"));
543 : // Timer tim;
544 : // tim.mark();
545 :
546 0 : AWProjectFT::initializeToSky(iimage,weight,vb);
547 : // cerr << "$$$$$ AWP initializesky " << tim.real()<< endl;
548 :
549 : // tim.mark();
550 0 : init(vb);
551 : //cerr << "$$$$$$$$$$$ init(vb) " << tim.real() << endl;
552 :
553 0 : }
554 0 : void AWProjectWBFTHPG::findConvFunction(const ImageInterface<Complex>& image,
555 : const VisBuffer2& vb)
556 : {
557 :
558 :
559 : // CoordinateSystem ftcoords;
560 : //We have to make sure awh is loaded
561 0 : if(awConvs_p.use_count()==0)
562 0 : throw(AipsError("Programmer's error:Convolution function has not been set"));
563 : /*
564 : if(avgPBReady_p){
565 : LatticeExprNode le( max( *avgPB_p ) );
566 : Float avgPB_max=le.getFloat();
567 :
568 : if(avgPB_max <= 0.0) avgPBReady_p = false;
569 : }
570 :
571 : if(!avgPBReady_p) makeSensitivityImage(vb,image,*avgPB_p);
572 :
573 :
574 :
575 : verifyShapes(avgPB_p->shape(), image.shape());
576 : */
577 :
578 :
579 :
580 0 : }
581 : }; // namespace refim
582 : }; // namespace casa
|
