Line data Source code
1 : // -*- C++ -*-
2 : //# PhaseGrad.cc: Implementation of the PhaseGrad class
3 : //# Copyright (C) 1997,1998,1999,2000,2001,2002,2003
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 :
30 : #include <synthesis/TransformMachines2/PhaseGrad.h>
31 : #include <synthesis/TransformMachines/SynthesisMath.h>
32 : #include <casacore/casa/Logging/LogIO.h>
33 : #include <casacore/casa/Logging/LogSink.h>
34 : #include <casacore/casa/Logging/LogOrigin.h>
35 :
36 : using namespace casacore;
37 : namespace casa{
38 : using namespace refim;
39 : //
40 : //----------------------------------------------------------------------
41 : //
42 0 : PhaseGrad& PhaseGrad::operator=(const PhaseGrad& other)
43 : {
44 0 : if(this!=&other)
45 : {
46 0 : field_phaseGrad_p = other.field_phaseGrad_p;
47 0 : antenna_phaseGrad_p = other.antenna_phaseGrad_p;
48 0 : cached_FieldOffset_p = other.cached_FieldOffset_p;
49 0 : needCFPhaseGrad_p=other.needCFPhaseGrad_p;
50 : }
51 0 : return *this;
52 : }
53 : //
54 : //----------------------------------------------------------------------
55 : //
56 3684117 : bool PhaseGrad::needsNewPhaseGrad(const VisBuffer2& vb,
57 : const int& /*row*/)
58 : {
59 3684117 : unsigned int nRow=vb.nRows();
60 3684117 : if (cached_FieldOffset_p.nelements() < nRow) cached_FieldOffset_p.resize(nRow,true);
61 :
62 7368058 : return ((field_phaseGrad_p.shape()[0] < maxCFShape_p[0]) ||
63 7368058 : (field_phaseGrad_p.shape()[1] < maxCFShape_p[1]));
64 : }
65 : //
66 : //----------------------------------------------------------------------
67 : //
68 : // bool PhaseGrad::ComputeFieldPointingGrad(const Vector<double>& pointingOffset,
69 : // const CountedPtr<CFBuffer>& cfb,
70 : // const Vector<int>&cfShape,
71 : // const Vector<int>& convOrigin,
72 : // const double& /*cfRefFreq*/,
73 : // const double& /*imRefFreq*/,
74 : // const int& spwID, const int& fieldId)
75 : // bool PhaseGrad::ComputeFieldPointingGrad(const CountedPtr<PointingOffsets>& pointingOffsets_p,
76 : // const CountedPtr<CFBuffer>& cfb,
77 : // const VisBuffer2& vb,
78 : // const int& row,
79 : // const pair<int,int> antGrp)
80 3684117 : bool PhaseGrad::ComputeFieldPointingGrad(const CountedPtr<PointingOffsets>& pointingOffsets_p,
81 : const VisBuffer2& vb,
82 : const int& row,
83 : const pair<int,int> antGrp)
84 :
85 : {
86 3684117 : if (!needCFPhaseGrad_p)
87 : {
88 0 : needsNewPOPG_p = false;
89 0 : return true;
90 : }
91 :
92 : //
93 : // Re-find the max. CF size if the CFB changed.
94 : //
95 : // CFBuffer *thisCFB = cfb.get();
96 : // if (thisCFB != cachedCFBPtr_p)
97 : // {
98 : // maxCFShape_p[0] = maxCFShape_p[1] = cfb->getMaxCFSize();
99 : // {
100 : // // LogIO log_l(LogOrigin("PhaseGrad","computeFieldPointingGrad"));
101 : // //cerr << "CFB changed: "<< thisCFB << " " << cachedCFBPtr_p << " " << vb.spectralWindows()(0) << " " << vb.fieldId()(0) << " " << maxCFShape_p << endl;
102 : // }
103 : // cachedCFBPtr_p = thisCFB;
104 : // }
105 : //
106 : // If the pointing or the max. CF size changed, recompute the phase gradient.
107 : //
108 11052351 : LogIO log_l(LogOrigin("PhaseGrad","computeAntennaPointingGrad"));
109 3684117 : Bool needsCFShape_l = (needsNewPhaseGrad(vb, row));
110 3684117 : if (needsCFShape_l || needsNewFieldPG_p || needsNewPOPG_p )
111 : {
112 456 : Vector<Vector<double> > pointingOffset;
113 228 : if(pointingOffsets_p->getDoPointing())
114 : {
115 111 : pointingOffset.assign(pointingOffsets_p->pullAntGridPointingOffsets());
116 : // cerr << "Pointing Offsets in PG: "<<pointingOffset<<endl;
117 : // cerr << "NeedsCFShape_l " << needsCFShape_l <<" NeedsNewFieldPG_p " << needsNewFieldPG_p << " NeedsNewPOPG_p " << needsNewPOPG_p <<" row " <<row << endl;
118 : }
119 : else
120 117 : pointingOffset = pointingOffsets_p->pullPointingOffsets();
121 228 : if(needsNewFieldPG_p)
122 : {
123 135 : log_l << "Computing Phase Grad for Field : " << vb.fieldId()(0) << " for spw :"
124 135 : << vb.spectralWindows()(0)
125 405 : << LogIO::POST;
126 : }
127 : // else if(needsCFShape_l)
128 : // {
129 : // log_l << "Computing Phase Grad for change of max CF size : " << maxCFShape_p[0]
130 : // << LogIO::POST;
131 : // }
132 : // else if(needsNewPOPG_p)
133 : // {
134 :
135 : // // cerr << "NeedsNewPOPG_p for row "<< row << " " << vb.fieldId()(0)<< " " << vb.spectralWindows()(0) << endl;
136 : // // log_l << "Computing Phase Grad for Antenna Pointing Offset for row : " << row << " " << pointingOffset[row][0] << " " << pointingOffset[row][1] << " spw :"
137 : // // << vb.spectralWindows()(0) << " and field " << vb.fieldId()(0)
138 : // // << LogIO::POST;
139 : // }
140 :
141 228 : int nx=maxCFShape_p(0), ny=maxCFShape_p(1);
142 : double grad;
143 228 : DComplex phx,phy;
144 912 : Vector<int> convOrigin = (maxCFShape_p[0]%2==0) ? maxCFShape_p/2 : maxCFShape_p/2+1 ;
145 :
146 228 : field_phaseGrad_p.resize(nx,ny);
147 : // cached_FieldOffset_p[row] = pointingOffset[row];
148 456 : vector<double> ant1PO_l, ant2PO_l, blPO_l;
149 228 : Vector<double> tmpblPO_l;
150 228 : ant1PO_l.resize(2,0);
151 228 : ant2PO_l.resize(2,0);
152 228 : blPO_l.resize(2,0);
153 :
154 228 : if(pointingOffsets_p->getDoPointing())
155 : {
156 111 : ant1PO_l[0] = (pointingOffset[0][antGrp.first]);
157 111 : ant1PO_l[1] = (pointingOffset[1][antGrp.first]);
158 111 : ant2PO_l[0] = (pointingOffset[0][antGrp.second]);
159 111 : ant2PO_l[1] = (pointingOffset[1][antGrp.second]);
160 :
161 : // cerr<<"Ant1: "<<antGrp.first << " "<< ant1PO_l[0] << " " << ant1PO_l[1] << endl;
162 : // cerr<<"Ant2: "<<antGrp.second << " "<< ant2PO_l[0] << " " << ant2PO_l[1] << endl;
163 111 : blPO_l[0] = (ant1PO_l[0] + ant2PO_l[0])/2;
164 111 : blPO_l[1] = (ant1PO_l[1] + ant2PO_l[1])/2;
165 111 : Vector<double> const blPO_l_casavec(blPO_l);
166 111 : tmpblPO_l = pointingOffsets_p->gradPerPixel(blPO_l_casavec);
167 111 : cached_FieldOffset_p[row](0) = tmpblPO_l[0];
168 111 : cached_FieldOffset_p[row](1) = tmpblPO_l[1];
169 : }
170 : else
171 : {
172 117 : cached_FieldOffset_p[row](0) = pointingOffset[row][0];
173 117 : cached_FieldOffset_p[row](1) = pointingOffset[row][1];
174 :
175 : // cached_FieldOffset_p[row](0) = pointingOffsets_p->gradPerPixel((ant1PO_l[0] + ant2PO_l[0])/2);
176 : // cached_FieldOffset_p[row](1) = pointingOffsets_p->gradPerPixel((ant1PO_l[1] + ant2PO_l[1])/2);
177 : }
178 :
179 : {
180 : //cerr << "Cached Field Offset is : " << cached_FieldOffset_p[row](0) << " " << cached_FieldOffset_p[row](1)<< " for row " << row << " field id " << vb.fieldId()(0) << " " << needCFPhaseGrad_p<< endl;
181 78404 : for(int ix=0;ix<nx;ix++)
182 : {
183 78176 : grad = (ix-convOrigin[0])*cached_FieldOffset_p[row](0);
184 : double sx,cx;
185 78176 : SINCOS(grad,sx,cx);
186 78176 : phx = Complex(cx,sx);
187 32681064 : for(int iy=0;iy<ny;iy++)
188 : {
189 32602888 : grad = (iy-convOrigin[1])*cached_FieldOffset_p[row](1);
190 : Double sy,cy;
191 32602888 : SINCOS(grad,sy,cy);
192 32602888 : phy = Complex(cy,sy);
193 32602888 : field_phaseGrad_p(ix,iy)=phx*phy;
194 : }
195 : }
196 : }
197 228 : needsNewPOPG_p = false;
198 228 : return true; // New phase gradient was computed
199 : }
200 3683889 : return false;
201 : }
202 : }
|
