From fb20401162488cdc56d8d26da61a802fc6b3af4d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Andr=C3=A9=20Offringa?= <offringa@astron.nl>
Date: Tue, 28 Sep 2010 19:29:51 +0000
Subject: [PATCH] Bug 1491: some test things
---
.../include/AOFlagger/imaging/model.h | 47 +++-
.../include/AOFlagger/imaging/observatorium.h | 264 +++++++++---------
.../include/AOFlagger/rfi/rfistatistics.h | 1 -
CEP/DP3/AOFlagger/src/gui/mswindow.cpp | 38 ++-
CEP/DP3/AOFlagger/src/imaging/model.cpp | 38 ++-
5 files changed, 230 insertions(+), 158 deletions(-)
diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/imaging/model.h b/CEP/DP3/AOFlagger/include/AOFlagger/imaging/model.h
index 0906e9904d9..b1d4535f1e5 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/imaging/model.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/imaging/model.h
@@ -26,6 +26,35 @@
#include <AOFlagger/msio/image2d.h>
#include <AOFlagger/msio/types.h>
+#include <AOFlagger/imaging/uvimager.h>
+
+template<typename T>
+struct OutputReceiver
+{
+};
+template<>
+struct OutputReceiver<class UVImager>
+{
+ UVImager *_imager;
+ void SetUVValue(size_t, double u, double v, double r, double i, double w)
+ {
+ _imager->SetUVValue(u, v, r, i, w);
+ _imager->SetUVValue(-u, -v, r, -i, w);
+ }
+ void SetY(size_t) { }
+};
+template<>
+struct OutputReceiver<class TimeFrequencyData>
+{
+ Image2DPtr _real, _imaginary;
+ size_t _y;
+ void SetUVValue(size_t x, double, double, double r, double i, double)
+ {
+ _real->SetValue(x, _y, r);
+ _imaginary->SetValue(x, _y, i);
+ }
+ void SetY(size_t y) { _y = y; }
+};
/**
@author A.R. Offringa <offringa@astro.rug.nl>
*/
@@ -33,6 +62,7 @@ class Model {
struct PointSource {
long double dec, ra, fluxIntensity, sqrtFluxIntensity;
};
+
public:
Model();
~Model();
@@ -48,8 +78,21 @@ class Model {
}
void SimulateAntenna(num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t frequency, num_t earthLattitude, num_t &r, num_t &i);
void SimulateUncoherentAntenna(num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t frequency, num_t earthLattitude, num_t &r, num_t &i, size_t index);
- void SimulateCorrelation(class UVImager &imager, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t frequency, double totalTime, double integrationTime);
- void SimulateObservation(class UVImager &imager, class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency);
+
+ template<typename T>
+ void SimulateCorrelation(struct OutputReceiver<T> &receiver, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t frequency, double totalTime, double integrationTime);
+
+ void SimulateObservation(class UVImager &imager, class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency)
+ {
+ OutputReceiver<UVImager> imagerOutputter;
+ imagerOutputter._imager = &imager;
+ SimulateObservation(imagerOutputter, observatorium, delayDirectionDEC, delayDirectionRA, frequency);
+ }
+
+ std::pair<TimeFrequencyData, TimeFrequencyMetaDataPtr> SimulateObservation(class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency);
+
+ template<typename T>
+ void SimulateObservation(struct OutputReceiver<T> &receiver, class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency);
static void GetUVPosition(num_t &u, num_t &v, num_t earthLattitudeAngle, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t waveLength);
static num_t GetWPosition(num_t delayDirectionDec, num_t delayDirectionRA, num_t frequency, num_t earthLattitudeAngle, num_t dx, num_t dy);
diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/imaging/observatorium.h b/CEP/DP3/AOFlagger/include/AOFlagger/imaging/observatorium.h
index 80e1459f8eb..6e019a21233 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/imaging/observatorium.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/imaging/observatorium.h
@@ -19,35 +19,39 @@ class Observatorium
_channelWidthHz = channelWidthHz;
}
double ChannelWidthHz() const { return _channelWidthHz; }
+ const struct BandInfo &BandInfo() const
+ {
+ return _bandInfo;
+ }
+ protected:
+ struct BandInfo &GetBandInfo() { return _bandInfo; }
private:
std::vector<AntennaInfo> _antennae;
double _channelWidthHz;
+ struct BandInfo _bandInfo;
};
struct WSRTObservatorium : public Observatorium
{
WSRTObservatorium()
{
- SetChannelWidthHz(10000.0);
AntennaInfo antennas[14];
for(unsigned i=0;i<14;++i)
+ {
WSRTCommon(antennas[i]);
- WSRT0(antennas[0]);
- WSRT1(antennas[1]);
- WSRT2(antennas[2]);
- WSRT3(antennas[3]);
- WSRT4(antennas[4]);
- WSRT5(antennas[5]);
- WSRT6(antennas[6]);
- WSRT7(antennas[7]);
- WSRT8(antennas[8]);
- WSRT9(antennas[9]);
- WSRTA(antennas[10]);
- WSRTB(antennas[11]);
- WSRTC(antennas[12]);
- WSRTD(antennas[13]);
- for(unsigned i=0;i<14;++i)
+ WSRTn(i, antennas[i]);
AddAntenna(antennas[i]);
+ }
+ }
+ WSRTObservatorium(size_t antenna1, size_t antenna2)
+ {
+ AntennaInfo antennas[2];
+ WSRTCommon(antennas[0]);
+ WSRTCommon(antennas[1]);
+ WSRTn(antenna1, antennas[0]);
+ WSRTn(antenna2, antennas[1]);
+ AddAntenna(antennas[0]);
+ AddAntenna(antennas[1]);
}
private:
@@ -57,117 +61,123 @@ struct WSRTObservatorium : public Observatorium
antenna.mount = "equatorial";
antenna.station = "WSRT";
}
- void WSRT0(AntennaInfo &antenna)
- {
- antenna.id = 0;
- antenna.name = "RT0";
- antenna.position.x = 3.82876e+06;
- antenna.position.y = 442449;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT1(AntennaInfo &antenna)
- {
- antenna.id = 1;
- antenna.name = "RT1";
- antenna.position.x = 3.82875e+06;
- antenna.position.y = 442592;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT2(AntennaInfo &antenna)
- {
- antenna.id = 2;
- antenna.name = "RT2";
- antenna.position.x = 3.82873e+06;
- antenna.position.y = 442735;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT3(AntennaInfo &antenna)
- {
- antenna.id = 3;
- antenna.name = "RT3";
- antenna.position.x = 3.82871e+06;
- antenna.position.y = 442878;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT4(AntennaInfo &antenna)
- {
- antenna.id = 4;
- antenna.name = "RT4";
- antenna.position.x = 3.8287e+06;
- antenna.position.y = 443021;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT5(AntennaInfo &antenna)
- {
- antenna.id = 5;
- antenna.name = "RT5";
- antenna.position.x = 3.82868e+06;
- antenna.position.y = 443164;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT6(AntennaInfo &antenna)
- {
- antenna.id = 6;
- antenna.name = "RT6";
- antenna.position.x = 3.82866e+06;
- antenna.position.y = 443307;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT7(AntennaInfo &antenna)
- {
- antenna.id = 7;
- antenna.name = "RT7";
- antenna.position.x = 3.82865e+06;
- antenna.position.y = 443450;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT8(AntennaInfo &antenna)
- {
- antenna.id = 8;
- antenna.name = "RT8";
- antenna.position.x = 3.82863e+06;
- antenna.position.y = 443593;
- antenna.position.z = 5.06492e+06;
- }
- void WSRT9(AntennaInfo &antenna)
- {
- antenna.id = 9;
- antenna.name = "RT9";
- antenna.position.x = 3.82861e+06;
- antenna.position.y = 443736;
- antenna.position.z = 5.06492e+06;
- }
- void WSRTA(AntennaInfo &antenna)
- {
- antenna.id = 10;
- antenna.name = "RTA";
- antenna.position.x = 3.8286e+06;
- antenna.position.y = 443832;
- antenna.position.z = 5.06492e+06;
- }
- void WSRTB(AntennaInfo &antenna)
- {
- antenna.id = 11;
- antenna.name = "RTB";
- antenna.position.x = 3.82859e+06;
- antenna.position.y = 443903;
- antenna.position.z = 5.06492e+06;
- }
- void WSRTC(AntennaInfo &antenna)
- {
- antenna.id = 12;
- antenna.name = "RTC";
- antenna.position.x = 3.82845e+06;
- antenna.position.y = 445119;
- antenna.position.z = 5.06492e+06;
- }
- void WSRTD(AntennaInfo &antenna)
- {
- antenna.id = 13;
- antenna.name = "RTD";
- antenna.position.x = 3.82845e+06;
- antenna.position.y = 445191;
- antenna.position.z = 5.06492e+06;
+ void WSRTn(size_t antennaIndex, AntennaInfo &antenna)
+ {
+ switch(antennaIndex)
+ {
+ case 0:
+ antenna.id = 0;
+ antenna.name = "RT0";
+ antenna.position.x = 3.82876e+06;
+ antenna.position.y = 442449;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 1:
+ antenna.id = 1;
+ antenna.name = "RT1";
+ antenna.position.x = 3.82875e+06;
+ antenna.position.y = 442592;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 2:
+ antenna.id = 2;
+ antenna.name = "RT2";
+ antenna.position.x = 3.82873e+06;
+ antenna.position.y = 442735;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 3:
+ antenna.id = 3;
+ antenna.name = "RT3";
+ antenna.position.x = 3.82871e+06;
+ antenna.position.y = 442878;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 4:
+ antenna.id = 4;
+ antenna.name = "RT4";
+ antenna.position.x = 3.8287e+06;
+ antenna.position.y = 443021;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 5:
+ antenna.id = 5;
+ antenna.name = "RT5";
+ antenna.position.x = 3.82868e+06;
+ antenna.position.y = 443164;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 6:
+ antenna.id = 6;
+ antenna.name = "RT6";
+ antenna.position.x = 3.82866e+06;
+ antenna.position.y = 443307;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 7:
+ antenna.id = 7;
+ antenna.name = "RT7";
+ antenna.position.x = 3.82865e+06;
+ antenna.position.y = 443450;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 8:
+ antenna.id = 8;
+ antenna.name = "RT8";
+ antenna.position.x = 3.82863e+06;
+ antenna.position.y = 443593;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 9:
+ antenna.id = 9;
+ antenna.name = "RT9";
+ antenna.position.x = 3.82861e+06;
+ antenna.position.y = 443736;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 10:
+ antenna.id = 10;
+ antenna.name = "RTA";
+ antenna.position.x = 3.8286e+06;
+ antenna.position.y = 443832;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 11:
+ antenna.id = 11;
+ antenna.name = "RTB";
+ antenna.position.x = 3.82859e+06;
+ antenna.position.y = 443903;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 12:
+ antenna.id = 12;
+ antenna.name = "RTC";
+ antenna.position.x = 3.82845e+06;
+ antenna.position.y = 445119;
+ antenna.position.z = 5.06492e+06;
+ break;
+ case 13:
+ antenna.id = 13;
+ antenna.name = "RTD";
+ antenna.position.x = 3.82845e+06;
+ antenna.position.y = 445191;
+ antenna.position.z = 5.06492e+06;
+ break;
+ }
+ }
+ void initBand()
+ {
+ SetChannelWidthHz(10000.0);
+ GetBandInfo().windowIndex = 0;
+ GetBandInfo().channelCount = 256;
+ for(size_t i=0;i<256;++i)
+ {
+ ChannelInfo channel;
+ channel.frequencyIndex = i;
+ channel.channelWidthHz = ChannelWidthHz();
+ channel.frequencyHz = 147000000.0 + ChannelWidthHz() * i;
+ GetBandInfo().channels[i] = channel;
+ }
}
};
diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/rfi/rfistatistics.h b/CEP/DP3/AOFlagger/include/AOFlagger/rfi/rfistatistics.h
index 1448ceb2f1b..c05bf8ea5ba 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/rfi/rfistatistics.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/rfi/rfistatistics.h
@@ -165,7 +165,6 @@ class RFIStatistics {
saveChannels(_crossChannels, _filePrefix + "counts-channels-cross.txt");
saveAmplitudes(_crossAmplitudes, _filePrefix + "counts-amplitudes-cross.txt");
}
- saveWithoutBaselines(_filePrefix);
}
void Add(const ChannelInfo &channel, bool autocorrelation);
void Add(const TimestepInfo ×tep, bool autocorrelation);
diff --git a/CEP/DP3/AOFlagger/src/gui/mswindow.cpp b/CEP/DP3/AOFlagger/src/gui/mswindow.cpp
index 2bbc88ae79b..1fe8f1377a2 100644
--- a/CEP/DP3/AOFlagger/src/gui/mswindow.cpp
+++ b/CEP/DP3/AOFlagger/src/gui/mswindow.cpp
@@ -1222,14 +1222,18 @@ void MSWindow::onTFWidgetMouseMoved(size_t x, size_t y)
{
Image2DCPtr image = _timeFrequencyWidget.Image();
num_t v = image->Value(x, y);
- UVW uvw = _timeFrequencyWidget.GetMetaData()->UVW()[x];
_statusbar.pop();
std::stringstream s;
- s << "x=" << x << ",y=" << y << ",value=" << v;
- const std::vector<double> × = _timeFrequencyWidget.GetMetaData()->ObservationTimes();
- s << " (t=" << Date::AipsMJDToString(times[x]) <<
- ", f=" << Frequency::ToString(_timeFrequencyWidget.GetMetaData()->Band().channels[y].frequencyHz)
- << ", uvw=" << uvw.u << "," << uvw.v << "," << uvw.w << ')';
+ s << "x=" << x << ",y=" << y << ",value=" << v;
+ const std::vector<double> × = _timeFrequencyWidget.GetMetaData()->ObservationTimes();
+ s << " (t=" << Date::AipsMJDToString(times[x]) <<
+ ", f=" << Frequency::ToString(_timeFrequencyWidget.GetMetaData()->Band().channels[y].frequencyHz);
+ if(_timeFrequencyWidget.GetMetaData()->HasUVW())
+ {
+ UVW uvw = _timeFrequencyWidget.GetMetaData()->UVW()[x];
+ s << ", uvw=" << uvw.u << "," << uvw.v << "," << uvw.w;
+ }
+ s << ')';
_statusbar.push(s.str(), 0);
}
}
@@ -1367,10 +1371,6 @@ void MSWindow::showError(const std::string &description)
void MSWindow::onSimulateCorrelation()
{
Model model;
- //model.AddSource(-M_PIn - 0.04,0.04,0.5);
- //model.AddSource(-M_PIn - 0.04075,0.04075,0.2);
- //model.AddSource(-M_PIn + 0.1,0.0,0.35);
- //model.AddSource(-M_PIn + .101,0.001,0.45);
model.loadUrsaMajor();
WSRTObservatorium wsrtObservatorium;
@@ -1381,25 +1381,19 @@ void MSWindow::onSimulateCorrelation()
void MSWindow::onSimulateDoubledBaselineCorrelation()
{
Model model;
- //model.AddSource(-M_PIn - 0.04,0.04,0.5);
- //model.AddSource(-M_PIn - 0.04075,0.04075,0.2);
- //model.AddSource(-M_PIn + 0.1,0.0,0.35);
- //model.AddSource(-M_PIn + .101,0.001,0.45);
model.loadUrsaMajor();
- WSRTObservatorium wsrtObservatorium;
- FourProductCorrelatorTester fpcTester(model, *_imagePlaneWindow->GetImager(), wsrtObservatorium);
- fpcTester.SimulateTwoProdObservation(-M_PIn-0.05, 0.05, 147000000.0); //TwoProd
- _imagePlaneWindow->Update();
+ WSRTObservatorium wsrtObservatorium(0,5);
+ std::pair<TimeFrequencyData, TimeFrequencyMetaDataPtr> pair = model.SimulateObservation(wsrtObservatorium, -M_PIn-0.05, 0.05, 147000000.0);
+ TimeFrequencyData data = pair.first;
+ TimeFrequencyMetaDataCPtr metaData = pair.second;
+ _timeFrequencyWidget.SetNewData(data, metaData);
+ _timeFrequencyWidget.Update();
}
void MSWindow::onSimulateFourProductCorrelation()
{
Model model;
- //model.AddSource(-M_PIn - 0.04,0.04,0.5);
- //model.AddSource(-M_PIn - 0.04075,0.04075,0.2);
- //model.AddSource(-M_PIn + 0.1,0.0,0.35);
- //model.AddSource(-M_PIn + .101,0.001,0.45);
model.loadUrsaMajor();
WSRTObservatorium wsrtObservatorium;
diff --git a/CEP/DP3/AOFlagger/src/imaging/model.cpp b/CEP/DP3/AOFlagger/src/imaging/model.cpp
index 18c304422d3..743c98c753f 100644
--- a/CEP/DP3/AOFlagger/src/imaging/model.cpp
+++ b/CEP/DP3/AOFlagger/src/imaging/model.cpp
@@ -20,7 +20,8 @@
#include <AOFlagger/imaging/model.h>
#include <AOFlagger/msio/image2d.h>
-#include <AOFlagger/imaging/uvimager.h>
+#include <AOFlagger/msio/timefrequencydata.h>
+
#include <AOFlagger/imaging/observatorium.h>
#include <AOFlagger/util/rng.h>
@@ -34,13 +35,15 @@ Model::~Model()
{
}
-void Model::SimulateObservation(UVImager &imager, Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency)
+template<typename T>
+void Model::SimulateObservation(struct OutputReceiver<T> &receiver, Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency)
{
size_t frequencySteps = 1;
for(size_t f=0;f<frequencySteps;++f)
{
double channelFrequency = frequency + observatorium.ChannelWidthHz() * f * 256/frequencySteps;
+ receiver.SetY(f);
for(size_t i=0;i<observatorium.AntennaCount();++i)
{
for(size_t j=i+1;j<observatorium.AntennaCount();++j)
@@ -53,13 +56,35 @@ void Model::SimulateObservation(UVImager &imager, Observatorium &observatorium,
dy = antenna1.position.y - antenna2.position.y,
dz = antenna1.position.z - antenna2.position.z;
- SimulateCorrelation(imager, delayDirectionDEC, delayDirectionRA, dx, dy, dz, channelFrequency, 12*60*60, 10.0);
+ SimulateCorrelation(receiver, delayDirectionDEC, delayDirectionRA, dx, dy, dz, channelFrequency, 12*60*60, 10.0);
}
}
}
}
-void Model::SimulateCorrelation(UVImager &imager, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t frequency, double totalTime, double integrationTime)
+template void Model::SimulateObservation(struct OutputReceiver<UVImager> &receiver, Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency);
+template void Model::SimulateObservation(struct OutputReceiver<TimeFrequencyData> &receiver, Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency);
+
+std::pair<TimeFrequencyData, TimeFrequencyMetaDataPtr> Model::SimulateObservation(class Observatorium &observatorium, num_t delayDirectionDEC, num_t delayDirectionRA, num_t frequency)
+{
+ OutputReceiver<TimeFrequencyData> tfOutputter;
+ tfOutputter._real = Image2D::CreateZeroImagePtr(12*60*60/10, 1);
+ tfOutputter._imaginary = Image2D::CreateZeroImagePtr(12*60*60/10, 1);
+ SimulateObservation(tfOutputter, observatorium, delayDirectionDEC, delayDirectionRA, frequency);
+
+ TimeFrequencyMetaData *metaData = new TimeFrequencyMetaData();
+ metaData->SetAntenna1(observatorium.GetAntenna(0));
+ metaData->SetAntenna2(observatorium.GetAntenna(1));
+ metaData->SetBand(observatorium.BandInfo());
+ std::vector<double> times;
+ for(num_t t=0.0;t<12*60*60;t+=10) times.push_back(t);
+ metaData->SetObservationTimes(times);
+
+ return std::pair<TimeFrequencyData, TimeFrequencyMetaDataPtr>(TimeFrequencyData(StokesIPolarisation, tfOutputter._real, tfOutputter._imaginary), TimeFrequencyMetaDataPtr(metaData));
+}
+
+template<typename T>
+void Model::SimulateCorrelation(struct OutputReceiver<T> &receiver, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t frequency, double totalTime, double integrationTime)
{
num_t wavelength = 1.0L / frequency;
size_t index = 0;
@@ -73,12 +98,13 @@ void Model::SimulateCorrelation(UVImager &imager, num_t delayDirectionDEC, num_t
num_t
r = r1 * r2 - (i1 * -i2),
i = r1 * -i2 + r2 * i1;
- imager.SetUVValue(u, v, r, i, 1.0);
- imager.SetUVValue(-u, -v, r, -i, 1.0);
+ receiver.SetUVValue(index, u, v, r, i, 1.0);
++index;
}
}
+template void Model::SimulateCorrelation(struct OutputReceiver<UVImager> &receiver, num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t dz, num_t frequency, double totalTime, double integrationTime);
+
void Model::SimulateAntenna(num_t delayDirectionDEC, num_t delayDirectionRA, num_t dx, num_t dy, num_t frequency, num_t earthLattitude, num_t &r, num_t &i)
{
r = 0.0;
--
GitLab