00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029 #ifndef SYNTHESIS_PBWPROJECTFT_H
00030 #define SYNTHESIS_PBWPROJECTFT_H
00031
00032 #include <synthesis/TransformMachines/FTMachine.h>
00033 #include <casa/Arrays/Matrix.h>
00034 #include <scimath/Mathematics/FFTServer.h>
00035 #include <msvis/MSVis/VisBuffer.h>
00036 #include <images/Images/ImageInterface.h>
00037 #include <casa/Containers/Block.h>
00038 #include <casa/Arrays/Array.h>
00039 #include <casa/Arrays/Vector.h>
00040 #include <casa/Arrays/Matrix.h>
00041 #include <scimath/Mathematics/ConvolveGridder.h>
00042 #include <lattices/Lattices/LatticeCache.h>
00043 #include <lattices/Lattices/ArrayLattice.h>
00044 #include <ms/MeasurementSets/MSColumns.h>
00045 #include <measures/Measures/Measure.h>
00046 #include <measures/Measures/MDirection.h>
00047 #include <measures/Measures/MPosition.h>
00048 #include <coordinates/Coordinates/DirectionCoordinate.h>
00049 #include <synthesis/TransformMachines/VPSkyJones.h>
00050 #include <synthesis/TransformMachines/VLACalcIlluminationConvFunc.h>
00051 #include <synthesis/TransformMachines/VLAIlluminationConvFunc.h>
00052 #include <synthesis/TransformMachines/Utils.h>
00053
00054
00055 #include <synthesis/MeasurementComponents/SolvableVisCal.h>
00056 #include <synthesis/MeasurementComponents/ConvFuncDiskCache.h>
00057
00058 namespace casa {
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136 class SolvableVisJones;
00137 class nPBWProjectFT : public FTMachine {
00138 public:
00139
00140
00141
00142
00143
00144
00145 nPBWProjectFT(Int nFacets, Long cachesize, String& cfCacheDirName,
00146 Bool applyPointingOffset=True,
00147 Bool doPBCorr=True,
00148 Int tilesize=16,
00149 Float paSteps=5.0, Float pbLimit=5e-2,
00150 Bool usezero=False);
00151
00152
00153
00154 nPBWProjectFT(const RecordInterface& stateRec);
00155
00156
00157 nPBWProjectFT(const nPBWProjectFT &other);
00158
00159
00160 nPBWProjectFT &operator=(const nPBWProjectFT &other);
00161
00162 ~nPBWProjectFT();
00163
00164
00165 void setEPJones(SolvableVisJones* ep_j) {epJ = ep_j;}
00166
00167 void setDOPBCorrection(Bool doit=True) {doPBCorrection=doit;};
00168
00169
00170
00171 virtual void initializeToVis(ImageInterface<Complex>& image,
00172 const VisBuffer& vb);
00173
00174
00175 virtual void initializeToVis(ImageInterface<Complex>& image,
00176 const VisBuffer& vb, Array<Complex>& griddedVis,
00177 Vector<Double>& uvscale);
00178
00179
00180
00181 virtual void finalizeToVis();
00182
00183
00184 virtual void initializeToSky(ImageInterface<Complex>& image, Matrix<Float>& weight,
00185 const VisBuffer& vb);
00186
00187
00188
00189
00190 virtual void finalizeToSky();
00191
00192 virtual void initVisBuffer(VisBuffer& vb, Type whichVBColumn);
00193 void initVisBuffer(VisBuffer& vb, Type whichVBColumn, Int row);
00194
00195
00196 void get(VisBuffer& vb, Int row=-1);
00197
00198
00199
00200
00201 void get(VisBuffer& vb, Cube<Complex>& degrid,
00202 Array<Complex>& griddedVis, Vector<Double>& scale,
00203 Int row=-1);
00204
00205 void get(VisBuffer& vb, Cube<Float>& pointingOffsets, Int row=-1,
00206 Type whichVBColumn=FTMachine::MODEL,Int Conj=0)
00207 {
00208 get(vb,vb,vb,pointingOffsets,row,whichVBColumn,whichVBColumn,Conj,0);
00209 }
00210
00211 void get(VisBuffer& vb, VisBuffer& gradAzVB,VisBuffer& gradElVB,
00212 Cube<Float>& pointingOffsets,Int row=-1,
00213 Type whichVBColumn=FTMachine::MODEL,
00214 Type whichGradVBColumn=FTMachine::MODEL,
00215 Int Conj=0, Int doGrad=1) ;
00216 void nget(VisBuffer& vb,
00217
00218 Array<Float>& l_off, Array<Float>& m_off,
00219 Cube<Complex>& Mout,
00220 Cube<Complex>& dMout1,
00221 Cube<Complex>& dMout2,
00222 Int Conj=0, Int doGrad=1);
00223
00224
00225
00226 void get(VisBuffer& vb, Cube<Complex>& degrid,
00227 Array<Complex>& griddedVis, Vector<Double>& scale,
00228 Cube<Float>& pointingOffsets,Int row=-1);
00229
00230
00231
00232
00233
00234
00235
00236 void put(const VisBuffer&,
00237 TempImage<Complex>&, Vector<Double>&, int,
00238 UVWMachine*, Bool)
00239 {
00240
00241 }
00242 void put(const VisBuffer& vb, Int row=-1, Bool dopsf=False,
00243 FTMachine::Type type=FTMachine::OBSERVED);
00244
00245
00246 void makeImage(FTMachine::Type type,
00247 VisSet& vs,
00248 ImageInterface<Complex>& image,
00249 Matrix<Float>& weight);
00250
00251
00252
00253 virtual ImageInterface<Complex>& getImage(Matrix<Float>&, Bool normalize=True);
00254 virtual void normalizeImage(Lattice<Complex>& ,
00255 const Matrix<Double>& ,
00256 Lattice<Float>& ,
00257 Bool )
00258 {throw(AipsError("nPBWProjectFT::normalizeImage() called"));}
00259
00260
00261 void getWeightImage(ImageInterface<Float>&, Matrix<Float>&);
00262
00263
00264 Bool toRecord(String& error, RecordInterface& outRec,
00265 Bool withImage=False, const String diskimage="");
00266 Bool fromRecord(String& error, const RecordInterface& inRec);
00267
00268
00269 Bool isFourier() {return True;}
00270
00271
00272 Bool changed(const VisBuffer& ) {return False;}
00273
00274 virtual Int findPointingOffsets(const VisBuffer&, Array<Float>&, Array<Float>&,
00275 Bool Evaluate=True);
00276 virtual Int findPointingOffsets(const VisBuffer&, Cube<Float>&,
00277 Array<Float>&, Array<Float>&,
00278 Bool Evaluate=True);
00279 virtual Double getVBPA(const VisBuffer& vb)
00280 {
00281
00282
00283 return getPA(vb);
00284 };
00285 MDirection::Convert makeCoordinateMachine(const VisBuffer&,
00286 const MDirection::Types&,
00287 const MDirection::Types&,
00288 MEpoch& last);
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309 virtual void makeSensitivityImage(Lattice<Complex>& wtImage,
00310 ImageInterface<Float>& sensitivityImage,
00311 const Matrix<Float>& sumWt=Matrix<Float>(),
00312 const Bool& doFFTNorm=True);
00313 virtual Bool makeAveragePB0(const VisBuffer& vb,
00314 const ImageInterface<Complex>& image,
00315 Int& polInUse,
00316 TempImage<Float>& avgPB);
00317
00318
00319
00320
00321
00322
00323 void makeConjPolMap(const VisBuffer& vb, const Vector<Int> cfPolMap, Vector<Int>& conjPolMap);
00324
00325 void makeCFPolMap(const VisBuffer& vb, Vector<Int>& polM);
00326
00327 void reset() {paChangeDetector.reset();}
00328
00329 void setPAIncrement(const Quantity &paIncrement);
00330
00331 Vector<Int>& getPolMap() {return polMap;};
00332 virtual String name() const { return "PBWProjectFT";};
00333 virtual Bool verifyAvgPB(ImageInterface<Float>& pb, ImageInterface<Float>& sky)
00334 {return verifyShapes(pb.shape(),sky.shape());}
00335 virtual Bool verifyAvgPB(ImageInterface<Float>& pb, ImageInterface<Complex>& sky)
00336 {return verifyShapes(pb.shape(),sky.shape());}
00337 virtual Bool verifyShapes(IPosition shape0, IPosition shape1);
00338 Bool findSupport(Array<Complex>& func, Float& threshold, Int& origin, Int& R);
00339 void makeAntiAliasingOp(Vector<Complex>& val, const Int len);
00340 void makeAntiAliasingCorrection(Vector<Complex>& correction,
00341 const Vector<Complex>& op,
00342 const Int nx);
00343 void applyAntiAliasingOp(ImageInterface<Complex>& cf,
00344 Vector<IPosition>& offset,
00345 Int op=0,
00346 Bool Square=False);
00347 void applyAntiAliasingOp(ImageInterface<Float>& cf,
00348 Vector<IPosition>& offset,
00349 Int op=0,
00350 Bool Square=False);
00351 void correctAntiAliasing(Lattice<Complex>& cf);
00352 virtual void setMiscInfo(const Int qualifier){(void)qualifier;};
00353 virtual void ComputeResiduals(VisBuffer&, Bool ) {};
00354
00355 protected:
00356
00357
00358 Float padding_p;
00359
00360 Int nint(Double val){return Int(floor(val+0.5));};
00361
00362
00363 Int makePBPolnCoords(
00364 CoordinateSystem& coord, const VisBuffer& vb);
00365
00366 Int locateConvFunction(Int Nw, Int polInUse, const VisBuffer& vb, Float &pa);
00367 void cacheConvFunction(Int which, Array<Complex>& cf, CoordinateSystem& coord);
00368
00369 void findConvFunction(const ImageInterface<Complex>& image,
00370 const VisBuffer& vb);
00371 void makeConvFunction(const ImageInterface<Complex>& image,
00372 const VisBuffer& vb, Float pa);
00373
00374 Int nWPlanes_p;
00375
00376
00377 Array<Complex>* getDataPointer(const IPosition&, Bool);
00378
00379 void ok();
00380
00381 void init();
00382
00383 Int getVisParams(const VisBuffer& vb);
00384
00385
00386 Bool recordOnGrid(const VisBuffer& vb, Int rownr) const;
00387
00388
00389 LatticeCache<Complex> * imageCache;
00390
00391
00392 Long cachesize;
00393 Int tilesize;
00394
00395
00396 ConvolveGridder<Double, Complex>* gridder;
00397
00398
00399 Bool isTiled;
00400
00401
00402
00403 CountedPtr<Lattice<Complex> > arrayLattice;
00404
00405
00406
00407
00408 CountedPtr<Lattice<Complex> > lattice;
00409
00410 Float maxAbsData;
00411
00412
00413 IPosition centerLoc, offsetLoc;
00414
00415
00416 Vector<Double> uvScale, uvOffset;
00417
00418
00419 Array<Complex> griddedData;
00420
00421
00422 MSPointingColumns* mspc;
00423
00424
00425 MSAntennaColumns* msac;
00426
00427 DirectionCoordinate directionCoord;
00428
00429 MDirection::Convert* pointingToImage;
00430
00431 Vector<Double> xyPos;
00432
00433 MDirection worldPosMeas;
00434
00435 Int priorCacheSize;
00436
00437
00438 Bool usezero_p;
00439
00440 Array<Complex> convFunc;
00441 Array<Complex> convWeights;
00442 CoordinateSystem convFuncCS_p;
00443 Int convSize;
00444
00445
00446
00447
00448
00449 Int convSampling;
00450 Cube<Int> convSupport, convWtSupport;
00451
00452
00453
00454
00455
00456 PtrBlock < Array<Complex> *> convFuncCache, convWeightsCache;
00457
00458
00459
00460
00461
00462 Int PAIndex;
00463
00464
00465
00466 Bool convFuncCacheReady;
00467
00468 Int wConvSize;
00469
00470 Int lastIndex_p;
00471
00472 Int getIndex(const ROMSPointingColumns& mspc, const Double& time,
00473 const Double& interval);
00474
00475 Bool getXYPos(const VisBuffer& vb, Int row);
00476
00477
00478
00479
00480
00481 TempImage<Float> avgPB;
00482
00483
00484
00485
00486 Int polInUse, bandID_p;
00487 Int maxConvSupport;
00488
00489
00490
00491
00492 Float pbLimit_p;
00493
00494
00495 SolvableVisJones *epJ;
00496 Double HPBW, Diameter_p, sigma;
00497 Int Nant_p;
00498 Int doPointing;
00499 Bool doPBCorrection;
00500 Bool makingPSF;
00501
00502 Unit Second, Radian, Day;
00503 Array<Float> l_offsets,m_offsets;
00504 Int noOfPASteps;
00505 Vector<Float> pbPeaks;
00506 Bool pbNormalized,resetPBs,rotateAperture_p;
00507 Vector<Float> paList;
00508 ConvFuncDiskCache cfCache;
00509 Double currentCFPA;
00510 ParAngleChangeDetector paChangeDetector;
00511 Vector<Int> cfStokes;
00512 Vector<Complex> Area;
00513 Double cfRefFreq_p;
00514 Bool avgPBSaved;
00515 Bool avgPBReady;
00516 Vector<Complex> antiAliasingOp,antiAliasingCorrection;
00517 Float lastPAUsedForWtImg;
00518
00519
00520
00521
00522 virtual void normalizeAvgPB();
00523 virtual void runFortranGet(Matrix<Double>& uvw,Vector<Double>& dphase,
00524 Cube<Complex>& visdata,
00525 IPosition& s,
00526
00527
00528
00529 Int& Conj,
00530 Cube<Int>& flags,Vector<Int>& rowFlags,
00531 Int& rownr,Vector<Double>& actualOffset,
00532 Array<Complex>* dataPtr,
00533 Int& aNx, Int& aNy, Int& npol, Int& nchan,
00534 VisBuffer& vb,Int& Nant_p, Int& scanNo,
00535 Double& sigma,
00536 Array<Float>& raoffsets,
00537 Array<Float>& decoffsets,
00538 Double area,
00539 Int& doGrad,Int paIndex);
00540 virtual void runFortranPut(Matrix<Double>& uvw,Vector<Double>& dphase,
00541 const Complex& visdata_p,
00542 IPosition& s,
00543
00544
00545
00546 Int& Conj,
00547 Cube<Int>& flags,Vector<Int>& rowFlags,
00548 const Matrix<Float>& weight,
00549 Int& rownr,Vector<Double>& actualOffset,
00550 Array<Complex>& dataPtr,
00551 Int& aNx, Int& aNy, Int& npol, Int& nchan,
00552 const VisBuffer& vb,Int& Nant_p, Int& scanNo,
00553 Double& sigma,
00554 Array<Float>& raoffsets,
00555 Array<Float>& decoffsets,
00556 Matrix<Double>& sumWeight,
00557 Double& area,
00558 Int& doGrad,
00559 Int& doPSF,Int paIndex);
00560 void runFortranGetGrad(Matrix<Double>& uvw,Vector<Double>& dphase,
00561 Cube<Complex>& visdata,
00562 IPosition& s,
00563 Cube<Complex>& gradVisAzData,
00564 Cube<Complex>& gradVisElData,
00565
00566 Int& Conj,
00567 Cube<Int>& flags,Vector<Int>& rowFlags,
00568 Int& rownr,Vector<Double>& actualOffset,
00569 Array<Complex>* dataPtr,
00570 Int& aNx, Int& aNy, Int& npol, Int& nchan,
00571 VisBuffer& vb,Int& Nant_p, Int& scanNo,
00572 Double& sigma,
00573 Array<Float>& l_off,
00574 Array<Float>& m_off,
00575 Double area,
00576 Int& doGrad,
00577 Int paIndex);
00578 };
00579
00580
00581 }
00582
00583 #endif