diff --git a/nf2ff/main.cpp b/nf2ff/main.cpp
new file mode 100644
index 0000000..1dc7cfe
--- /dev/null
+++ b/nf2ff/main.cpp
@@ -0,0 +1,44 @@
+/*
+*	Copyright (C) 2012 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 <iostream>
+
+#include "nf2ff.h"
+
+using namespace std;
+
+int main(int argc, char *argv[])
+{
+	cout << " ---------------------------------------------------------------------- " << endl;
+	cout << " | nf2ff, near-field to far-field transformation for openEMS " << endl;
+	cout << " | (C) 2012 Thorsten Liebig <thorsten.liebig@gmx.de>  GPL license" << endl;
+	cout << " ---------------------------------------------------------------------- " << endl;
+
+	if (argc<=1)
+	{
+		cout << " Usage: nf2ff <nf2ff-xml-file>" << endl << endl;
+		cout << endl;
+		exit(-1);
+	}
+
+	if (argc>=2)
+	{
+		return !nf2ff::AnalyseXMLFile(argv[argc-1]);
+	}
+
+	return 0;
+}
diff --git a/nf2ff/nf2ff.cpp b/nf2ff/nf2ff.cpp
new file mode 100644
index 0000000..17b2f7c
--- /dev/null
+++ b/nf2ff/nf2ff.cpp
@@ -0,0 +1,365 @@
+/*
+*	Copyright (C) 2012 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 "nf2ff.h"
+#include "nf2ff_calc.h"
+#include "../tools/array_ops.h"
+#include "../tools/useful.h"
+#include "../tools/hdf5_file_reader.h"
+#include "../tools/hdf5_file_writer.h"
+#include <hdf5.h>
+#include <boost/algorithm/string.hpp>
+#include <stdio.h>
+#include <stdlib.h>
+#include <vector>
+#include <cmath>
+#include <complex>
+#include <iostream>
+#include <sstream>
+
+//external libs
+#include "tinyxml.h"
+
+nf2ff::nf2ff(vector<float> freq, vector<float> theta, vector<float> phi, unsigned int numThreads)
+{
+	m_freq = freq;
+
+	m_numTheta = theta.size();
+	m_theta = new float[m_numTheta];
+	for (size_t n=0;n<m_numTheta;++n)
+		 m_theta[n]=theta.at(n);
+
+	m_numPhi = phi.size();
+	m_phi = new float[m_numPhi];
+	for (size_t n=0;n<m_numPhi;++n)
+		 m_phi[n]=phi.at(n);
+
+	m_nf2ff.resize(freq.size(),NULL);
+	for (size_t fn=0;fn<freq.size();++fn)
+	{
+		m_nf2ff.at(fn) = new nf2ff_calc(freq.at(fn),theta, phi);
+		if (numThreads)
+			m_nf2ff.at(fn)->SetNumThreads(numThreads);
+	}
+
+	m_radius = 1;
+	m_Verbose = 0;
+}
+
+nf2ff::~nf2ff()
+{
+	m_freq.clear();
+	for (size_t fn=0;fn<m_nf2ff.size();++fn)
+		delete m_nf2ff.at(fn);
+	m_nf2ff.clear();
+
+	delete[] m_phi;
+	m_phi = NULL;
+	delete[] m_theta;
+	m_theta = NULL;
+}
+
+bool nf2ff::AnalyseXMLNode(TiXmlElement* ti_nf2ff)
+{
+	if (ti_nf2ff==NULL)
+		return false;
+
+	unsigned int numThreads=0;
+	int ihelp=0;
+	if (ti_nf2ff->QueryIntAttribute("NumThreads",&ihelp) == TIXML_SUCCESS)
+	{
+		numThreads = ihelp;
+		cerr << "nf2ff: Set number of threads to: " << numThreads << endl;
+	}
+	int Verbose=0;
+	if (ti_nf2ff->QueryIntAttribute("Verbose",&Verbose) == TIXML_SUCCESS)
+		cerr << "nf2ff: Set verbose level to " << Verbose << endl;
+	else
+		Verbose = 0;
+
+	const char* attr = NULL;
+	attr = ti_nf2ff->Attribute("freq");
+	if (attr==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read frequency inforamtions ... " << endl;
+		return false;
+	}
+	vector<float> freq = SplitString2Float(attr);
+	attr = ti_nf2ff->Attribute("Outfile");
+	if (attr==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read frequency inforamtions ... " << endl;
+		return false;
+	}
+	string outfile = string(attr);
+	if (outfile.empty())
+	{
+		cerr << "nf2ff::AnalyseXMLNode: outfile is empty, skipping nf2ff... " << endl;
+		return false;
+	}
+
+	TiXmlElement* ti_theta = ti_nf2ff->FirstChildElement("theta");
+	if (ti_theta==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read theta values ... " << endl;
+		return false;
+	}
+	TiXmlNode* ti_theta_node = ti_theta->FirstChild();
+	if (ti_theta_node==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read theta text child ... " << endl;
+		return false;
+	}
+	TiXmlText* ti_theta_text = ti_theta_node->ToText();
+	if (ti_theta_text==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read theta text values ... " << endl;
+		return false;
+	}
+	vector<float> theta = SplitString2Float(ti_theta_text->Value());
+
+	TiXmlElement* ti_phi = ti_nf2ff->FirstChildElement("phi");
+	if (ti_phi==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read phi values ... " << endl;
+		return false;
+	}
+	TiXmlNode* ti_phi_node = ti_phi->FirstChild();
+	if (ti_phi_node==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read phi text child ... " << endl;
+		return false;
+	}
+	TiXmlText* ti_phi_text = ti_phi_node->ToText();
+	if (ti_phi_text==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLNode: Can't read phi text values ... " << endl;
+		return false;
+	}
+	vector<float> phi = SplitString2Float(ti_phi_text->Value());
+
+	nf2ff* l_nf2ff = new nf2ff(freq,theta,phi,numThreads);
+	l_nf2ff->SetVerboseLevel(Verbose);
+
+	TiXmlElement* ti_Planes = ti_nf2ff->FirstChildElement();
+	string E_name;
+	string H_name;
+	while (ti_Planes!=NULL)
+	{
+		E_name = string(ti_Planes->Attribute("E_Field"));
+		H_name = string(ti_Planes->Attribute("H_Field"));
+		if ((!E_name.empty()) && (!H_name.empty()))
+		{
+			if (l_nf2ff->AnalyseFile(E_name,H_name)==false)
+			{
+				cerr << "nf2ff::AnalyseXMLNode: Error, analysing Plane ... " << endl;
+				return false;
+			}
+		}
+		else
+		{
+			cerr << "nf2ff::AnalyseXMLNode: Error, invalid plane entry ... " << endl;
+			return false;
+		}
+		ti_Planes = ti_Planes->NextSiblingElement("Planes");
+	}
+	l_nf2ff->Write2HDF5(outfile);
+	delete l_nf2ff;
+	return true;
+}
+
+bool nf2ff::AnalyseXMLFile(string filename)
+{
+	TiXmlDocument doc(filename.c_str());
+	if (!doc.LoadFile())
+	{
+		cerr << "nf2ff::AnalyseXMLFile: Error loading xml-file failed!!! File: " << filename << endl;
+		return false;
+	}
+	TiXmlElement* ti_nf2ff = doc.FirstChildElement("nf2ff");
+	if (ti_nf2ff==NULL)
+	{
+		cerr << "nf2ff::AnalyseXMLFile: Can't read nf2ff ... " << endl;
+		return false;
+	}
+
+	return AnalyseXMLNode(ti_nf2ff);
+}
+
+bool nf2ff::AnalyseFile(string E_Field_file, string H_Field_file)
+{
+	HDF5_File_Reader E_file(E_Field_file);
+	HDF5_File_Reader H_file(H_Field_file);
+
+	if (m_Verbose>0)
+		cerr << "nf2ff: Reading planes: " << E_Field_file << " & " << E_Field_file << endl;
+
+	// read E-mesh
+	float* E_lines[3]={NULL,NULL,NULL};
+	unsigned int E_numLines[3];
+	int E_meshType;
+	if (E_file.ReadMesh(E_lines, E_numLines, E_meshType) == false)
+	{
+		cerr << "nf2ff::AnalyseFile: Error reading  E-field mesh..." << endl;
+		return false;
+	}
+
+	// read H-mesh
+	float* H_lines[3]={NULL,NULL,NULL};
+	unsigned int H_numLines[3];
+	int H_meshType;
+	if (H_file.ReadMesh(H_lines, H_numLines, H_meshType) == false)
+	{
+		cerr << "nf2ff::AnalyseFile: Error reading H-Field mesh..." << endl;
+		return false;
+	}
+
+	// compare E/H meshs
+	if (E_meshType!=H_meshType)
+	{
+		cerr << "nf2ff::AnalyseFile: Error mesh types don't agree" << endl;
+		return false;
+	}
+	if ((E_numLines[0]!=H_numLines[0]) || (E_numLines[1]!=H_numLines[1]) || (E_numLines[2]!=H_numLines[2]))
+	{
+		cerr << "nf2ff::AnalyseFile: Error mesh dimensions don't agree" << endl;
+		return false;
+	}
+	for (int n=0;n<3;++n)
+		for (unsigned int m=0;m<E_numLines[n];++m)
+			if (E_lines[n][m]!=H_lines[n][m])
+			{
+				cerr << "nf2ff::AnalyseFile: Error mesh lines don't agree" << endl;
+				return false;
+			}
+
+	if (m_Verbose>1)
+		cerr << "nf2ff: Data-Size: " << E_numLines[0] << "x" << E_numLines[1] << "x"  << E_numLines[2] << endl;
+	if (m_Verbose>1)
+		cerr << "nf2ff: calculate dft..." << endl;
+
+	unsigned int data_size[4];
+	vector<complex<float>****> E_fd_data;
+	E_file.CalcFDVectorData(m_freq,E_fd_data,data_size);
+
+	vector<complex<float>****> H_fd_data;
+	H_file.CalcFDVectorData(m_freq,H_fd_data,data_size);
+
+	if (m_Verbose>0)
+		cerr << "nf2ff: Analysing far-field for " <<  m_nf2ff.size() << " frequencies.  " << endl;
+
+	for (size_t fn=0;fn<m_nf2ff.size();++fn)
+	{
+		if (m_Verbose>1)
+			cerr << "nf2ff: f = " << m_freq.at(fn) << "Hz (" << fn+1 << "/" << m_nf2ff.size() << ") ...";
+		m_nf2ff.at(fn)->AddPlane(E_lines, E_numLines, E_fd_data.at(fn), H_fd_data.at(fn));
+		if (m_Verbose>1)
+			cerr << " done." << endl;
+	}
+	for (int n=0;n<3;++n)
+	{
+		delete[] H_lines[n];
+		delete[] E_lines[n];
+	}
+	return true;
+}
+
+bool nf2ff::Write2HDF5(string filename)
+{
+	HDF5_File_Writer hdf_file(filename);
+
+	//write mesh information
+	hdf_file.SetCurrentGroup("/Mesh");
+	size_t meshsize[1]={m_numTheta};
+	if (hdf_file.WriteData(string("theta"),m_theta,1,meshsize)==false)
+		return false;
+	meshsize[0]=m_numPhi;
+	if (hdf_file.WriteData(string("phi"),m_phi,1,meshsize)==false)
+		return false;
+	meshsize[0]=1;
+	float rad[1]={m_radius};
+	if (hdf_file.WriteData(string("r"),rad,1,meshsize)==false)
+		return false;
+
+	float attr_value = 2;
+	hdf_file.WriteAtrribute("/Mesh", "MeshType", &attr_value, 1);
+
+	//write field data
+	size_t dim = 2;
+	size_t pos = 0;
+	size_t datasize[2]={m_numPhi,m_numTheta};
+	size_t size = datasize[0]*datasize[1];
+	float* buffer = new float[size];
+	complex<float>** field_data;
+	string field_names[2]={"E_theta", "E_phi"};
+	for (int n=0;n<2;++n)
+	{
+		hdf_file.SetCurrentGroup("/nf2ff/" + field_names[n] + "/FD");
+		for (size_t fn=0;fn<m_freq.size();++fn)
+		{
+			stringstream ss;
+			ss << "f" << fn;
+			pos = 0;
+			if (n==0)
+				field_data = GetETheta(fn);
+			else
+				field_data = GetEPhi(fn);
+			for (size_t j=0;j<m_numPhi;++j)
+				for (size_t i=0;i<m_numTheta;++i)
+				{
+					buffer[pos++]=real(field_data[i][j]);
+				}
+			if (hdf_file.WriteData(ss.str() + "_real",buffer,dim,datasize)==false)
+			{
+				delete[] buffer;
+				cerr << "nf2ff::Write2HDF5: Error writing field data" << endl;
+				return false;
+			}
+
+			pos = 0;
+			for (size_t j=0;j<m_numPhi;++j)
+				for (size_t i=0;i<m_numTheta;++i)
+				{
+					buffer[pos++]=imag(field_data[i][j]);
+				}
+			if (hdf_file.WriteData(ss.str() + "_imag",buffer,dim,datasize)==false)
+			{
+				delete[] buffer;
+				cerr << "nf2ff::Write2HDF5: Error writing field data" << endl;
+				return false;
+			}
+		}
+	}
+	delete[] buffer;
+
+	//write frequency attribute
+	hdf_file.WriteAtrribute("/nf2ff", "Frequency",m_freq);
+
+	buffer = new float[m_freq.size()];
+	//write radiated power attribute
+	for (size_t fn=0;fn<m_freq.size();++fn)
+		buffer[fn] = GetRadPower(fn);
+	hdf_file.WriteAtrribute("/nf2ff", "Prad",buffer,m_freq.size());
+
+	//write max directivity attribute
+	for (size_t fn=0;fn<m_freq.size();++fn)
+		buffer[fn] = GetMaxDirectivity(fn);
+	hdf_file.WriteAtrribute("/nf2ff", "Dmax",buffer,m_freq.size());
+
+	delete[] buffer;
+	return true;
+}
diff --git a/nf2ff/nf2ff.h b/nf2ff/nf2ff.h
new file mode 100644
index 0000000..60b10de
--- /dev/null
+++ b/nf2ff/nf2ff.h
@@ -0,0 +1,65 @@
+/*
+*	Copyright (C) 2012 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 NF2FF_H
+#define NF2FF_H
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <vector>
+#include <cmath>
+#include <complex>
+#include "nf2ff_calc.h"
+
+using namespace std;
+
+class TiXmlElement;
+
+class nf2ff
+{
+public:
+	nf2ff(vector<float> freq, vector<float> theta, vector<float> phi, unsigned int numThreads=0);
+	~nf2ff();
+
+	bool AnalyseFile(string E_Field_file, string H_Field_file);
+
+	float GetRadPower(size_t f_idx) const {return m_nf2ff.at(f_idx)->GetRadPower();}
+	float GetMaxDirectivity(size_t f_idx) const {return m_nf2ff.at(f_idx)->GetMaxDirectivity();}
+
+	complex<float>** GetETheta(size_t f_idx) const {return m_nf2ff.at(f_idx)->GetETheta();}
+	complex<float>** GetEPhi(size_t f_idx) const {return m_nf2ff.at(f_idx)->GetEPhi();}
+
+	//! Write results to a hdf5 file
+	bool Write2HDF5(string filename);
+
+	void SetVerboseLevel(int level) {m_Verbose=level;}
+
+	static bool AnalyseXMLNode(TiXmlElement* ti_nf2ff);
+	static bool AnalyseXMLFile(string filename);
+
+protected:
+	vector<float> m_freq;
+	unsigned int m_numTheta;
+	unsigned int m_numPhi;
+	float* m_theta;
+	float* m_phi;
+	float m_radius;
+	int m_Verbose;
+	vector<nf2ff_calc*> m_nf2ff;
+};
+
+#endif // NF2FF_H
diff --git a/nf2ff/nf2ff.pro b/nf2ff/nf2ff.pro
new file mode 100644
index 0000000..f8b6dc7
--- /dev/null
+++ b/nf2ff/nf2ff.pro
@@ -0,0 +1,60 @@
+TARGET = nf2ff
+CONFIG += console
+CONFIG -= app_bundle qt
+TEMPLATE = app
+OBJECTS_DIR = obj
+INCLUDEPATH += .
+INCLUDEPATH += ../../tinyxml
+CONFIG += debug_and_release
+
+win32 {
+    INCLUDEPATH += ../../hdf5/include ../../hdf5/include/cpp ../../boost/include/boost-1_42
+	LIBS +=  ../../hdf5/lib/hdf5.lib
+    LIBS += ../../boost/lib/libboost_thread-mgw44-mt.lib
+	LIBS += ../../tinyxml/release/libtinyxml2.a
+}
+!win32 {
+    LIBS += -lboost_thread
+	LIBS += -lhdf5
+	LIBS += ../../tinyxml/libtinyxml.so
+}
+QMAKE_LFLAGS += \'-Wl,-rpath,\$$ORIGIN/../../tinyxml\'
+
+TOOLSPATH = ../tools
+
+#### SOURCES ################################################################
+SOURCES += main.cpp \
+	nf2ff.cpp \
+	nf2ff_calc.cpp
+
+# tools
+ SOURCES += $$TOOLSPATH/global.cpp \
+	$$TOOLSPATH/useful.cpp \
+	$$TOOLSPATH/array_ops.cpp \
+	$$TOOLSPATH/hdf5_file_reader.cpp \
+	$$TOOLSPATH/hdf5_file_writer.cpp
+
+#### HEADERS ################################################################
+HEADERS += nf2ff.h \
+	nf2ff_calc.h
+
+# tools
+HEADERS += $$TOOLSPATH/constants.h \
+	$$TOOLSPATH/array_ops.h \
+	$$TOOLSPATH/global.h \
+	$$TOOLSPATH/useful.h \
+	$$TOOLSPATH/aligned_allocator.h \
+	$$TOOLSPATH/hdf5_file_reader.h \
+	$$TOOLSPATH/hdf5_file_writer.h
+
+QMAKE_CXXFLAGS_RELEASE = -O3 \
+	-g \
+	-march=native
+QMAKE_CXXFLAGS_DEBUG = -O0 \
+	-g \
+	-march=native
+
+# add git revision
+# QMAKE_CXXFLAGS += -DGIT_VERSION=\\\"`git describe --tags`\\\"
+
+
diff --git a/nf2ff/nf2ff_calc.cpp b/nf2ff/nf2ff_calc.cpp
new file mode 100644
index 0000000..3f04f11
--- /dev/null
+++ b/nf2ff/nf2ff_calc.cpp
@@ -0,0 +1,353 @@
+/*
+*	Copyright (C) 2012 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 "nf2ff_calc.h"
+#include "../tools/array_ops.h"
+#include "../tools/useful.h"
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <vector>
+#include <cmath>
+#include <complex>
+#include <iostream>
+#include <sstream>
+
+nf2ff_calc_thread::nf2ff_calc_thread(nf2ff_calc* nfc, unsigned int start, unsigned int stop, unsigned int threadID, nf2ff_data &data)
+{
+	m_nf_calc = nfc;
+	m_start = start;
+	m_stop = stop;
+	m_threadID = threadID;
+	m_data = data;
+}
+
+void nf2ff_calc_thread::operator()()
+{
+	m_nf_calc->m_Barrier->wait(); // start
+
+	int ny = m_data.ny;
+	int nP = (ny+1)%3;
+	int nPP = (ny+2)%3;
+
+	unsigned int* numLines = m_data.numLines;
+	float* normDir = m_data.normDir;
+	float **lines = m_data.lines;
+	float* edge_length_P = m_data.edge_length_P;
+	float* edge_length_PP = m_data.edge_length_PP;
+
+	unsigned int pos[3];
+	unsigned int pos_t=0;
+	unsigned int num_t=m_stop-m_start+1;
+
+
+	complex<float>**** Js=m_data.Js;
+	complex<float>**** Ms=m_data.Ms;
+	complex<float>**** E_field=m_data.E_field;
+	complex<float>**** H_field=m_data.H_field;
+
+	// calc Js and Ms (eq. 8.15a/b)
+	pos[ny]=0;
+	for (pos_t=0; pos_t<num_t; ++pos_t)
+	{
+		pos[nP] = m_start+pos_t;
+		for (pos[nPP]=0; pos[nPP]<numLines[nPP]; ++pos[nPP])
+			{
+				// Js =  n x H
+				Js[0][pos[0]][pos[1]][pos[2]] = normDir[1]*H_field[2][pos[0]][pos[1]][pos[2]] - normDir[2]*H_field[1][pos[0]][pos[1]][pos[2]];
+				Js[1][pos[0]][pos[1]][pos[2]] = normDir[2]*H_field[0][pos[0]][pos[1]][pos[2]] - normDir[0]*H_field[2][pos[0]][pos[1]][pos[2]];
+				Js[2][pos[0]][pos[1]][pos[2]] = normDir[0]*H_field[1][pos[0]][pos[1]][pos[2]] - normDir[1]*H_field[0][pos[0]][pos[1]][pos[2]];
+
+				// Ms = -n x E
+				Ms[0][pos[0]][pos[1]][pos[2]] = normDir[2]*E_field[1][pos[0]][pos[1]][pos[2]] - normDir[1]*E_field[2][pos[0]][pos[1]][pos[2]];
+				Ms[1][pos[0]][pos[1]][pos[2]] = normDir[0]*E_field[2][pos[0]][pos[1]][pos[2]] - normDir[2]*E_field[0][pos[0]][pos[1]][pos[2]];
+				Ms[2][pos[0]][pos[1]][pos[2]] = normDir[1]*E_field[0][pos[0]][pos[1]][pos[2]] - normDir[0]*E_field[1][pos[0]][pos[1]][pos[2]];
+			}
+	}
+
+	complex<float>** m_Nt=m_data.m_Nt;
+	complex<float>** m_Np=m_data.m_Np;
+	complex<float>** m_Lt=m_data.m_Lt;
+	complex<float>** m_Lp=m_data.m_Lp;
+
+	// calc local Nt,Np,Lt and Lp
+	float area;
+	float cosT_cosP,cosP_sinT;
+	float cosT_sinP,sinT_sinP;
+	float sinT,sinP;
+	float cosP,cosT;
+	float r_cos_psi;
+	float k = 2*M_PI*m_nf_calc->m_freq/__C0__;
+	complex<float> exp_jkr;
+	complex<float> _I_(0,1);
+	for (unsigned int tn=0;tn<m_nf_calc->m_numTheta;++tn)
+		for (unsigned int pn=0;pn<m_nf_calc->m_numPhi;++pn)
+		{
+			sinT = sin(m_nf_calc->m_theta[tn]);
+			sinP = sin(m_nf_calc->m_phi[pn]);
+			cosT = cos(m_nf_calc->m_theta[tn]);
+			cosP = cos(m_nf_calc->m_phi[pn]);
+			cosT_cosP = cosT*cosP;
+			cosT_sinP = cosT*sinP;
+			cosP_sinT = cosP*sinT;
+			sinT_sinP = sinP*sinT;
+
+			for (pos_t=0; pos_t<num_t; ++pos_t)
+			{
+				pos[nP] = m_start+pos_t;
+				for (pos[nPP]=0; pos[nPP]<numLines[nPP]; ++pos[nPP])
+					{
+						r_cos_psi = lines[0][pos[0]]*cosP_sinT + lines[1][pos[1]]*sinT_sinP + lines[2][pos[2]]*cosT;
+						exp_jkr = exp(_I_*k*r_cos_psi);
+						area = edge_length_P[pos[nP]]*edge_length_PP[pos[nPP]];
+						m_Nt[tn][pn] += area*exp_jkr*(Js[0][pos[0]][pos[1]][pos[2]]*cosT_cosP + Js[1][pos[0]][pos[1]][pos[2]]*cosT_sinP \
+													- Js[2][pos[0]][pos[1]][pos[2]]*sinT);
+						m_Np[tn][pn] += area*exp_jkr*(Js[1][pos[0]][pos[1]][pos[2]]*cosP - Js[0][pos[0]][pos[1]][pos[2]]*sinP);
+
+						m_Lt[tn][pn] += area*exp_jkr*(Ms[0][pos[0]][pos[1]][pos[2]]*cosT_cosP + Ms[1][pos[0]][pos[1]][pos[2]]*cosT_sinP \
+													- Ms[2][pos[0]][pos[1]][pos[2]]*sinT);
+						m_Lp[tn][pn] += area*exp_jkr*(Ms[1][pos[0]][pos[1]][pos[2]]*cosP - Ms[0][pos[0]][pos[1]][pos[2]]*sinP);
+					}
+
+			}
+		}
+
+	m_nf_calc->m_Barrier->wait(); //combine all thread local Nt,Np,Lt and Lp
+
+	m_nf_calc->m_Barrier->wait(); //wait for termination
+}
+
+
+/***********************************************************************/
+
+
+nf2ff_calc::nf2ff_calc(float freq, vector<float> theta, vector<float> phi)
+{
+	m_freq = freq;
+
+	m_numTheta = theta.size();
+	m_theta = new float[m_numTheta];
+	for (size_t n=0;n<m_numTheta;++n)
+		 m_theta[n]=theta.at(n);
+
+	m_numPhi = phi.size();
+	m_phi = new float[m_numPhi];
+	for (size_t n=0;n<m_numPhi;++n)
+		 m_phi[n]=phi.at(n);
+
+	unsigned int numLines[2] = {m_numTheta, m_numPhi};
+	m_E_theta = Create2DArray<std::complex<float> >(numLines);
+	m_E_phi = Create2DArray<std::complex<float> >(numLines);
+	m_H_theta = Create2DArray<std::complex<float> >(numLines);
+	m_H_phi = Create2DArray<std::complex<float> >(numLines);
+	m_P_rad = Create2DArray<float>(numLines);
+
+	m_centerCoord[0]=m_centerCoord[1]=m_centerCoord[2]=0;
+
+	m_radPower = 0;
+	m_maxDir = 0;
+	m_radius = 1;
+
+	m_Barrier = NULL;
+	m_numThreads = boost::thread::hardware_concurrency();
+}
+
+nf2ff_calc::~nf2ff_calc()
+{
+	delete[] m_phi;
+	m_phi = NULL;
+	delete[] m_theta;
+	m_theta = NULL;
+
+	unsigned int numLines[2] = {m_numTheta, m_numPhi};
+	Delete2DArray(m_E_theta,numLines);
+	m_E_theta = NULL;
+	Delete2DArray(m_E_phi,numLines);
+	m_E_phi = NULL;
+	Delete2DArray(m_H_theta,numLines);
+	m_H_theta = NULL;
+	Delete2DArray(m_H_phi,numLines);
+	m_H_phi = NULL;
+	Delete2DArray(m_P_rad,numLines);
+	m_P_rad = NULL;
+
+	delete m_Barrier;
+	m_Barrier = NULL;
+}
+
+bool nf2ff_calc::AddPlane(float **lines, unsigned int* numLines, complex<float>**** E_field, complex<float>**** H_field)
+{
+	//find normal direction
+	int ny = -1;
+	int nP,nPP;
+	for (int n=0;n<3;++n)
+	{
+		nP = (n+1)%3;
+		nPP = (n+2)%3;
+		if ((numLines[n]==1) && (numLines[nP]>2) && (numLines[nPP]>2))
+			ny=n;
+	}
+	nP = (ny+1)%3;
+	nPP = (ny+2)%3;
+	if (ny<0)
+	{
+		cerr << "nf2ff_calc::AddPlane: Error can't determine normal direction..." << endl;
+		return false;
+	}
+
+	complex<float>**** Js = Create_N_3DArray<complex<float> >(numLines);
+	complex<float>**** Ms = Create_N_3DArray<complex<float> >(numLines);
+
+	float normDir[3]= {0,0,0};
+	if (lines[ny][0]>=m_centerCoord[ny])
+		normDir[ny]=1;
+	else
+		normDir[ny]=-1;
+	unsigned int pos[3];
+
+	float edge_length_P[numLines[nP]];
+	for (unsigned int n=1;n<numLines[nP]-1;++n)
+		edge_length_P[n]=0.5*(lines[nP][n+1]-lines[nP][n-1]);
+	edge_length_P[0]=0.5*(lines[nP][1]-lines[nP][0]);
+	edge_length_P[numLines[nP]-1]=0.5*(lines[nP][numLines[nP]-1]-lines[nP][numLines[nP]-2]);
+
+	float edge_length_PP[numLines[nPP]];
+	for (unsigned int n=1;n<numLines[nPP]-1;++n)
+		edge_length_PP[n]=0.5*(lines[nPP][n+1]-lines[nPP][n-1]);
+	edge_length_PP[0]=0.5*(lines[nPP][1]-lines[nPP][0]);
+	edge_length_PP[numLines[nPP]-1]=0.5*(lines[nPP][numLines[nPP]-1]-lines[nPP][numLines[nPP]-2]);
+
+	complex<float> power = 0;
+	float area;
+	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])
+			{
+				area = edge_length_P[pos[nP]]*edge_length_PP[pos[nPP]];
+				power = (E_field[nP][pos[0]][pos[1]][pos[2]]*conj(H_field[nPP][pos[0]][pos[1]][pos[2]]) \
+							- E_field[nPP][pos[0]][pos[1]][pos[2]]*conj(H_field[nP][pos[0]][pos[1]][pos[2]]));
+				m_radPower += 0.5*area*real(power)*normDir[ny];
+			}
+
+	unsigned int numAngles[2] = {m_numTheta, m_numPhi};
+
+	// setup multi-threading jobs
+	vector<unsigned int> jpt = AssignJobs2Threads(numLines[nP], m_numThreads, true);
+	m_numThreads = jpt.size();
+	nf2ff_data thread_data[m_numThreads];
+	m_Barrier = new boost::barrier(m_numThreads+1); // numThread workers + 1 controller
+	unsigned int start=0;
+	unsigned int stop=jpt.at(0)-1;
+	for (unsigned int n=0; n<m_numThreads; n++)
+	{
+		thread_data[n].ny=ny;
+		thread_data[n].normDir=normDir;
+		thread_data[n].numLines=numLines;
+		thread_data[n].lines=lines;
+		thread_data[n].edge_length_P=edge_length_P;
+		thread_data[n].edge_length_PP=edge_length_PP;
+		thread_data[n].E_field=E_field;
+		thread_data[n].H_field=H_field;
+		thread_data[n].Js=Js;
+		thread_data[n].Ms=Ms;
+		thread_data[n].m_Nt=Create2DArray<complex<float> >(numAngles);
+		thread_data[n].m_Np=Create2DArray<complex<float> >(numAngles);
+		thread_data[n].m_Lt=Create2DArray<complex<float> >(numAngles);
+		thread_data[n].m_Lp=Create2DArray<complex<float> >(numAngles);
+
+		boost::thread *t = new boost::thread( nf2ff_calc_thread(this,start,stop,n,thread_data[n]) );
+
+		m_thread_group.add_thread( t );
+
+		start = stop+1;
+		if (n<m_numThreads-1)
+			stop = start + jpt.at(n+1)-1;
+	}
+	//all threads a running and waiting for the barrier
+
+	m_Barrier->wait(); //start
+
+	// threads: calc Js and Ms (eq. 8.15a/b)
+	// threads calc their local Nt,Np,Lt and Lp
+
+	m_Barrier->wait(); //combine all thread local Nt,Np,Lt and Lp
+
+	//cleanup E- & H-Fields
+	Delete_N_3DArray(E_field,numLines);
+	Delete_N_3DArray(H_field,numLines);
+
+	complex<float>** Nt = Create2DArray<complex<float> >(numAngles);
+	complex<float>** Np = Create2DArray<complex<float> >(numAngles);
+	complex<float>** Lt = Create2DArray<complex<float> >(numAngles);
+	complex<float>** Lp = Create2DArray<complex<float> >(numAngles);
+
+	for (unsigned int n=0; n<m_numThreads; n++)
+	{
+		for (unsigned int tn=0;tn<m_numTheta;++tn)
+			for (unsigned int pn=0;pn<m_numPhi;++pn)
+			{
+				Nt[tn][pn] += thread_data[n].m_Nt[tn][pn];
+				Np[tn][pn] += thread_data[n].m_Np[tn][pn];
+				Lt[tn][pn] += thread_data[n].m_Lt[tn][pn];
+				Lp[tn][pn] += thread_data[n].m_Lp[tn][pn];
+			}
+		Delete2DArray(thread_data[n].m_Nt,numAngles);
+		Delete2DArray(thread_data[n].m_Np,numAngles);
+		Delete2DArray(thread_data[n].m_Lt,numAngles);
+		Delete2DArray(thread_data[n].m_Lp,numAngles);
+	}
+
+	m_Barrier->wait(); //wait for termination
+	m_thread_group.join_all(); // wait for termination
+	delete m_Barrier;
+	m_Barrier = NULL;
+
+	//cleanup Js & Ms
+	Delete_N_3DArray(Js,numLines);
+	Delete_N_3DArray(Ms,numLines);
+
+	// calc equations 8.23a/b and 8.24a/b
+	float k = 2*M_PI*m_freq/__C0__;
+	complex<float> factor(0,k/4.0/M_PI/m_radius);
+	complex<float> f_exp(0,-1*k*m_radius);
+	factor *= exp(f_exp);
+	complex<float> Z0 = __Z0__;
+	float P_max = 0;
+	for (unsigned int tn=0;tn<m_numTheta;++tn)
+		for (unsigned int pn=0;pn<m_numPhi;++pn)
+		{
+			m_E_theta[tn][pn] -= factor*(Lp[tn][pn] + Z0*Nt[tn][pn]);
+			m_E_phi[tn][pn] += factor*(Lt[tn][pn] - Z0*Np[tn][pn]);
+
+			m_H_theta[tn][pn] += factor*(Np[tn][pn] - Lt[tn][pn]/Z0);
+			m_H_phi[tn][pn] -= factor*(Nt[tn][pn] + Lp[tn][pn]/Z0);
+
+			m_P_rad[tn][pn] = m_radius*m_radius/(2*__Z0__) * abs((m_E_theta[tn][pn]*conj(m_E_theta[tn][pn])+m_E_phi[tn][pn]*conj(m_E_phi[tn][pn])));
+			if (m_P_rad[tn][pn]>P_max)
+				P_max = m_P_rad[tn][pn];
+		}
+
+	//cleanup Nx and Lx
+	Delete2DArray(Nt,numAngles);
+	Delete2DArray(Np,numAngles);
+	Delete2DArray(Lt,numAngles);
+	Delete2DArray(Lp,numAngles);
+
+	m_maxDir = 4*M_PI*P_max / m_radPower;
+
+	return true;
+}
diff --git a/nf2ff/nf2ff_calc.h b/nf2ff/nf2ff_calc.h
new file mode 100644
index 0000000..fbfbb6a
--- /dev/null
+++ b/nf2ff/nf2ff_calc.h
@@ -0,0 +1,114 @@
+/*
+*	Copyright (C) 2012 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 NF2FF_CALC_H
+#define NF2FF_CALC_H
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <vector>
+#include <cmath>
+#include <complex>
+#include <boost/thread.hpp>
+
+using namespace std;
+
+class nf2ff_calc;
+
+// data structure to exchange data between thread-controller and worker-threads
+typedef struct
+{
+	//local working data IN
+	int ny;
+	float* normDir;
+	unsigned int* numLines;
+	float **lines;
+	float* edge_length_P;
+	float* edge_length_PP;
+
+	complex<float>**** E_field;
+	complex<float>**** H_field;
+	complex<float>**** Js;
+	complex<float>**** Ms;
+
+	//local working data OUT
+	complex<float>** m_Nt;
+	complex<float>** m_Np;
+	complex<float>** m_Lt;
+	complex<float>** m_Lp;
+
+} nf2ff_data;
+
+class nf2ff_calc_thread
+{
+public:
+	nf2ff_calc_thread(nf2ff_calc* nfc, unsigned int start, unsigned int stop, unsigned int threadID, nf2ff_data &data);
+	void operator()();
+
+protected:
+	unsigned int m_start, m_stop, m_threadID;
+	nf2ff_calc *m_nf_calc;
+
+	nf2ff_data m_data;
+};
+
+class nf2ff_calc
+{
+	// allow full data access to nf2ff_calc_thread class
+	friend class nf2ff_calc_thread;
+public:
+	nf2ff_calc(float freq, vector<float> theta, vector<float> phi);
+	~nf2ff_calc();
+
+	float GetRadPower() const {return m_radPower;}
+	float GetMaxDirectivity() const {return m_maxDir;}
+
+	complex<float>** GetETheta() const {return m_E_theta;}
+	complex<float>** GetEPhi() const {return m_E_phi;}
+
+	unsigned int GetNumThreads() const {return m_numThreads;}
+	void SetNumThreads(unsigned int n) {m_numThreads=n;}
+
+	bool AddPlane(float **lines, unsigned int* numLines, complex<float>**** E_field, complex<float>**** H_field);
+
+protected:
+	float m_freq;
+	float m_radius;
+
+	float m_radPower;
+	float m_maxDir;
+
+	complex<float>** m_E_theta;
+	complex<float>** m_E_phi;
+	complex<float>** m_H_theta;
+	complex<float>** m_H_phi;
+	float** m_P_rad;
+
+	float m_centerCoord[3];
+	unsigned int m_numTheta;
+	unsigned int m_numPhi;
+	float* m_theta;
+	float* m_phi;
+
+	//boost multi-threading
+	unsigned int m_numThreads;
+	boost::thread_group m_thread_group;
+	boost::barrier *m_Barrier;
+};
+
+
+#endif // NF2FF_CALC_H