diff --git a/CEP/BB/BBSKernel/include/BBSKernel/MNS/MeqStatUVW.h b/CEP/BB/BBSKernel/include/BBSKernel/MNS/MeqStatUVW.h index dbd18c787e4a4b040b404fb15147d5152205f27a..b33612bdd3c9cfd6eff3688d4a352b05f3c7639b 100644 --- a/CEP/BB/BBSKernel/include/BBSKernel/MNS/MeqStatUVW.h +++ b/CEP/BB/BBSKernel/include/BBSKernel/MNS/MeqStatUVW.h @@ -20,7 +20,7 @@ //# //# $Id$ -#if !defined(MNS_MEQSTATUVW_H) +#ifndef MNS_MEQSTATUVW_H #define MNS_MEQSTATUVW_H // \file @@ -31,7 +31,11 @@ #include <BBSKernel/MNS/MeqRequest.h> #include <BBSKernel/MNS/MeqStation.h> #include <Common/lofar_map.h> -#include <measures/Measures/MeasFrame.h> +#include <Common/lofar_vector.h> +#include <utility> +#include <measures/Measures/MDirection.h> +#include <measures/Measures/MPosition.h> + namespace LOFAR { @@ -44,8 +48,6 @@ namespace BBS //# Forward declarations class MeqStation; -class MeqPhaseRef; - class MeqStatUVW { @@ -53,7 +55,8 @@ public: // The default constructor. MeqStatUVW() {}; - MeqStatUVW (MeqStation*, const MeqPhaseRef*); +//# MeqStatUVW (MeqStation*, const MeqPhaseRef*); + MeqStatUVW(MeqStation* station, const pair<double, double> &phaseCenter, const vector<double> &arrayPosition); void calculate (const MeqRequest&); @@ -65,25 +68,17 @@ public: const MeqResult& getW (const MeqRequest& request) { if (request.getId() != itsLastReqId) calculate(request); return itsW; } - // Clear the map of UVW coordinates and reset the last request id. - void clear(); - // Set the UVW coordinate for the given time. void set (double time, double u, double v, double w); - const MeqStation *getStation() const - { - return itsStation; - } - private: struct MeqTime { MeqTime (double time) :itsTime(time) {} bool operator< (const MeqTime& other) const - { return itsTime < other.itsTime-0.000001; } - //operator== (const MeqTime& other) - // { return itsTime > other.itsTime-0.001 && itsTime < other.itsTime+0.001; } + { return itsTime < other.itsTime-0.000001; } + //# operator== (const MeqTime& other) + //# { return itsTime > other.itsTime-0.001 && itsTime < other.itsTime+0.001; } double itsTime; }; @@ -96,14 +91,17 @@ private: double itsW; }; - MeqStation* itsStation; - const MeqPhaseRef* itsPhaseRef; - MeqResult itsU; - MeqResult itsV; - MeqResult itsW; - MeqRequestId itsLastReqId; - casa::MeasFrame itsFrame; - map<MeqTime,MeqUVW> itsUVW; // association of time and UVW coordinates + MeqStation *itsStation; + casa::MDirection itsPhaseCenter; + casa::MPosition itsArrayPosition; + MeqResult itsU; + MeqResult itsV; + MeqResult itsW; + MeqRequestId itsLastReqId; + map<MeqTime,MeqUVW> itsUVW; // association of time and UVW coordinates + +//# const MeqPhaseRef* itsPhaseRef; +//# casa::MeasFrame itsFrame; }; // @} diff --git a/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.cc b/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.cc index 2e23a65faf5547371a13470bfb04f30be7fb9192..71f64a9b9ef7d5b2cbe6b44e612c5125c9e0257e 100644 --- a/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.cc +++ b/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.cc @@ -23,17 +23,17 @@ #include <lofar_config.h> #include <BBSKernel/MNS/MeqStatUVW.h> -#include <BBSKernel/MNS/MeqPhaseRef.h> #include <BBSKernel/MNS/MeqRequest.h> #include <BBSKernel/MNS/MeqExpr.h> #include <BBSKernel/MNS/MeqMatrixTmp.h> #include <Common/LofarLogger.h> +#include <Common/lofar_iomanip.h> + #include <measures/Measures/MBaseline.h> -#include <measures/Measures/MPosition.h> #include <measures/Measures/MEpoch.h> #include <measures/Measures/MeasConvert.h> +#include <measures/Measures/MeasFrame.h> #include <casa/Quanta/MVuvw.h> -#include <iomanip> using namespace casa; @@ -42,6 +42,10 @@ namespace LOFAR namespace BBS { +/* +//# Old constructor, uses MeqPhaseRef class which contains a hard-coded +//# reference to the WSRT. + MeqStatUVW::MeqStatUVW (MeqStation* station, const MeqPhaseRef* phaseRef) : itsStation (station), @@ -53,25 +57,32 @@ MeqStatUVW::MeqStatUVW (MeqStation* station, { itsFrame.set (itsPhaseRef->direction()); itsFrame.set (itsPhaseRef->earthPosition()); + LOG_TRACE_FLOW("Phase ref earth position: " << itsPhaseRef->earthPosition()); } +*/ -void MeqStatUVW::clear() + +MeqStatUVW::MeqStatUVW(MeqStation* station, const pair<double, double> &phaseCenter, const vector<double> &arrayPosition) +: itsStation (station), + itsPhaseCenter(MVDirection(phaseCenter.first, phaseCenter.second), MDirection::J2000), + itsArrayPosition(MVPosition(arrayPosition[0], arrayPosition[1], arrayPosition[2]), MPosition::ITRF), + itsU (0), + itsV (0), + itsW (0), + itsLastReqId (InitMeqRequestId) { - itsLastReqId = InitMeqRequestId; - itsUVW.clear(); - itsU.clear(); - itsV.clear(); - itsW.clear(); } + void MeqStatUVW::calculate (const MeqRequest& request) { - // Make sure the MeqResult/Matrix objects have the correct size. + //# Make sure the MeqResult/Matrix objects have the correct size. int ntime = request.ny(); double* uptr = itsU.setDoubleFormat (1, ntime); double* vptr = itsV.setDoubleFormat (1, ntime); double* wptr = itsW.setDoubleFormat (1, ntime); - // Use the UVW coordinates if already calculated for a time. + + //# Use the UVW coordinates if already calculated for a time. int ndone = 0; for (int i=0; i<ntime; ++i) { MeqTime meqtime(request.y(i)); @@ -83,54 +94,75 @@ void MeqStatUVW::calculate (const MeqRequest& request) ndone++; } } - // If all found we can exit. + //# If all found we can exit. if (ndone == ntime) { itsLastReqId = request.getId(); return; } - // Calculate the other UVW coordinates using the AIPS++ code. - // First form the time-independent stuff. + + //# Calculate the other UVW coordinates using the AIPS++ code. + //# First form the time-independent stuff. ASSERTSTR (itsStation, "UVW coordinates cannot be calculated"); + + //# Get station coordinates in ITRF coordinates MeqResult posx = itsStation->getPosX().getResult (request); MeqResult posy = itsStation->getPosY().getResult (request); MeqResult posz = itsStation->getPosZ().getResult (request); + LOG_TRACE_FLOW_STR ("posx" << posx.getValue()); LOG_TRACE_FLOW_STR ("posy" << posy.getValue()); LOG_TRACE_FLOW_STR ("posz" << posz.getValue()); - // Get position relative to center to keep values small. - const MVPosition& mvcpos = itsPhaseRef->earthPosition().getValue(); + + //# Get position relative to center to keep values small. +//# const MVPosition& mvcpos = itsPhaseRef->earthPosition().getValue(); + const MVPosition& mvcpos = itsArrayPosition.getValue(); MVPosition mvpos(posx.getValue().getDouble() - mvcpos(0), - posy.getValue().getDouble() - mvcpos(1), - posz.getValue().getDouble() - mvcpos(2)); + posy.getValue().getDouble() - mvcpos(1), + posz.getValue().getDouble() - mvcpos(2)); + MVBaseline mvbl(mvpos); MBaseline mbl(mvbl, MBaseline::ITRF); LOG_TRACE_FLOW_STR ("mbl " << mbl); + Quantum<double> qepoch(0, "s"); +//# ?? qepoch.setValue(time(0)) ?? MEpoch mepoch(qepoch, MEpoch::UTC); - itsFrame.set (mepoch); - mbl.getRefPtr()->set(itsFrame); // attach frame +//# itsFrame.set (mepoch); + MeasFrame frame(itsArrayPosition); + frame.set(itsPhaseCenter); + frame.set(mepoch); + +//# mbl.getRefPtr()->set(itsFrame); // attach frame + mbl.getRefPtr()->set(frame); // attach frame + MBaseline::Convert mcvt(mbl, MBaseline::J2000); + for (int i=0; i<ntime; ++i) { double time = request.y(i); map<MeqTime,MeqUVW>::iterator pos = itsUVW.find (MeqTime(time)); if (pos == itsUVW.end()) { - qepoch.setValue (time); + qepoch.setValue(time); + mepoch.set(qepoch); + frame.set(mepoch); + LOG_TRACE_FLOW_STR ("frame " << mbl.getRefPtr()->getFrame()); + const MVBaseline& bas2000 = mcvt().getValue(); LOG_TRACE_FLOW_STR (bas2000); - LOG_TRACE_FLOW_STR (itsPhaseRef->direction().getValue()); - MVuvw uvw2000 (bas2000, itsPhaseRef->direction().getValue()); +//# LOG_TRACE_FLOW_STR (itsPhaseRef->direction().getValue()); + LOG_TRACE_FLOW_STR (itsPhaseCenter.getValue()); +//# MVuvw uvw2000 (bas2000, itsPhaseRef->direction().getValue()); + MVuvw uvw2000 (bas2000, itsPhaseCenter.getValue()); + const Vector<double>& xyz = uvw2000.getValue(); + LOG_TRACE_FLOW_STR (xyz(0) << ' ' << xyz(1) << ' ' << xyz(2)); *uptr++ = xyz(0); *vptr++ = xyz(1); *wptr++ = xyz(2); // Save the UVW coordinates in the map. itsUVW[MeqTime(time)] = MeqUVW(xyz(0), xyz(1), xyz(2)); -/// if (itsStation->getName() == "SR9" || itsStation->getName() == "SR12") { -/// cout << "UVW="<<itsStation->getName()<<' ' << std::setprecision(12)<<xyz(0) << ' '<< std::setprecision(12)<<xyz(1)<<' '<< std::setprecision(12)<<xyz(2) << ' '<<std::setprecision(12)<<time<<' '<<mvcpos(0)<<' '<<mvcpos(1)<<' '<<mvcpos(2)<<' ' <<posz.getValue()<<' '<<posy.getValue()<<' '<<posz.getValue()<<endl; -/// } - } + } } LOG_TRACE_FLOW_STR ('U' << itsU.getValue()); LOG_TRACE_FLOW_STR ('V' << itsV.getValue()); @@ -138,6 +170,7 @@ void MeqStatUVW::calculate (const MeqRequest& request) itsLastReqId = request.getId(); } + void MeqStatUVW::set (double time, double u, double v, double w) { itsUVW[MeqTime(time)] = MeqUVW(u,v,w); diff --git a/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.h b/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.h index dbd18c787e4a4b040b404fb15147d5152205f27a..b33612bdd3c9cfd6eff3688d4a352b05f3c7639b 100644 --- a/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.h +++ b/CEP/BB/BBSKernel/src/MNS/MeqStatUVW.h @@ -20,7 +20,7 @@ //# //# $Id$ -#if !defined(MNS_MEQSTATUVW_H) +#ifndef MNS_MEQSTATUVW_H #define MNS_MEQSTATUVW_H // \file @@ -31,7 +31,11 @@ #include <BBSKernel/MNS/MeqRequest.h> #include <BBSKernel/MNS/MeqStation.h> #include <Common/lofar_map.h> -#include <measures/Measures/MeasFrame.h> +#include <Common/lofar_vector.h> +#include <utility> +#include <measures/Measures/MDirection.h> +#include <measures/Measures/MPosition.h> + namespace LOFAR { @@ -44,8 +48,6 @@ namespace BBS //# Forward declarations class MeqStation; -class MeqPhaseRef; - class MeqStatUVW { @@ -53,7 +55,8 @@ public: // The default constructor. MeqStatUVW() {}; - MeqStatUVW (MeqStation*, const MeqPhaseRef*); +//# MeqStatUVW (MeqStation*, const MeqPhaseRef*); + MeqStatUVW(MeqStation* station, const pair<double, double> &phaseCenter, const vector<double> &arrayPosition); void calculate (const MeqRequest&); @@ -65,25 +68,17 @@ public: const MeqResult& getW (const MeqRequest& request) { if (request.getId() != itsLastReqId) calculate(request); return itsW; } - // Clear the map of UVW coordinates and reset the last request id. - void clear(); - // Set the UVW coordinate for the given time. void set (double time, double u, double v, double w); - const MeqStation *getStation() const - { - return itsStation; - } - private: struct MeqTime { MeqTime (double time) :itsTime(time) {} bool operator< (const MeqTime& other) const - { return itsTime < other.itsTime-0.000001; } - //operator== (const MeqTime& other) - // { return itsTime > other.itsTime-0.001 && itsTime < other.itsTime+0.001; } + { return itsTime < other.itsTime-0.000001; } + //# operator== (const MeqTime& other) + //# { return itsTime > other.itsTime-0.001 && itsTime < other.itsTime+0.001; } double itsTime; }; @@ -96,14 +91,17 @@ private: double itsW; }; - MeqStation* itsStation; - const MeqPhaseRef* itsPhaseRef; - MeqResult itsU; - MeqResult itsV; - MeqResult itsW; - MeqRequestId itsLastReqId; - casa::MeasFrame itsFrame; - map<MeqTime,MeqUVW> itsUVW; // association of time and UVW coordinates + MeqStation *itsStation; + casa::MDirection itsPhaseCenter; + casa::MPosition itsArrayPosition; + MeqResult itsU; + MeqResult itsV; + MeqResult itsW; + MeqRequestId itsLastReqId; + map<MeqTime,MeqUVW> itsUVW; // association of time and UVW coordinates + +//# const MeqPhaseRef* itsPhaseRef; +//# casa::MeasFrame itsFrame; }; // @}