From f76b9ed9b6e5fb9a78cfa0e34c20a1afec43ff99 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Andr=C3=A9=20Offringa?= <offringa@astron.nl>
Date: Tue, 28 Feb 2012 00:07:26 +0000
Subject: [PATCH] Task #1892: variable step size for dnds
---
.../AOFlagger/gui/quality/histogrampage.h | 1 +
.../include/AOFlagger/quality/loghistogram.h | 19 +++++++++++++++++++
.../src/gui/quality/histogrampage.cpp | 15 +++++++++++----
3 files changed, 31 insertions(+), 4 deletions(-)
diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/gui/quality/histogrampage.h b/CEP/DP3/AOFlagger/include/AOFlagger/gui/quality/histogrampage.h
index 68d46a520b7..63cfa4ba21d 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/gui/quality/histogrampage.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/gui/quality/histogrampage.h
@@ -107,6 +107,7 @@ class HistogramPage : public Gtk::HBox {
Gtk::Frame _functionFrame;
Gtk::VBox _functionBox;
Gtk::RadioButton _nsButton, _dndsButton;
+ Gtk::Entry _deltaSEntry;
Gtk::Button _plotPropertiesButton, _dataExportButton;
diff --git a/CEP/DP3/AOFlagger/include/AOFlagger/quality/loghistogram.h b/CEP/DP3/AOFlagger/include/AOFlagger/quality/loghistogram.h
index 3daf1efc9fc..2f23c823a0d 100644
--- a/CEP/DP3/AOFlagger/include/AOFlagger/quality/loghistogram.h
+++ b/CEP/DP3/AOFlagger/include/AOFlagger/quality/loghistogram.h
@@ -291,6 +291,25 @@ class LogHistogram : public Serializable
return (double) sqrtl(sqErrorSum/(xsqErrorSum * (long double) (n-2)));
}
+ double NormalizedSlopeStdDevBySampling(double startAmplitude, double endAmplitude, double slope, double stepFactor) const
+ {
+ long double sum = 0.0;
+ unsigned long n = 0;
+ if(stepFactor <= 1.0001) stepFactor = 1.0001;
+ while(startAmplitude < endAmplitude)
+ {
+ const double stepEnd = startAmplitude * stepFactor;
+ double sampledSlope = NormalizedSlope(startAmplitude, stepEnd);
+ double sampleError = sampledSlope - slope;
+ sum += sampleError * sampleError;
+ ++n;
+
+ startAmplitude = stepEnd;
+ }
+
+ return (double) sqrtl(sum / ((long double) n*n - n));
+ }
+
double PowerLawUpperLimit(double constrainingAmplitude, double exponent, double factor) const
{
const double count = NormalizedCountAbove(constrainingAmplitude);
diff --git a/CEP/DP3/AOFlagger/src/gui/quality/histogrampage.cpp b/CEP/DP3/AOFlagger/src/gui/quality/histogrampage.cpp
index 4654bb8e4fd..8fc2cd56fba 100644
--- a/CEP/DP3/AOFlagger/src/gui/quality/histogrampage.cpp
+++ b/CEP/DP3/AOFlagger/src/gui/quality/histogrampage.cpp
@@ -49,6 +49,7 @@ HistogramPage::HistogramPage() :
_functionFrame("Function"),
_nsButton("N(S)"),
_dndsButton("dN(S)/dS"),
+ _deltaSEntry(),
_plotPropertiesButton("Properties"),
_dataExportButton("Data"),
_slopeFrame("Slope"),
@@ -123,6 +124,9 @@ HistogramPage::HistogramPage() :
_dndsButton.signal_clicked().connect(sigc::mem_fun(*this, &HistogramPage::updatePlot));
_dndsButton.set_group(group);
_nsButton.set_active(true);
+ _functionBox.pack_start(_deltaSEntry, Gtk::PACK_SHRINK);
+ _deltaSEntry.set_text("2");
+ _deltaSEntry.signal_activate().connect(sigc::mem_fun(*this, &HistogramPage::updatePlot));
_functionFrame.add(_functionBox);
_sideBox.pack_start(_functionFrame, Gtk::PACK_SHRINK);
@@ -341,13 +345,15 @@ void HistogramPage::plotFit(const LogHistogram &histogram, const std::string &ti
void HistogramPage::addHistogramToPlot(const LogHistogram &histogram)
{
const bool derivative = _dndsButton.get_active();
+ double deltaS = atof(_deltaSEntry.get_text().c_str());
+ if(deltaS <= 1.0001) deltaS = 1.0001;
for(LogHistogram::iterator i=histogram.begin();i!=histogram.end();++i)
{
if(derivative)
{
const double x = i.value();
const double logx = log10(x);
- const double cslope = histogram.NormalizedSlope(x*0.5, x*2.0);
+ const double cslope = histogram.NormalizedSlope(x/deltaS, x*deltaS);
if(std::isfinite(logx) && std::isfinite(cslope))
_plot.PushDataPoint(logx, cslope);
} else {
@@ -463,8 +469,6 @@ void HistogramPage::onDataExportClicked()
void HistogramPage::updateSlopeFrame(const LogHistogram &histogram)
{
std::stringstream str;
- str << "Slopes:";
-
addSlopeText(str, histogram, true);
_slopeTextView.get_buffer()->set_text(str.str());
@@ -472,6 +476,8 @@ void HistogramPage::updateSlopeFrame(const LogHistogram &histogram)
void HistogramPage::addSlopeText(std::stringstream &str, const LogHistogram &histogram, bool updateRange)
{
+ double deltaS = atof(_deltaSEntry.get_text().c_str());
+ if(deltaS <= 1.0001) deltaS = 1.0001;
double minRange, maxRange;
if(_slopeAutoRangeButton.get_active())
{
@@ -494,11 +500,12 @@ void HistogramPage::addSlopeText(std::stringstream &str, const LogHistogram &his
slope = histogram.NormalizedSlope(minRange, maxRange),
offset = histogram.NormalizedSlopeOffset(minRange, maxRange, slope),
error = histogram.NormalizedSlopeStdDev(minRange, maxRange, slope, offset),
+ errorB = histogram.NormalizedSlopeStdDevBySampling(minRange, maxRange, slope, deltaS),
upperLimit = histogram.PowerLawUpperLimit(minRange, slope, pow10(offset)),
lowerLimit = histogram.PowerLawLowerLimit(minRange, slope, pow10(offset), rfiRatio),
lowerError = fabs(lowerLimit - histogram.PowerLawLowerLimit(minRange, slope - error, pow10(offset), rfiRatio)),
lowerLimit2 = histogram.PowerLawLowerLimit2(minRange, slope, pow10(offset), rfiRatio);
- str << '\n' << slope << "±" << error << "\n["
+ str << '\n' << slope << "±" << error << "\n/±" << errorB << "\n["
<< log10(lowerLimit) << "±" << lowerError << ';' << log10(upperLimit) << ']' << '\n'
<< log10(lowerLimit2);
}
--
GitLab