diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/gui/mswindow.h b/CEP/DP3/AOFlagger/include/AOFlagger/gui/mswindow.h index 61c1be8c88b8f18a8f073f378d9bd9083706ced4..b9a927eb80dd4a5fc27d8916ff7170627e3672a7 100644 --- a/CEP/DP3/AOFlagger/include/AOFlagger/gui/mswindow.h +++ b/CEP/DP3/AOFlagger/include/AOFlagger/gui/mswindow.h @@ -163,6 +163,7 @@ class MSWindow : public Gtk::Window { void onOpenTestSetBAligned() { openTestSet(25); } void onOpenTestSetGaussianBroadband() { openTestSet(26); } void onOpenTestSetSinusoidalBroadband() { openTestSet(27); } + void onOpenTestSetSlewedGaussianBroadband() { openTestSet(28); } void onGaussianTestSets() { _gaussianTestSets = 1; } void onRayleighTestSets() { _gaussianTestSets = 0; } void onZeroTestSets() { _gaussianTestSets = 2; } diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/strategy/algorithms/mitigationtester.h b/CEP/DP3/AOFlagger/include/AOFlagger/strategy/algorithms/mitigationtester.h index 8d91e065b9e557b29d00ef5921bf73f3bfb82988..0b34ec34c38fd7972632be76376ac109c024adb0 100644 --- a/CEP/DP3/AOFlagger/include/AOFlagger/strategy/algorithms/mitigationtester.h +++ b/CEP/DP3/AOFlagger/include/AOFlagger/strategy/algorithms/mitigationtester.h @@ -42,6 +42,20 @@ class MitigationTester{ ~MitigationTester(); void GenerateNoise(size_t scanCount, size_t frequencyCount, bool independentComplex = false, double sigma=1.0, enum NoiseType noiseType=GaussianPartialProduct); + + static double shapeLevel(enum BroadbandShape shape, double x) + { + switch(shape) + { + default: + case UniformShape: + return 1.0; + case GaussianShape: + return exp(-x*x*3.0*3.0); + case SinusoidalShape: + return (1.0 + cos(x*M_PI*2.0*1.5)) * 0.5; + } + } static class Image2D *CreateRayleighData(unsigned width, unsigned height); static class Image2D *CreateGaussianData(unsigned width, unsigned height); @@ -65,6 +79,7 @@ class MitigationTester{ } static void AddBroadbandLine(Image2DPtr data, Mask2DPtr rfi, double lineStrength, size_t startTime, size_t duration, double frequencyRatio, double frequencyOffsetRatio); static void AddBroadbandLinePos(Image2DPtr data, Mask2DPtr rfi, double lineStrength, size_t startTime, size_t duration, unsigned frequencyStart, double frequencyEnd, enum BroadbandShape shape); + static void AddSlewedBroadbandLinePos(Image2DPtr data, Mask2DPtr rfi, double lineStrength, double slewrate, size_t startTime, size_t duration, unsigned frequencyStart, double frequencyEnd, enum BroadbandShape shape); static void AddRfiPos(Image2DPtr data, Mask2DPtr rfi, double lineStrength, size_t startTime, size_t duration, unsigned frequencyPos); void AddRFI(size_t &rfiCount); @@ -91,9 +106,14 @@ class MitigationTester{ { AddBroadbandToTestSet(image, rfi, 1.0, 1.0, false, SinusoidalShape); } + static void AddSlewedGaussianBroadbandToTestSet(Image2DPtr image, Mask2DPtr rfi) + { + AddSlewedBroadbandToTestSet(image, rfi, 1.0); + } private: static void AddBroadbandToTestSet(Image2DPtr image, Mask2DPtr rfi, long double length, double strength=1.0, bool align=false, enum BroadbandShape shape=UniformShape); + static void AddSlewedBroadbandToTestSet(Image2DPtr image, Mask2DPtr rfi, long double length, double strength=1.0, double slewrate=0.02, enum BroadbandShape shape=GaussianShape); static void AddVarBroadbandToTestSet(Image2DPtr image, Mask2DPtr rfi); static void AddModelData(Image2DPtr image, unsigned sources); static void SubtractBackground(Image2DPtr image); diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/test/experiments/rankoperatorrocexperiment.h b/CEP/DP3/AOFlagger/include/AOFlagger/test/experiments/rankoperatorrocexperiment.h index ecf22b0b847218d45ee1b67aa8ab6bdc16475044..d968096e899d24c3a42e657f9617b6636305eb6f 100644 --- a/CEP/DP3/AOFlagger/include/AOFlagger/test/experiments/rankoperatorrocexperiment.h +++ b/CEP/DP3/AOFlagger/include/AOFlagger/test/experiments/rankoperatorrocexperiment.h @@ -49,6 +49,7 @@ class RankOperatorROCExperiment : public UnitTest { { AddTest(RankOperatorROCGaussian(), "Constructing rank operator & dilation ROC curve for Gaussian broadband RFI"); AddTest(RankOperatorROCSinusoidal(), "Constructing rank operator & dilation ROC curve for sinusoidal broadband RFI"); + AddTest(RankOperatorROCSGaussian(), "Constructing rank operator & dilation ROC curve for slewed Gaussian RFI"); } private: @@ -75,10 +76,14 @@ class RankOperatorROCExperiment : public UnitTest { { void operator()(); }; + struct RankOperatorROCSGaussian : public Asserter + { + void operator()(); + }; static const unsigned _repeatCount; - enum TestType { GaussianBroadband, SinusoidalBroadband}; + enum TestType { GaussianBroadband, SinusoidalBroadband, SlewedGaussianBroadband }; static void TestNoisePerformance(size_t totalRFI, double totalRFISum, const std::string &testname); static void evaluateIterationResults(Mask2DPtr result, Mask2DPtr maskGroundTruth, Image2DPtr fuzzyGroundTruth, const size_t totalRFI, struct EvaluationResult &evaluationResult); @@ -169,6 +174,11 @@ inline void RankOperatorROCExperiment::RankOperatorROCSinusoidal::operator()() executeTest(SinusoidalBroadband); } +inline void RankOperatorROCExperiment::RankOperatorROCSGaussian::operator()() +{ + executeTest(SlewedGaussianBroadband); +} + void RankOperatorROCExperiment::evaluateIterationResults(Mask2DPtr result, Mask2DPtr maskGroundTruth, Image2DPtr fuzzyGroundTruth, const size_t totalRFI, EvaluationResult &evaluationResult) { size_t totalPositives = result->GetCount<true>(); @@ -218,6 +228,9 @@ void RankOperatorROCExperiment::executeTest(enum TestType testType) case SinusoidalBroadband: testname = "sinusoidal"; break; + case SlewedGaussianBroadband: + testname = "slewed_gaussian"; + break; } EvaluationResult rocResults[ETA_STEPS+1], dilResults[DIL_STEPS+1]; @@ -251,6 +264,23 @@ void RankOperatorROCExperiment::executeTest(enum TestType testType) groundTruth = Image2D::CreateCopy(TimeFrequencyData(SinglePolarisation, realTruth, imagTruth).GetSingleImage()); data = TimeFrequencyData(SinglePolarisation, realImage, imagImage); } break; + case SlewedGaussianBroadband: + { + Image2DPtr + realTruth = Image2D::CreateZeroImagePtr(width, height), + imagTruth = Image2D::CreateZeroImagePtr(width, height); + realImage = MitigationTester::CreateTestSet(2, mask, width, height), + imagImage = MitigationTester::CreateTestSet(2, mask, width, height); + rfiLessData = TimeFrequencyData(SinglePolarisation, realImage, imagImage); + rfiLessData.Trim(0, 0, 180, height); + rfiLessImage = rfiLessData.GetSingleImage(); + MitigationTester::AddSlewedGaussianBroadbandToTestSet(realImage, mask); + MitigationTester::AddSlewedGaussianBroadbandToTestSet(imagImage, mask); + MitigationTester::AddSlewedGaussianBroadbandToTestSet(realTruth, mask); + MitigationTester::AddSlewedGaussianBroadbandToTestSet(imagTruth, mask); + groundTruth = Image2D::CreateCopy(TimeFrequencyData(SinglePolarisation, realTruth, imagTruth).GetSingleImage()); + data = TimeFrequencyData(SinglePolarisation, realImage, imagImage); + } break; case SinusoidalBroadband: { Image2DPtr @@ -271,6 +301,7 @@ void RankOperatorROCExperiment::executeTest(enum TestType testType) data = TimeFrequencyData(SinglePolarisation, realImage, imagImage); } break; } + std::cout << '_' << std::flush; data.Trim(0, 0, 180, height); groundTruth->SetTrim(0, 0, 180, height); @@ -301,6 +332,8 @@ void RankOperatorROCExperiment::executeTest(enum TestType testType) evaluateIterationResults(resultMask, mask, groundTruth, totalRFI, rocResults[i]); } + std::cout << '.' << std::flush; + for(size_t i=0;i<DIL_STEPS;++i) { const size_t dilSize = i * MAX_DILATION / DIL_STEPS; @@ -315,7 +348,7 @@ void RankOperatorROCExperiment::executeTest(enum TestType testType) //grTotalRFISum += totalRFISum; std::cout << '.' << std::flush; - } + } // next repetition const std::string rankOperatorFilename(std::string("rank-operator-roc-") + testname + ".txt"); std::ofstream rankOperatorFile(rankOperatorFilename.c_str()); diff --git a/CEP/DP3/AOFlagger/src/gui/mswindow.cpp b/CEP/DP3/AOFlagger/src/gui/mswindow.cpp index 29ceafdd567196a8521bc2633012adfd1e90d016..0583f06489e79ee3aa51cbe2209a59baeb40e521 100644 --- a/CEP/DP3/AOFlagger/src/gui/mswindow.cpp +++ b/CEP/DP3/AOFlagger/src/gui/mswindow.cpp @@ -566,6 +566,8 @@ void MSWindow::createToolbar() sigc::mem_fun(*this, &MSWindow::onOpenTestSetGaussianBroadband)); _actionGroup->add( Gtk::Action::create("OpenTestSetSinusoidalBroadband", "Sinusoidal broadband"), sigc::mem_fun(*this, &MSWindow::onOpenTestSetSinusoidalBroadband)); + _actionGroup->add( Gtk::Action::create("OpenTestSetSlewedGaussianBroadband", "Slewed Gaussian"), + sigc::mem_fun(*this, &MSWindow::onOpenTestSetSlewedGaussianBroadband)); _actionGroup->add( Gtk::Action::create("AddTestModification", "Test modify") ); _actionGroup->add( Gtk::Action::create("AddStaticFringe", "Static fringe"), sigc::mem_fun(*this, &MSWindow::onAddStaticFringe) ); @@ -802,6 +804,7 @@ void MSWindow::createToolbar() " <menuitem action='OpenTestSetBAligned'/>" " <menuitem action='OpenTestSetGaussianBroadband'/>" " <menuitem action='OpenTestSetSinusoidalBroadband'/>" + " <menuitem action='OpenTestSetSlewedGaussianBroadband'/>" " </menu>" " <menu action='AddTestModification'>" " <menuitem action='AddStaticFringe'/>" diff --git a/CEP/DP3/AOFlagger/src/strategy/algorithms/mitigationtester.cpp b/CEP/DP3/AOFlagger/src/strategy/algorithms/mitigationtester.cpp index 4ddcc20650f84fa01d87a356c33b9ce6d05c3553..d8855b2aa7676c5f0e380bf5e115217626a82047 100644 --- a/CEP/DP3/AOFlagger/src/strategy/algorithms/mitigationtester.cpp +++ b/CEP/DP3/AOFlagger/src/strategy/algorithms/mitigationtester.cpp @@ -93,23 +93,10 @@ void MitigationTester::AddBroadbandLinePos(Image2DPtr data, Mask2DPtr rfi, doubl { const double s = (frequencyEnd-frequencyStart); for(size_t f=frequencyStart;f<frequencyEnd;++f) { + // x will run from -1 to 1 + const double x = (double) ((f-frequencyStart)*2)/s-1.0; + double factor = shapeLevel(shape, x); for(size_t t=startTime;t<startTime+duration;++t) { - // x will run from -1 to 1 - const double x = (double) ((f-frequencyStart)*2)/s-1.0; - double factor; - switch(shape) - { - default: - case UniformShape: - factor = 1.0; - break; - case GaussianShape: - factor = exp(-x*x*3.0*3.0); - break; - case SinusoidalShape: - factor = (1.0 + cos(x*M_PI*2.0*1.5)) * 0.5; - break; - } data->AddValue(t, f, lineStrength * factor); if(lineStrength > 0.0) rfi->SetValue(t, f, true); @@ -117,6 +104,31 @@ void MitigationTester::AddBroadbandLinePos(Image2DPtr data, Mask2DPtr rfi, doubl } } +void MitigationTester::AddSlewedBroadbandLinePos(Image2DPtr data, Mask2DPtr rfi, double lineStrength, double slewrate, size_t startTime, size_t duration, unsigned frequencyStart, double frequencyEnd, enum BroadbandShape shape) +{ + const double s = (frequencyEnd-frequencyStart); + for(size_t f=frequencyStart;f<frequencyEnd;++f) { + // x will run from -1 to 1 + const double x = (double) ((f-frequencyStart)*2)/s-1.0; + double factor = shapeLevel(shape, x); + double slew = slewrate * (double) f; + size_t slewInt = (size_t) slew; + double slewRest = slew - slewInt; + + data->AddValue(startTime+slewInt, f, lineStrength * factor * (1.0 - slewRest)); + if(lineStrength > 0.0) + rfi->SetValue(startTime+slewInt, f, true); + for(size_t t=startTime+1;t<startTime+duration;++t) { + data->AddValue(t+slewInt, f, lineStrength * factor); + if(lineStrength > 0.0) + rfi->SetValue(t+slewInt, f, true); + } + data->AddValue(startTime+duration+slewInt, f, lineStrength * factor * slewRest); + if(lineStrength > 0.0) + rfi->SetValue(startTime+duration+slewInt, f, true); + } +} + void MitigationTester::AddRfiPos(Image2DPtr data, Mask2DPtr rfi, double lineStrength, size_t startTime, size_t duration, unsigned frequencyPos) { for(size_t t=startTime;t<startTime+duration;++t) { @@ -454,6 +466,10 @@ Image2DPtr MitigationTester::CreateTestSet(int number, Mask2DPtr rfi, unsigned w image = Image2DPtr(CreateNoise(width, height, gaussianNoise)); AddBroadbandToTestSet(image, rfi, 1.0, 1.0, false, SinusoidalShape); } break; + case 28: { // Several slewed Gaussian broadband lines + image = Image2DPtr(CreateNoise(width, height, gaussianNoise)); + AddSlewedBroadbandToTestSet(image, rfi, 1.0); + } break; } return image; } @@ -514,6 +530,25 @@ void MitigationTester::AddBroadbandToTestSet(Image2DPtr image, Mask2DPtr rfi, lo } } +void MitigationTester::AddSlewedBroadbandToTestSet(Image2DPtr image, Mask2DPtr rfi, long double length, double strength, double slewrate, enum BroadbandShape shape) +{ + size_t frequencyCount = image->Height(); + unsigned step = image->Width()/11; + unsigned fStart = (unsigned) ((0.5 - length/2.0) * frequencyCount); + unsigned fEnd = (unsigned) ((0.5 + length/2.0) * frequencyCount); + AddSlewedBroadbandLinePos(image, rfi, 3.0*strength, slewrate, step*1, 3, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 2.5*strength, slewrate, step*2, 3, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 2.0*strength, slewrate, step*3, 3, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 1.8*strength, slewrate, step*4, 3, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 1.6*strength, slewrate, step*5, 3, fStart, fEnd, shape); + + AddSlewedBroadbandLinePos(image, rfi, 3.0*strength, slewrate, step*6, 1, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 2.5*strength, slewrate, step*7, 1, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 2.0*strength, slewrate, step*8, 1, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 1.8*strength, slewrate, step*9, 1, fStart, fEnd, shape); + AddSlewedBroadbandLinePos(image, rfi, 1.6*strength, slewrate, step*10, 1, fStart, fEnd, shape); +} + void MitigationTester::AddVarBroadbandToTestSet(Image2DPtr image, Mask2DPtr rfi) { // The "randomness" should be reproducable randomness, so calling