Bug 1491: adding time convolution operation

.gitattributes

@@ -232,6 +232,7 @@ CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/strategy.h -text
 CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/strategyiterator.h -text
 CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/svdaction.h -text
 CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/thresholdaction.h -text
+CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/timeconvolutionaction.h -text
 CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/timeselectionaction.h -text
 CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/types.h -text
 CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/writeflagsaction.h -text

CEP/DP3/AOFlagger/include/AOFlagger/rfi/strategy/timeconvolutionaction.h

+/***************************************************************************
+ *   Copyright (C) 2008 by A.R. Offringa   *
+ *   offringa@astro.rug.nl   *
+ *                                                                         *
+ *   This program is free software; you can redistribute it and/or modify  *
+ *   it under the terms of the GNU General Public License as published by  *
+ *   the Free Software Foundation; either version 2 of the License, or     *
+ *   (at your option) any later version.                                   *
+ *                                                                         *
+ *   This program is distributed in the hope that it will be useful,       *
+ *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
+ *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
+ *   GNU General Public License for more details.                          *
+ *                                                                         *
+ *   You should have received a copy of the GNU General Public License     *
+ *   along with this program; if not, write to the                         *
+ *   Free Software Foundation, Inc.,                                       *
+ *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
+ ***************************************************************************/
+#ifndef RFI_TIME_CONVOLUTION_ACTION
+#define RFI_TIME_CONVOLUTION_ACTION
+
+#include "artifactset.h"
+#include "actionblock.h"
+
+#include <AOFlagger/rfi/thresholdtools.h>
+
+namespace rfiStrategy {
+
+	class TimeConvolutionAction : public Action
+	{
+		public:
+			TimeConvolutionAction() : Action(), _sincSize(100.0)
+			{
+			}
+			virtual std::string Description()
+			{
+				return "Time convolution";
+			}
+			virtual ActionType Type() const { return AdapterType; }
+			virtual void Perform(ArtifactSet &artifacts, class ProgressListener &)
+			{
+				TimeFrequencyData data = artifacts.ContaminatedData();
+				Image2DCPtr image = data.GetSingleImage();
+				num_t *row = new num_t[image->Width()];
+				Image2DPtr newImage = Image2D::CreateEmptyImagePtr(image->Width(), image->Height());
+				for(unsigned y=0;y<image->Height();++y)
+				{
+					for(unsigned x=0;x<image->Width();++x)
+						row[x] = image->Value(x, y);
+					ThresholdTools::OneDimensionalSincConvolution(row, image->Width(), _sincSize);
+					for(unsigned x=0;x<image->Width();++x)
+						newImage->SetValue(x, y, row[x]);
+				}
+				delete[] row;
+				
+				TimeFrequencyData newRevisedData = TimeFrequencyData(data.PhaseRepresentation(), data.Polarisation(), newImage);
+				newRevisedData.SetMask(artifacts.RevisedData());
+
+				TimeFrequencyData *contaminatedData =
+					TimeFrequencyData::CreateTFDataFromDiff(artifacts.ContaminatedData(), newRevisedData);
+				contaminatedData->SetMask(artifacts.ContaminatedData());
+				artifacts.SetRevisedData(newRevisedData);
+				artifacts.SetContaminatedData(*contaminatedData);
+				delete contaminatedData;
+			}
+		private:
+			double _sincSize;
+	};
+
+} // namespace
+
+#endif // RFI_TIME_CONVOLUTION_ACTION

CEP/DP3/AOFlagger/include/AOFlagger/rfi/thresholdtools.h

@@ -46,6 +46,7 @@ class ThresholdTools{
 		static void CountMaskLengths(Mask2DCPtr mask, int *lengths, size_t lengthsSize);
 		static void OneDimensionalConvolution(num_t *data, unsigned dataSize, const num_t *kernel, unsigned kernelSize);
 		static void OneDimensionalGausConvolution(num_t *data, unsigned dataSize, num_t variance);
+		static void OneDimensionalSincConvolution(num_t *data, unsigned dataSize, num_t stretchFactor);
 		static void FilterConnectedSamples(Mask2DPtr mask, size_t minConnectedSampleArea, bool eightConnected=true);
 		static void FilterConnectedSample(Mask2DPtr mask, unsigned x, unsigned y, size_t minConnectedSampleArea, bool eightConnected=true);
 		static void UnrollPhase(Image2DPtr image);

CEP/DP3/AOFlagger/src/imaging/model.cpp

@@ -273,5 +273,5 @@ void Model::loadUrsaMajorDistortingSource()
 	double cr = 0.05;
 	double fluxoffset = 0.0;
 
-	AddSource(cd + s*rs*160, cr + s*300, 4.0 + fluxoffset); // Dubhe
+	AddSource(cd + s*rs*320, cr + s*600, 4.0 + fluxoffset); // Dubhe
 }

CEP/DP3/AOFlagger/src/rfi/strategy/actionfactory.cpp

@@ -43,6 +43,7 @@
 #include <AOFlagger/rfi/strategy/statisticalflagaction.h>
 #include <AOFlagger/rfi/strategy/svdaction.h>
 #include <AOFlagger/rfi/strategy/thresholdaction.h>
+#include <AOFlagger/rfi/strategy/timeconvolutionaction.h>
 #include <AOFlagger/rfi/strategy/timeselectionaction.h>
 #include <AOFlagger/rfi/strategy/writeflagsaction.h>
 
@@ -74,6 +75,7 @@ const std::vector<std::string> ActionFactory::GetActionList()
 	list.push_back("Spatial composition");
 	list.push_back("Statistical flagging");
 	list.push_back("Threshold");
+	list.push_back("Time convolution");
 	list.push_back("Time selection");
 	list.push_back("Write flags");
 	return list;
@@ -127,6 +129,8 @@ Action *ActionFactory::CreateAction(const std::string &action)
 		return new StatisticalFlagAction();
 	else if(action == "Threshold")
 		return new ThresholdAction();
+	else if(action == "Time convolution")
+		return new TimeConvolutionAction();
 	else if(action == "Time selection")
 		return new TimeSelectionAction();
 	else if(action == "Write flags")

CEP/DP3/AOFlagger/src/rfi/thresholdtools.cpp

@@ -424,6 +424,21 @@ void ThresholdTools::OneDimensionalGausConvolution(num_t *data, unsigned dataSiz
 	delete[] kernel;
 }
 
+void ThresholdTools::OneDimensionalSincConvolution(num_t *data, unsigned dataSize, num_t stretchFactor)
+{
+	num_t *kernel = new num_t[dataSize];
+	for(unsigned i=0;i<dataSize;++i)
+	{
+		num_t x = (((num_t) i-(num_t) dataSize/2.0) * 1.0 / stretchFactor);
+		if(x!=0.0)
+			kernel[i] = sinn(x) / x;
+		else
+			kernel[i] = 1.0;
+	}
+	OneDimensionalConvolution(data, dataSize, kernel, dataSize);
+	delete[] kernel;
+}
+
 num_t ThresholdTools::Mode(Image2DCPtr image, Mask2DCPtr mask)
 {
 	num_t mode = 0.0;