From 0973f80680f120e951f2bcfef1461997b2035cea Mon Sep 17 00:00:00 2001
From: Thorsten Liebig <thorsten.liebig@gmx.de>
Date: Sun, 19 Dec 2010 20:41:08 +0100
Subject: [PATCH] new ProcessFieldsFD class + fixes

changes:
- new field processing class producing frequency domain field-dumps
- Processing needs dualMesh and dualTime information
- fix in TD field dumps using always dual-Time

Signed-off-by: Thorsten Liebig <thorsten.liebig@gmx.de>
---
 Common/processcurrent.cpp   |   2 -
 Common/processfields.cpp    |  77 ++++++++++++++
 Common/processfields.h      |   4 +
 Common/processfields_fd.cpp | 205 ++++++++++++++++++++++++++++++++++++
 Common/processfields_fd.h   |  41 ++++++++
 Common/processfields_td.cpp |   2 +-
 Common/processing.cpp       |   1 +
 Common/processing.h         |   6 ++
 Common/processintegral.cpp  |   3 +-
 Common/processintegral.h    |   3 -
 Common/processmodematch.cpp |   7 --
 openEMS.pro                 |   8 +-
 openems.cpp                 |  64 +++++++----
 13 files changed, 384 insertions(+), 39 deletions(-)
 create mode 100644 Common/processfields_fd.cpp
 create mode 100644 Common/processfields_fd.h

diff --git a/Common/processcurrent.cpp b/Common/processcurrent.cpp
index 923ff62..434e8d1 100644
--- a/Common/processcurrent.cpp
+++ b/Common/processcurrent.cpp
@@ -22,8 +22,6 @@
 
 ProcessCurrent::ProcessCurrent(Engine_Interface_Base* eng_if) : ProcessIntegral(eng_if)
 {
-	m_TimeShift = Op->GetTimestep()/2.0;
-	m_dualMesh = true;
 }
 
 ProcessCurrent::~ProcessCurrent()
diff --git a/Common/processfields.cpp b/Common/processfields.cpp
index 1346758..f9f25ff 100644
--- a/Common/processfields.cpp
+++ b/Common/processfields.cpp
@@ -28,6 +28,7 @@ ProcessFields::ProcessFields(Engine_Interface_Base* eng_if) : Processing(eng_if)
 	m_fileType = VTK_FILETYPE;
 	SetSubSampling(1);
 	SetPrecision(6);
+	m_dualTime = false;
 
 	for (int n=0; n<3; ++n)
 	{
@@ -84,6 +85,8 @@ void ProcessFields::InitProcess()
 
 		group = new H5::Group( file->createGroup( "/FieldData" ));
 		delete group;
+		group = new H5::Group( file->createGroup( "/FieldData/FD" ));
+		delete group;
 		delete file;
 	}
 }
@@ -461,6 +464,76 @@ bool ProcessFields::DumpVectorArray2HDF5(string filename, string name, FDTD_FLOA
 	return true;
 }
 
+bool ProcessFields::DumpVectorArray2HDF5(string filename, string name, std::complex<float> const* const* const* const* array, unsigned int const* numLines, float weight, float frequency)
+{
+	const H5std_string FILE_NAME(filename);
+	const H5std_string DATASET_NAME_RE( name + "_real");
+	const H5std_string DATASET_NAME_IM( name + "_imag");
+
+	H5::H5File file( FILE_NAME, H5F_ACC_RDWR );
+
+	H5::Group group( file.openGroup( "/FieldData/FD" ));
+
+	hsize_t t_dimsf[] = {1};
+	H5::DataSpace t_dataspace( 1, t_dimsf );
+
+	hsize_t dimsf[4];              // dataset dimensions
+	dimsf[0] = 3;
+	dimsf[1] = numLines[2];
+	dimsf[2] = numLines[1];
+	dimsf[3] = numLines[0];
+
+	H5::DataSpace dataspace( 4, dimsf );
+	H5::FloatType datatype( H5::PredType::NATIVE_FLOAT );
+
+	//create and write real part
+	H5::DataSet dataset = group.createDataSet( DATASET_NAME_RE, datatype, dataspace );
+	H5::Attribute attr = dataset.createAttribute("frequency",H5::PredType::NATIVE_FLOAT,t_dataspace);
+	attr.write( H5::PredType::NATIVE_FLOAT , &frequency);
+	// I have not the slightest idea why this array-copy action is necessary...  but it's the only way hdf5 does what it is supposed to do anyway!!
+	// at least it is save in case FDTD_FLOAT was defined as double...
+	// why does hdf5 write the dimensions backwards??? or matlab???
+	float hdf5array[3][numLines[2]][numLines[1]][numLines[0]];
+	for (int n=0; n<3; ++n)
+	{
+		for (unsigned int i=0; i<numLines[0]; ++i)
+		{
+			for (unsigned int j=0; j<numLines[1]; ++j)
+			{
+				for (unsigned int k=0; k<numLines[2]; ++k)
+				{
+					hdf5array[n][k][j][i] = array[n][i][j][k].real() * weight;
+				}
+			}
+		}
+	}
+	dataset.write( hdf5array, H5::PredType::NATIVE_FLOAT );
+
+	//create and write imaginary part
+	dataset = group.createDataSet( DATASET_NAME_IM, datatype, dataspace );
+	attr = dataset.createAttribute("frequency",H5::PredType::NATIVE_FLOAT,t_dataspace);
+	attr.write( H5::PredType::NATIVE_FLOAT , &frequency);
+	// I have not the slightest idea why this array-copy action is necessary...  but it's the only way hdf5 does what it is supposed to do anyway!!
+	// at least it is save in case FDTD_FLOAT was defined as double...
+	// why does hdf5 write the dimensions backwards??? or matlab???
+	for (int n=0; n<3; ++n)
+	{
+		for (unsigned int i=0; i<numLines[0]; ++i)
+		{
+			for (unsigned int j=0; j<numLines[1]; ++j)
+			{
+				for (unsigned int k=0; k<numLines[2]; ++k)
+				{
+					hdf5array[n][k][j][i] = array[n][i][j][k].imag() * weight;
+				}
+			}
+		}
+	}
+	dataset.write( hdf5array, H5::PredType::NATIVE_FLOAT );
+
+	return true;
+}
+
 FDTD_FLOAT**** ProcessFields::CalcField()
 {
 	unsigned int pos[3];
@@ -486,6 +559,7 @@ FDTD_FLOAT**** ProcessFields::CalcField()
 				}
 			}
 		}
+		return field;
 	}
 
 	if (m_DumpType==H_FIELD_DUMP)
@@ -506,7 +580,10 @@ FDTD_FLOAT**** ProcessFields::CalcField()
 				}
 			}
 		}
+		return field;
 	}
+
+	cerr << "ProcessFields::CalcField(): Error, unknown dump type..." << endl;
 	return field;
 }
 
diff --git a/Common/processfields.h b/Common/processfields.h
index 75a750c..9e8a404 100644
--- a/Common/processfields.h
+++ b/Common/processfields.h
@@ -71,8 +71,12 @@ public:
 	static bool DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT const* const* const* array, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12, string header_info = string(), MeshType meshT = CARTESIAN_MESH, double discLines_scaling = 1);
 	static bool DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT const* const* const* const* array, unsigned int numFields, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12, string header_info = string(), MeshType meshT = CARTESIAN_MESH, double discLines_scaling = 1);
 
+	//! Dump a time-domain vector dump to an HDF5 file
 	static bool DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines, float time=0);
 
+	//! Dump a frequency-domain complex-vector dump to an HDF5 file
+	static bool DumpVectorArray2HDF5(string filename, string name, std::complex<float> const* const* const* const* array, unsigned int const* numLines, float weight, float frequency);
+
 	double CalcTotalEnergy() const;
 
 	void SetFileType(FileType fileType) {m_fileType=fileType;}
diff --git a/Common/processfields_fd.cpp b/Common/processfields_fd.cpp
new file mode 100644
index 0000000..36bb88b
--- /dev/null
+++ b/Common/processfields_fd.cpp
@@ -0,0 +1,205 @@
+/*
+*	Copyright (C) 2010 Thorsten Liebig (Thorsten.Liebig@gmx.de)
+*
+*	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 3 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, see <http://www.gnu.org/licenses/>.
+*/
+
+#include "processfields_fd.h"
+#include "Common/operator_base.h"
+#include <iomanip>
+#include <sstream>
+#include <string>
+
+ProcessFieldsFD::ProcessFieldsFD(Engine_Interface_Base* eng_if) : ProcessFields(eng_if)
+{
+}
+
+ProcessFieldsFD::~ProcessFieldsFD()
+{
+	for (size_t n = 0; n<m_FD_Fields.size(); ++n)
+	{
+		Delete_N_3DArray(m_FD_Fields.at(n),numLines);
+	}
+	m_FD_Fields.clear();
+}
+
+void ProcessFieldsFD::InitProcess()
+{
+	if (Enabled==false) return;
+
+	if (m_FD_Samples.size()==0)
+	{
+		cerr << "ProcessFieldsFD::InitProcess: No frequencies found... skipping this dump!" << endl;
+		Enabled=false;
+		return;
+	}
+
+	//setup the hdf5 file
+	ProcessFields::InitProcess();
+
+	//create data structures...
+	for (size_t n = 0; n<m_FD_Samples.size(); ++n)
+	{
+		std::complex<float>**** field_fd = Create_N_3DArray<std::complex<float> >(numLines);
+		m_FD_Fields.push_back(field_fd);
+	}
+}
+
+int ProcessFieldsFD::Process()
+{
+	if (Enabled==false) return -1;
+	if (CheckTimestep()==false) return GetNextInterval();
+
+	if ((m_FD_Interval==0) || (m_Eng_Interface->GetNumberOfTimesteps()%m_FD_Interval!=0))
+		return GetNextInterval();
+
+	FDTD_FLOAT**** field_td = CalcField();
+	std::complex<float>**** field_fd = NULL;
+
+	double T = m_Eng_Interface->GetTime(m_dualTime);
+	unsigned int pos[3];
+	for (size_t n = 0; n<m_FD_Samples.size(); ++n)
+	{
+		std::complex<float> exp_jwt = std::exp( (std::complex<float>)(-2.0 * _I * M_PI * m_FD_Samples.at(n) * T) );
+		field_fd = m_FD_Fields.at(n);
+		for (pos[0]=0; pos[0]<numLines[0]; ++pos[0])
+		{
+			for (pos[1]=0; pos[1]<numLines[1]; ++pos[1])
+			{
+				for (pos[2]=0; pos[2]<numLines[2]; ++pos[2])
+				{
+					field_fd[0][pos[0]][pos[1]][pos[2]] += field_td[0][pos[0]][pos[1]][pos[2]] * exp_jwt;
+					field_fd[1][pos[0]][pos[1]][pos[2]] += field_td[1][pos[0]][pos[1]][pos[2]] * exp_jwt;
+					field_fd[2][pos[0]][pos[1]][pos[2]] += field_td[2][pos[0]][pos[1]][pos[2]] * exp_jwt;
+				}
+			}
+		}
+	}
+	++m_FD_SampleCount;
+	return GetNextInterval();
+}
+
+void ProcessFieldsFD::PostProcess()
+{
+	DumpFDData();
+}
+
+void ProcessFieldsFD::DumpFDData()
+{
+#ifdef OUTPUT_IN_DRAWINGUNITS
+	double discLines_scaling = 1;
+#else
+	double discLines_scaling = Op->GetGridDelta();
+#endif
+
+	if (m_fileType==VTK_FILETYPE)
+	{
+		unsigned int pos[3];
+		FDTD_FLOAT**** field = Create_N_3DArray<float>(numLines);
+		std::complex<float>**** field_fd = NULL;
+		double angle=0;
+		int Nr_Ph = 21;
+
+		for (size_t n = 0; n<m_FD_Samples.size(); ++n)
+		{
+			//dump multiple phase to vtk-files
+			for (int p=0; p<Nr_Ph; ++p)
+			{
+				angle = 2.0 * M_PI * p / Nr_Ph;
+				std::complex<float> exp_jwt = std::exp( (std::complex<float>)( _I * angle) );
+				field_fd = m_FD_Fields.at(n);
+				for (pos[0]=0; pos[0]<numLines[0]; ++pos[0])
+				{
+					for (pos[1]=0; pos[1]<numLines[1]; ++pos[1])
+					{
+						for (pos[2]=0; pos[2]<numLines[2]; ++pos[2])
+						{
+							field[0][pos[0]][pos[1]][pos[2]] = real(field_fd[0][pos[0]][pos[1]][pos[2]] * exp_jwt)/(float)m_FD_SampleCount;
+							field[1][pos[0]][pos[1]][pos[2]] = real(field_fd[1][pos[0]][pos[1]][pos[2]] * exp_jwt)/(float)m_FD_SampleCount;
+							field[2][pos[0]][pos[1]][pos[2]] = real(field_fd[2][pos[0]][pos[1]][pos[2]] * exp_jwt)/(float)m_FD_SampleCount;
+						}
+					}
+				}
+				stringstream ss;
+				ss << m_filename << fixed << "_f=" << m_FD_Samples.at(n) << "_p=" << std::setw( 3 ) << std::setfill( '0' ) <<(int)(angle * 180 / M_PI) << ".vtk";
+				ofstream file(ss.str().c_str());
+				if (file.is_open()==false)
+					cerr << "ProcessFieldsFD::DumpFDData: can't open file '" << ss.str() << "' for writing... abort! " << endl;
+				DumpVectorArray2VTK(file,GetFieldNameByType(m_DumpType),field,discLines,numLines,m_precision,string("Interpolation: ")+m_Eng_Interface->GetInterpolationTypeString(), m_Mesh_Type, discLines_scaling);
+				file.close();
+			}
+
+			{
+				//dump magnitude to vtk-files
+				for (pos[0]=0; pos[0]<numLines[0]; ++pos[0])
+				{
+					for (pos[1]=0; pos[1]<numLines[1]; ++pos[1])
+					{
+						for (pos[2]=0; pos[2]<numLines[2]; ++pos[2])
+						{
+							field[0][pos[0]][pos[1]][pos[2]] = abs(field_fd[0][pos[0]][pos[1]][pos[2]])/(float)m_FD_SampleCount;
+							field[1][pos[0]][pos[1]][pos[2]] = abs(field_fd[1][pos[0]][pos[1]][pos[2]])/(float)m_FD_SampleCount;
+							field[2][pos[0]][pos[1]][pos[2]] = abs(field_fd[2][pos[0]][pos[1]][pos[2]])/(float)m_FD_SampleCount;
+						}
+					}
+				}
+				stringstream ss;
+				ss << m_filename << fixed << "_f=" << m_FD_Samples.at(n) << "_abs" << ".vtk";
+				ofstream file(ss.str().c_str());
+				if (file.is_open()==false)
+					cerr << "ProcessFieldsFD::DumpFDData: can't open file '" << ss.str() << "' for writing... abort! " << endl;
+				DumpVectorArray2VTK(file,GetFieldNameByType(m_DumpType),field,discLines,numLines,m_precision,string("Interpolation: ")+m_Eng_Interface->GetInterpolationTypeString(), m_Mesh_Type, discLines_scaling);
+				file.close();
+			}
+
+			{
+				//dump phase to vtk-files
+				for (pos[0]=0; pos[0]<numLines[0]; ++pos[0])
+				{
+					for (pos[1]=0; pos[1]<numLines[1]; ++pos[1])
+					{
+						for (pos[2]=0; pos[2]<numLines[2]; ++pos[2])
+						{
+							field[0][pos[0]][pos[1]][pos[2]] = arg(field_fd[0][pos[0]][pos[1]][pos[2]]);
+							field[1][pos[0]][pos[1]][pos[2]] = arg(field_fd[1][pos[0]][pos[1]][pos[2]]);
+							field[2][pos[0]][pos[1]][pos[2]] = arg(field_fd[2][pos[0]][pos[1]][pos[2]]);
+						}
+					}
+				}
+				stringstream ss;
+				ss << m_filename << fixed << "_f=" << m_FD_Samples.at(n) << "_arg" << ".vtk";
+				ofstream file(ss.str().c_str());
+				if (file.is_open()==false)
+					cerr << "ProcessFieldsFD::DumpFDData: can't open file '" << ss.str() << "' for writing... abort! " << endl;
+				DumpVectorArray2VTK(file,GetFieldNameByType(m_DumpType),field,discLines,numLines,m_precision,string("Interpolation: ")+m_Eng_Interface->GetInterpolationTypeString(), m_Mesh_Type, discLines_scaling);
+				file.close();
+			}
+		}
+		Delete_N_3DArray(field,numLines);
+		return;
+	}
+
+	if (m_fileType==HDF5_FILETYPE)
+	{
+		for (size_t n = 0; n<m_FD_Samples.size(); ++n)
+		{
+			stringstream ss;
+			ss << "f" << n;
+			DumpVectorArray2HDF5(m_filename.c_str(), ss.str(), m_FD_Fields.at(n),numLines,1.0/(float)m_FD_SampleCount,m_FD_Samples.at(n));
+		}
+		return;
+	}
+
+	cerr << "ProcessFieldsTD::Process: unknown File-Type" << endl;
+}
diff --git a/Common/processfields_fd.h b/Common/processfields_fd.h
new file mode 100644
index 0000000..3a47474
--- /dev/null
+++ b/Common/processfields_fd.h
@@ -0,0 +1,41 @@
+/*
+*	Copyright (C) 2010 Thorsten Liebig (Thorsten.Liebig@gmx.de)
+*
+*	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 3 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, see <http://www.gnu.org/licenses/>.
+*/
+
+#ifndef PROCESSFIELDS_FD_H
+#define PROCESSFIELDS_FD_H
+
+#include "processfields.h"
+
+class ProcessFieldsFD : public ProcessFields
+{
+public:
+	ProcessFieldsFD(Engine_Interface_Base* eng_if);
+	virtual ~ProcessFieldsFD();
+
+	virtual void InitProcess();
+
+	virtual int Process();
+	virtual void PostProcess();
+
+protected:
+	void DumpFDData();
+
+	//! frequency domain field storage
+	vector<std::complex<float>****> m_FD_Fields;
+};
+
+#endif // PROCESSFIELDS_FD_H
diff --git a/Common/processfields_td.cpp b/Common/processfields_td.cpp
index 107afa6..e9d1ad9 100644
--- a/Common/processfields_td.cpp
+++ b/Common/processfields_td.cpp
@@ -69,7 +69,7 @@ int ProcessFieldsTD::Process()
 	{
 		stringstream ss;
 		ss << std::setw( pad_length ) << std::setfill( '0' ) << m_Eng_Interface->GetNumberOfTimesteps();
-		DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),field,numLines,(0.5+m_Eng_Interface->GetNumberOfTimesteps())*Op->GetTimestep());
+		DumpVectorArray2HDF5(filename.c_str(),string( ss.str() ),field,numLines,m_Eng_Interface->GetTime(m_dualTime));
 	}
 	else
 		cerr << "ProcessFieldsTD::Process: unknown File-Type" << endl;
diff --git a/Common/processing.cpp b/Common/processing.cpp
index 9d18349..9697990 100644
--- a/Common/processing.cpp
+++ b/Common/processing.cpp
@@ -37,6 +37,7 @@ Processing::Processing(Engine_Interface_Base* eng_if)
 	m_weight=1;
 	m_Flush = false;
 	m_dualMesh = false;
+	m_dualTime = false;
 	m_Mesh_Type = CARTESIAN_MESH;
 }
 
diff --git a/Common/processing.h b/Common/processing.h
index c262420..1807b77 100644
--- a/Common/processing.h
+++ b/Common/processing.h
@@ -97,6 +97,9 @@ public:
 	//! Dump probe geometry to file
 	virtual void DumpBox2File( string vtkfilenameprefix, bool dualMesh) const;
 
+	virtual void SetDualMesh(bool val) {m_dualMesh=val;}
+	virtual void SetDualTime(bool val) {m_dualTime=val;}
+
 protected:
 	Engine_Interface_Base* m_Eng_Interface;
 	const Operator_Base* Op;
@@ -130,6 +133,9 @@ protected:
 	//! define if given coords are on main or dualMesh (default is false)
 	bool m_dualMesh;
 
+	//! define if given processing uses the dual time concept (default is false);
+	bool m_dualTime;
+
 	//! define/store snapped start/stop coords as mesh index
 	unsigned int start[3];
 	unsigned int stop[3];
diff --git a/Common/processintegral.cpp b/Common/processintegral.cpp
index deeb926..52c4510 100644
--- a/Common/processintegral.cpp
+++ b/Common/processintegral.cpp
@@ -20,7 +20,6 @@
 
 ProcessIntegral::ProcessIntegral(Engine_Interface_Base* eng_if)  : Processing(eng_if)
 {
-	m_TimeShift = 0.0;
 	m_Results=NULL;
 }
 
@@ -56,7 +55,7 @@ int ProcessIntegral::Process()
 	int NrInt = GetNumberOfIntegrals();
 	double integral = m_Results[0] * m_weight;
 
-	double time = m_Eng_Interface->GetTime(false) + m_TimeShift;
+	double time = m_Eng_Interface->GetTime(m_dualTime);
 
 	if (ProcessInterval)
 	{
diff --git a/Common/processintegral.h b/Common/processintegral.h
index c6fda5e..b44fbeb 100644
--- a/Common/processintegral.h
+++ b/Common/processintegral.h
@@ -53,9 +53,6 @@ public:
 protected:
 	ProcessIntegral(Engine_Interface_Base* eng_if);
 
-	//! timeshift to be used in TD and FD data, e.g. 0.5*dT in case of current based parameter
-	double m_TimeShift;
-
 	vector<double> TD_Values;
 	vector<double_complex> FD_Values;
 
diff --git a/Common/processmodematch.cpp b/Common/processmodematch.cpp
index feb9034..1cfc68c 100644
--- a/Common/processmodematch.cpp
+++ b/Common/processmodematch.cpp
@@ -27,8 +27,6 @@ ProcessModeMatch::ProcessModeMatch(Engine_Interface_Base* eng_if) : ProcessInteg
 		m_ModeParser[n] = new CSFunctionParser();
 		m_ModeDist[n] = NULL;
 	}
-	m_dualMesh = false;
-
 	delete[] m_Results;
 	m_Results = new double[2];
 }
@@ -47,11 +45,6 @@ void ProcessModeMatch::InitProcess()
 {
 	if (!Enabled) return;
 
-	if (m_ModeFieldType==1)
-	{
-		m_TimeShift = Op->GetTimestep()/2.0;
-	}
-
 	if (m_Eng_Interface==NULL)
 	{
 		cerr << "ProcessModeMatch::InitProcess: Error, Engine_Interface is NULL, abort mode mathcing..." << endl;
diff --git a/openEMS.pro b/openEMS.pro
index a85eaa6..24a64e7 100644
--- a/openEMS.pro
+++ b/openEMS.pro
@@ -64,7 +64,8 @@ SOURCES += main.cpp \
     FDTD/excitation.cpp \
     FDTD/operator_cylindermultigrid.cpp \
     FDTD/engine_cylindermultigrid.cpp \
-    FDTD/engine_interface_fdtd.cpp
+    FDTD/engine_interface_fdtd.cpp \
+    Common/processfields_fd.cpp
 
 # FDTD/extensions source files
 SOURCES += FDTD/extensions/engine_extension.cpp \
@@ -116,9 +117,10 @@ HEADERS += tools/ErrorMsg.h \
     FDTD/operator_cylindermultigrid.h \
     FDTD/engine_cylindermultigrid.h \
     tools/aligned_allocator.h \
-    FDTD/engine_interface_fdtd.h  
+    FDTD/engine_interface_fdtd.h \
+    Common/processfields_fd.h
 
-# FDTD/extensions header files   
+# FDTD/extensions header files
 HEADERS += FDTD/operator_extension.h \
     FDTD/extensions/engine_extension.h \
     FDTD/extensions/engine_ext_mur_abc.h \
diff --git a/openems.cpp b/openems.cpp
index b57565d..b880131 100644
--- a/openems.cpp
+++ b/openems.cpp
@@ -33,6 +33,7 @@
 #include "Common/process_hfield.h"
 #include "Common/processmodematch.h"
 #include "Common/processfields_td.h"
+#include "Common/processfields_fd.h"
 #include <sys/time.h>
 #include <time.h>
 #include <H5Cpp.h> // only for H5get_libversion()
@@ -488,6 +489,11 @@ int openEMS::SetupFDTD(const char* file)
 				}
 				if (CylinderCoords)
 					proc->SetMeshType(Processing::CYLINDRICAL_MESH);
+				if ((pb->GetProbeType()==1) || (pb->GetProbeType()==3) || (pb->GetProbeType()==11))
+				{
+					proc->SetDualTime(true);
+					proc->SetDualMesh(true);
+				}
 				proc->SetProcessInterval(Nyquist/m_OverSampling);
 				proc->AddFrequency(pb->GetFDSamples());
 				proc->SetName(pb->GetName());
@@ -505,15 +511,11 @@ int openEMS::SetupFDTD(const char* file)
 	vector<CSProperties*> DumpProps = CSX.GetPropertyByType(CSProperties::DUMPBOX);
 	for (size_t i=0; i<DumpProps.size(); ++i)
 	{
-		ProcessFieldsTD* ProcTD = new ProcessFieldsTD(new Engine_Interface_FDTD(FDTD_Op,FDTD_Eng));
-		ProcTD->SetEnable(Enable_Dumps);
-		ProcTD->SetProcessInterval(Nyquist/m_OverSampling);
+		ProcessFields* ProcField=NULL;
 
 		//only looking for one prim atm
 		CSPrimitives* prim = DumpProps.at(i)->GetPrimitive(0);
-		if (prim==NULL)
-			delete 	ProcTD;
-		else
+		if (prim!=NULL)
 		{
 			bool acc;
 			double bnd[6] = {0,0,0,0,0,0};
@@ -527,22 +529,42 @@ int openEMS::SetupFDTD(const char* file)
 			CSPropDumpBox* db = DumpProps.at(i)->ToDumpBox();
 			if (db)
 			{
-				ProcTD->SetDumpType((ProcessFields::DumpType)db->GetDumpType());
-				ProcTD->SetDumpMode((Engine_Interface_Base::InterpolationType)db->GetDumpMode());
-				ProcTD->SetFileType((ProcessFields::FileType)db->GetFileType());
-				if (CylinderCoords)
-					ProcTD->SetMeshType(Processing::CYLINDRICAL_MESH);
-				for (int n=0; n<3; ++n)
-					ProcTD->SetSubSampling(db->GetSubSampling(n),n);
-				ProcTD->SetFilePattern(db->GetName());
-				ProcTD->SetFileName(db->GetName());
-				ProcTD->DefineStartStopCoord(start,stop);
-				ProcTD->InitProcess();
-				PA->AddProcessing(ProcTD);
-				prim->SetPrimitiveUsed(true);
+				if ((db->GetDumpType()>=0) && (db->GetDumpType()<2))
+					ProcField = new ProcessFieldsTD(new Engine_Interface_FDTD(FDTD_Op,FDTD_Eng));
+				else if  ((db->GetDumpType()>=10) && (db->GetDumpType()<12))
+					ProcField = new ProcessFieldsFD(new Engine_Interface_FDTD(FDTD_Op,FDTD_Eng));
+				else
+					cerr << "openEMS::SetupFDTD: unknown dump box type... skipping!" << endl;
+				if (ProcField)
+				{
+					ProcField->SetEnable(Enable_Dumps);
+					ProcField->SetProcessInterval(Nyquist/m_OverSampling);
+					if ((db->GetDumpType()==1) || (db->GetDumpType()==11))
+					{
+						ProcField->SetDualTime(true);
+						ProcField->SetDualMesh(true);
+					}
+					if ((db->GetDumpType()==10) || (db->GetDumpType()==11))
+						ProcField->AddFrequency(db->GetFDSamples());
+					if (db->GetDumpType()>=10)
+						ProcField->SetDumpType((ProcessFields::DumpType)(db->GetDumpType()-10));
+					else
+						ProcField->SetDumpType((ProcessFields::DumpType)db->GetDumpType());
+
+					ProcField->SetDumpMode((Engine_Interface_Base::InterpolationType)db->GetDumpMode());
+					ProcField->SetFileType((ProcessFields::FileType)db->GetFileType());
+					if (CylinderCoords)
+						ProcField->SetMeshType(Processing::CYLINDRICAL_MESH);
+					for (int n=0; n<3; ++n)
+						ProcField->SetSubSampling(db->GetSubSampling(n),n);
+					ProcField->SetFilePattern(db->GetName());
+					ProcField->SetFileName(db->GetName());
+					ProcField->DefineStartStopCoord(start,stop);
+					ProcField->InitProcess();
+					PA->AddProcessing(ProcField);
+					prim->SetPrimitiveUsed(true);
+				}
 			}
-			else
-				delete 	ProcTD;
 		}
 	}