diff --git a/FDTD/engine_interface_cylindrical_fdtd.cpp b/FDTD/engine_interface_cylindrical_fdtd.cpp
index df73e3f..d19df87 100644
--- a/FDTD/engine_interface_cylindrical_fdtd.cpp
+++ b/FDTD/engine_interface_cylindrical_fdtd.cpp
@@ -37,10 +37,30 @@ double* Engine_Interface_Cylindrical_FDTD::GetHField(const unsigned int* pos, do
 		return Engine_Interface_FDTD::GetHField(pos, out);
 
 	unsigned int iPos[] = {pos[0],pos[1],pos[2]};
+	if (iPos[1]==0)
+		iPos[1]=m_Op->GetNumberOfLines(1)-1;
 
-	if (pos[1]==m_Op->GetNumberOfLines(1)-1)
-		iPos[1]=0;
-	return Engine_Interface_FDTD::GetHField(iPos, out);
+	int nP,nPP;
+	for (int n=0; n<3; ++n)
+	{
+		nP = (n+1)%3;
+		nPP = (n+2)%3;
+		if ((iPos[0]==m_Op->GetNumberOfLines(0)-1) || (iPos[2]==m_Op->GetNumberOfLines(2)-1) || (iPos[nP]==0) || (iPos[nPP]==0))
+		{
+			out[n] = 0;
+			continue;
+		}
+		out[n]=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
+		--iPos[nP];
+		out[n]+=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
+		--iPos[nPP];
+		out[n]+=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
+		++iPos[nP];
+		out[n]+=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
+		++iPos[nPP];
+		out[n]/=4;
+	}
+	return out;
 }
 
 double* Engine_Interface_Cylindrical_FDTD::GetRawInterpolatedField(const unsigned int* pos, double* out, int type) const
diff --git a/FDTD/operator.cpp b/FDTD/operator.cpp
index 7064a2c..3fae3c8 100644
--- a/FDTD/operator.cpp
+++ b/FDTD/operator.cpp
@@ -205,49 +205,44 @@ double Operator::GetNodeArea(int ny, const int pos[3], bool dualMesh) const
 	return GetNodeArea(ny, uiPos, dualMesh);
 }
 
-bool Operator::SnapToMesh(const double* dcoord, unsigned int* uicoord, bool lower, bool* inside) const
+unsigned int Operator::SnapToMeshLine(int ny, double coord, bool &inside, bool dualMesh) const
+{
+	inside = false;
+	if ((ny<0) || (ny>2))
+		return 0;
+	if (coord<GetDiscLine(ny,0))
+		return 0;
+	unsigned int numLines = GetNumberOfLines(ny);
+	if (coord>GetDiscLine(ny,numLines-1))
+		return numLines-1;
+	inside=true;
+	if (dualMesh==false)
+	{
+		for (unsigned int n=0;n<numLines;++n)
+		{
+			if (coord<=GetDiscLine(ny,n,true))
+				return n;
+		}
+	}
+	else
+	{
+		for (unsigned int n=1;n<numLines;++n)
+		{
+			if (coord<=GetDiscLine(ny,n,false))
+				return n-1;
+		}
+	}
+	//should not happen
+	return 0;
+}
+
+bool Operator::SnapToMesh(const double* dcoord, unsigned int* uicoord, bool dualMesh, bool* inside) const
 {
 	bool ok=true;
-	unsigned int numLines[3];
 	for (int n=0; n<3; ++n)
 	{
-		numLines[n] = GetNumberOfLines(n);
-		if (inside) //set defaults
-			inside[n] = true;
-		uicoord[n]=0;
-		if (dcoord[n]<discLines[n][0])
-		{
-			ok=false;
-			uicoord[n]=0;
-			if (inside) inside[n] = false;
-		}
-		else if (dcoord[n]==discLines[n][0])
-			uicoord[n]=0;
-		else if (dcoord[n]>discLines[n][numLines[n]-1])
-		{
-			ok=false;
-			uicoord[n]=numLines[n]-1;
-			if (lower) uicoord[n]=numLines[n]-2;
-			if (inside) inside[n] = false;
-		}
-		else if (dcoord[n]==discLines[n][numLines[n]-1])
-		{
-			uicoord[n]=numLines[n]-1;
-			if (lower) uicoord[n]=numLines[n]-2;
-		}
-		else
-			for (unsigned int i=1; i<numLines[n]; ++i)
-			{
-				if (dcoord[n]<discLines[n][i])
-				{
-					if (fabs(dcoord[n]-discLines[n][i])<(fabs(dcoord[n]-discLines[n][i-1])))
-						uicoord[n]=i;
-					else
-						uicoord[n]=i-1;
-					if (lower) uicoord[n]=i-1;
-					i = numLines[n];
-				}
-			}
+		uicoord[n] = SnapToMeshLine(n,dcoord[n],inside[n],dualMesh);
+		ok &= inside[n];
 	}
 //	cerr << "Operator::SnapToMesh Wish: " << dcoord[0] << " " << dcoord[1] << " " << dcoord[2] << endl;
 //	cerr << "Operator::SnapToMesh Found: " << discLines[0][uicoord[0]] << " " << discLines[1][uicoord[1]] << " " << discLines[2][uicoord[2]] << endl;
diff --git a/FDTD/operator.h b/FDTD/operator.h
index 7021818..676b6ef 100644
--- a/FDTD/operator.h
+++ b/FDTD/operator.h
@@ -118,7 +118,9 @@ public:
 	*/
 	virtual double GetEdgeArea(int ny, const unsigned int pos[3], bool dualMesh = false) const {return GetNodeArea(ny,pos,dualMesh);}
 
-	virtual bool SnapToMesh(const double* coord, unsigned int* uicoord, bool lower=false, bool* inside=NULL) const;
+	virtual unsigned int SnapToMeshLine(int ny, double coord, bool &inside, bool dualMesh=false) const;
+
+	virtual bool SnapToMesh(const double* coord, unsigned int* uicoord, bool dualMesh=false, bool* inside=NULL) const;
 
 	virtual void AddExtension(Operator_Extension* op_ext);
 	virtual size_t GetNumberOfExtentions() const {return m_Op_exts.size();}