diff --git a/FDTD/engine.h b/FDTD/engine.h
index fb7f19e..90d58be 100644
--- a/FDTD/engine.h
+++ b/FDTD/engine.h
@@ -34,8 +34,8 @@ public:
 
 	virtual unsigned int GetNumberOfTimesteps() {return numTS;};
 
-	virtual FDTD_FLOAT**** GetVoltages() {return volt;};
-	virtual FDTD_FLOAT**** GetCurrents() {return curr;};
+	inline virtual FDTD_FLOAT GetVolt( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return volt[n][x][y][z]; }
+	inline virtual FDTD_FLOAT GetCurr( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return curr[n][x][y][z]; }
 
 protected:
 	Engine(const Operator* op);
diff --git a/FDTD/processcurrent.cpp b/FDTD/processcurrent.cpp
index eec845c..8176c9d 100644
--- a/FDTD/processcurrent.cpp
+++ b/FDTD/processcurrent.cpp
@@ -34,7 +34,6 @@ void ProcessCurrent::DefineStartStopCoord(double* dstart, double* dstop)
 
 int ProcessCurrent::Process()
 {
-	FDTD_FLOAT**** curr = Eng->GetCurrents();
 	if (Enabled==false) return -1;
 	if (CheckTimestep()==false) return GetNextInterval();
 	FDTD_FLOAT current=0;
@@ -56,27 +55,27 @@ int ProcessCurrent::Process()
 	//x-current
 	if (m_start_inside[1] && m_start_inside[2])
 		for (unsigned int i=start[0];i<stop[0];++i)
-			current+=curr[0][i][start[1]][start[2]];
+			current+=Eng->GetCurr(0,i,start[1],start[2]);
 	//y-current
 	if (m_stop_inside[0] && m_start_inside[2])
 		for (unsigned int i=start[1];i<stop[1];++i)
-			current+=curr[1][stop[0]][i][start[2]];
+			current+=Eng->GetCurr(1,stop[0],i,start[2]);
 	//z-current
 	if (m_stop_inside[0] && m_stop_inside[1])
 		for (unsigned int i=start[2];i<stop[2];++i)
-			current+=curr[2][stop[0]][stop[1]][i];
+			current+=Eng->GetCurr(2,stop[0],stop[1],i);
 	//x-current
 	if (m_stop_inside[1] && m_stop_inside[2])
 		for (unsigned int i=start[0];i<stop[0];++i)
-			current-=curr[0][i][stop[1]][stop[2]];
+			current-=Eng->GetCurr(0,i,stop[1],stop[2]);
 	//y-current
 	if (m_start_inside[0] && m_stop_inside[2])
 		for (unsigned int i=start[1];i<stop[1];++i)
-			current-=curr[1][start[0]][i][stop[2]];
+			current-=Eng->GetCurr(1,start[0],i,stop[2]);
 	//z-current
 	if (m_start_inside[0] && m_start_inside[1])
 		for (unsigned int i=start[2];i<stop[2];++i)
-			current-=curr[2][start[0]][start[1]][i];
+			current-=Eng->GetCurr(2,start[0],start[1],i);
 
 //	cerr << "ts: " << Eng->numTS << " i: " << current << endl;
 	v_current.push_back(current);
diff --git a/FDTD/processfields.cpp b/FDTD/processfields.cpp
index 75a8257..4ec4fad 100644
--- a/FDTD/processfields.cpp
+++ b/FDTD/processfields.cpp
@@ -167,13 +167,12 @@ void ProcessFields::DefineStartStopCoord(double* dstart, double* dstop)
 	}
 }
 
-
-double ProcessFields::CalcTotalEnergy()
+double ProcessFields::CalcTotalEnergy() const
 {
-	FDTD_FLOAT**** volt = Eng->GetVoltages();
-	FDTD_FLOAT**** curr = Eng->GetCurrents();
+	if (!Eng)
+		return 0;
+
 	double energy=0;
-	if (Eng==NULL) return 0.0;
 	unsigned int pos[3];
 	for (pos[0]=0;pos[0]<Op->GetNumberOfLines(0);++pos[0])
 	{
@@ -181,12 +180,12 @@ double ProcessFields::CalcTotalEnergy()
 		{
 			for (pos[2]=0;pos[2]<Op->GetNumberOfLines(2);++pos[2])
 			{
-				energy+=fabs(volt[0][pos[0]][pos[1]][pos[2]] * curr[1][pos[0]][pos[1]][pos[2]]);
-				energy+=fabs(volt[0][pos[0]][pos[1]][pos[2]] * curr[2][pos[0]][pos[1]][pos[2]]);
-				energy+=fabs(volt[1][pos[0]][pos[1]][pos[2]] * curr[0][pos[0]][pos[1]][pos[2]]);
-				energy+=fabs(volt[1][pos[0]][pos[1]][pos[2]] * curr[2][pos[0]][pos[1]][pos[2]]);
-				energy+=fabs(volt[2][pos[0]][pos[1]][pos[2]] * curr[0][pos[0]][pos[1]][pos[2]]);
-				energy+=fabs(volt[2][pos[0]][pos[1]][pos[2]] * curr[1][pos[0]][pos[1]][pos[2]]);
+				energy+=fabs(Eng->GetVolt(0,pos[0],pos[1],pos[2]) * Eng->GetCurr(1,pos[0],pos[1],pos[2]));
+				energy+=fabs(Eng->GetVolt(0,pos[0],pos[1],pos[2]) * Eng->GetCurr(2,pos[0],pos[1],pos[2]));
+				energy+=fabs(Eng->GetVolt(1,pos[0],pos[1],pos[2]) * Eng->GetCurr(0,pos[0],pos[1],pos[2]));
+				energy+=fabs(Eng->GetVolt(1,pos[0],pos[1],pos[2]) * Eng->GetCurr(2,pos[0],pos[1],pos[2]));
+				energy+=fabs(Eng->GetVolt(2,pos[0],pos[1],pos[2]) * Eng->GetCurr(0,pos[0],pos[1],pos[2]));
+				energy+=fabs(Eng->GetVolt(2,pos[0],pos[1],pos[2]) * Eng->GetCurr(1,pos[0],pos[1],pos[2]));
 			}
 		}
 	}
@@ -205,7 +204,7 @@ void  ProcessFields::SetSubSampling(unsigned int subSampleRate, int dir)
 	else subSample[dir]=subSampleRate;
 }
 
-void ProcessFields::WriteVTKHeader(ofstream &file, double** discLines, unsigned int* numLines, unsigned int precision)
+void ProcessFields::WriteVTKHeader(ofstream &file, double const* const* discLines, unsigned int const* numLines, unsigned int precision)
 {
 	file << "# vtk DataFile Version 2.0" << endl;
 	file << "Rectilinear Grid openEMS_ProcessFields" << endl;
@@ -227,7 +226,7 @@ void ProcessFields::WriteVTKHeader(ofstream &file, double** discLines, unsigned
 	file << "POINT_DATA " << numLines[0]*numLines[1]*numLines[2] << endl;
 }
 
-void ProcessFields::WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT**** array, unsigned int* numLines, unsigned int precision)
+void ProcessFields::WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines, unsigned int precision)
 {
 	file << "VECTORS " << name << " float " << endl;
 
@@ -250,14 +249,14 @@ void ProcessFields::WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT*
 }
 
 
-bool ProcessFields::DumpVectorArray2VTK(ofstream &file, string name, FDTD_FLOAT**** array, double** discLines, unsigned int* numLines, unsigned int precision)
+bool ProcessFields::DumpVectorArray2VTK(ofstream &file, string name, FDTD_FLOAT const* const* const* const* array, double const* const* discLines, unsigned int const* numLines, unsigned int precision)
 {
 	WriteVTKHeader(file, discLines, numLines, precision);
 	WriteVTKVectorArray(file, name, array, numLines, precision);
 	return true;
 }
 
-bool ProcessFields::DumpMultiVectorArray2VTK(ofstream &file, string names[], FDTD_FLOAT**** array[], unsigned int numFields, double** discLines, unsigned int* numLines, unsigned int precision)
+bool ProcessFields::DumpMultiVectorArray2VTK(ofstream &file, string names[], FDTD_FLOAT const* const* const* const* const* array, unsigned int numFields, double const* const* discLines, unsigned int const* numLines, unsigned int precision)
 {
 	WriteVTKHeader(file, discLines, numLines, precision);
 	for (unsigned int n=0;n<numFields;++n)
@@ -268,7 +267,7 @@ bool ProcessFields::DumpMultiVectorArray2VTK(ofstream &file, string names[], FDT
 	return true;
 }
 
-void ProcessFields::WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT*** array, unsigned int* numLines, unsigned int precision)
+void ProcessFields::WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT const* const* const* array, unsigned int const* numLines, unsigned int precision)
 {
 	file << "SCALARS " << name << " float " << 1 << endl;
 	file << "LOOKUP_TABLE default" << endl;
@@ -289,14 +288,14 @@ void ProcessFields::WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT*
 	}
 }
 
-bool ProcessFields::DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT*** array, double** discLines, unsigned int* numLines, unsigned int precision)
+bool ProcessFields::DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT const* const* const* array, double const* const* discLines, unsigned int const* numLines, unsigned int precision)
 {
 	WriteVTKHeader(file, discLines, numLines, precision);
 	WriteVTKScalarArray(file, name, array, numLines, precision);
 	return true;
 }
 
-bool ProcessFields::DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT*** array[], unsigned int numFields, double** discLines, unsigned int* numLines, unsigned int precision)
+bool ProcessFields::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)
 {
 	WriteVTKHeader(file, discLines, numLines);
 	for (unsigned int n=0;n<numFields;++n)
@@ -307,7 +306,7 @@ bool ProcessFields::DumpMultiScalarArray2VTK(ofstream &file, string names[], FDT
 	return true;
 }
 
-bool ProcessFields::DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT**** array, unsigned int* numLines)
+bool ProcessFields::DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines)
 {
 	const H5std_string FILE_NAME(filename);
 	const H5std_string DATASET_NAME( name );
diff --git a/FDTD/processfields.h b/FDTD/processfields.h
index 71611f4..a44f68f 100644
--- a/FDTD/processfields.h
+++ b/FDTD/processfields.h
@@ -52,22 +52,22 @@ public:
 	//! Set dump type: 0 for E-fields, 1 for H-fields, 2 for D-fields, 3 for B-fields, 4 for J-fields, etc...
 	void SetDumpType(DumpType type) {m_DumpType=type;}
 
-	static bool DumpVectorArray2VTK(ofstream &file, string name, FDTD_FLOAT**** array, double** discLines, unsigned int* numLines, unsigned int precision=12);
-	static bool DumpMultiVectorArray2VTK(ofstream &file, string names[], FDTD_FLOAT**** array[], unsigned int numFields, double** discLines, unsigned int* numLines, unsigned int precision=12);
-	static bool DumpScalarArray2VTK(ofstream &file, string name, FDTD_FLOAT*** array, double** discLines, unsigned int* numLines, unsigned int precision=12);
-	static bool DumpMultiScalarArray2VTK(ofstream &file, string names[], FDTD_FLOAT*** array[], unsigned int numFields, double** discLines, unsigned int* numLines, unsigned int precision=12);
+	static bool DumpVectorArray2VTK(ofstream &file, string name, FDTD_FLOAT const* const* const * const* array, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12);
+	static bool DumpMultiVectorArray2VTK(ofstream &file, string names[], FDTD_FLOAT const* const* const* const* const* array, unsigned int numFields, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12);
+	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);
+	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);
 
-	static bool DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT**** array, unsigned int* numLines);
+	static bool DumpVectorArray2HDF5(string filename, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines);
 
-	double CalcTotalEnergy();
+	double CalcTotalEnergy() const;
 
 	void SetFileType(FileType fileType) {m_fileType=fileType;}
 
 //	virtual void Process();
 protected:
-	static void WriteVTKHeader(ofstream &file, double** discLines, unsigned int* numLines, unsigned int precision=12);
-	static void WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT**** array, unsigned int* numLines, unsigned int precision=12);
-	static void WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT*** array, unsigned int* numLines, unsigned int precision=12);
+	static void WriteVTKHeader(ofstream &file, double const* const* discLines, unsigned int const* numLines, unsigned int precision=12);
+	static void WriteVTKVectorArray(ofstream &file, string name, FDTD_FLOAT const* const* const* const* array, unsigned int const* numLines, unsigned int precision=12);
+	static void WriteVTKScalarArray(ofstream &file, string name, FDTD_FLOAT const* const* const* array, unsigned int const* numLines, unsigned int precision=12);
 
 	static string GetFieldNameByType(DumpType type);
 	DumpMode m_DumpMode;
diff --git a/FDTD/processfields_td.cpp b/FDTD/processfields_td.cpp
index ffc6c00..f0d1604 100644
--- a/FDTD/processfields_td.cpp
+++ b/FDTD/processfields_td.cpp
@@ -31,9 +31,6 @@ ProcessFieldsTD::~ProcessFieldsTD()
 
 void ProcessFieldsTD::DumpCellInterpol(string filename)
 {
-	FDTD_FLOAT**** volt = Eng->GetVoltages();
-	FDTD_FLOAT**** curr = Eng->GetCurrents();
-
 	if (m_DumpType==E_FIELD_DUMP)
 	{
 		//create array
@@ -55,21 +52,21 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
 					delta  = Op->GetMeshDelta(0,OpPos); //Op->discLines[0][OpPos[0]+1] - Op->discLines[0][OpPos[0]];
 					if (delta)
 					{
-						E_T[0][pos[0]][pos[1]][pos[2]] = volt[0][OpPos[0]][OpPos[1]][OpPos[2]] + volt[0][OpPos[0]][OpPos[1]+1][OpPos[2]] + volt[0][OpPos[0]][OpPos[1]][OpPos[2]+1] + volt[0][OpPos[0]][OpPos[1]+1][OpPos[2]+1];
+						E_T[0][pos[0]][pos[1]][pos[2]] = Eng->GetVolt(0,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetVolt(0,OpPos[0],OpPos[1]+1,OpPos[2]) + Eng->GetVolt(0,OpPos[0],OpPos[1],OpPos[2]+1) + Eng->GetVolt(0,OpPos[0],OpPos[1]+1,OpPos[2]+1);
 						E_T[0][pos[0]][pos[1]][pos[2]] /= (4*delta);//*Op->gridDelta);
 					}
 					//in y
 					delta  = Op->GetMeshDelta(1,OpPos); //Op->discLines[1][OpPos[1]+1] - Op->discLines[1][OpPos[1]];
 					if (delta)
 					{
-						E_T[1][pos[0]][pos[1]][pos[2]] = volt[1][OpPos[0]][OpPos[1]][OpPos[2]] + volt[1][OpPos[0]+1][OpPos[1]][OpPos[2]] + volt[1][OpPos[0]][OpPos[1]][OpPos[2]+1] + volt[1][OpPos[0]+1][OpPos[1]][OpPos[2]+1];
+						E_T[1][pos[0]][pos[1]][pos[2]] = Eng->GetVolt(1,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetVolt(1,OpPos[0]+1,OpPos[1],OpPos[2]) + Eng->GetVolt(1,OpPos[0],OpPos[1],OpPos[2]+1) + Eng->GetVolt(1,OpPos[0]+1,OpPos[1],OpPos[2]+1);
 						E_T[1][pos[0]][pos[1]][pos[2]] /= (4*delta);//*Op->gridDelta);
 					}
 					//in z
 					delta  = Op->GetMeshDelta(2,OpPos); //Op->discLines[2][OpPos[2]+1] - Op->discLines[2][OpPos[2]];
 					if (delta)
 					{
-						E_T[2][pos[0]][pos[1]][pos[2]] = volt[2][OpPos[0]][OpPos[1]][OpPos[2]] + volt[2][OpPos[0]][OpPos[1]+1][OpPos[2]] + volt[2][OpPos[0]+1][OpPos[1]][OpPos[2]] + volt[2][OpPos[0]+1][OpPos[1]+1][OpPos[2]];
+						E_T[2][pos[0]][pos[1]][pos[2]] = Eng->GetVolt(2,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetVolt(2,OpPos[0],OpPos[1]+1,OpPos[2]) + Eng->GetVolt(2,OpPos[0]+1,OpPos[1],OpPos[2]) + Eng->GetVolt(2,OpPos[0]+1,OpPos[1]+1,OpPos[2]);
 						E_T[2][pos[0]][pos[1]][pos[2]] /= (4*delta);//*Op->gridDelta);
 					}
 				}
@@ -118,7 +115,7 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
 					delta = Op->GetDiscLine(0,OpPos[0],true);
 					if (delta)
 					{
-						H_T[0][pos[0]][pos[1]][pos[2]] = curr[0][OpPos[0]][OpPos[1]][OpPos[2]] + curr[0][OpPos[0]+1][OpPos[1]][OpPos[2]];
+						H_T[0][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(0,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(0,OpPos[0]+1,OpPos[1],OpPos[2]);
 						H_T[0][pos[0]][pos[1]][pos[2]] /= (2*delta);//*Op->gridDelta);
 					}
 					//in y
@@ -126,7 +123,7 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
 					delta = Op->GetDiscLine(1,OpPos[1],true);
 					if (delta)
 					{
-						H_T[1][pos[0]][pos[1]][pos[2]] = curr[1][OpPos[0]][OpPos[1]][OpPos[2]] + curr[1][OpPos[0]][OpPos[1]+1][OpPos[2]];
+						H_T[1][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(1,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(1,OpPos[0],OpPos[1]+1,OpPos[2]);
 						H_T[1][pos[0]][pos[1]][pos[2]] /= (2*delta);//*Op->gridDelta);
 					}
 					//in z
@@ -134,7 +131,7 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
 					delta = Op->GetDiscLine(2,OpPos[2],true);
 					if (delta)
 					{
-						H_T[2][pos[0]][pos[1]][pos[2]] = curr[2][OpPos[0]][OpPos[1]][OpPos[2]] + curr[2][OpPos[0]][OpPos[1]][OpPos[2]+1];
+						H_T[2][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(2,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(2,OpPos[0],OpPos[1],OpPos[2]+1);
 						H_T[2][pos[0]][pos[1]][pos[2]] /= (2*delta);//*Op->gridDelta);
 					}
 				}
@@ -162,9 +159,6 @@ void ProcessFieldsTD::DumpCellInterpol(string filename)
 
 void ProcessFieldsTD::DumpNoInterpol(string filename)
 {
-	FDTD_FLOAT**** volt = Eng->GetVoltages();
-	FDTD_FLOAT**** curr = Eng->GetCurrents();
-
 	unsigned int pos[3];
 	unsigned int OpPos[3];
 	double delta[3];
@@ -185,11 +179,11 @@ void ProcessFieldsTD::DumpNoInterpol(string filename)
 					OpPos[2]=start[2]+pos[2]*subSample[2];
 					delta[2]=Op->GetMeshDelta(2,OpPos);//fabs(Op->MainOp->GetIndexDelta(2,OpPos[2]));
 					if (delta[0])
-						E_T[0][pos[0]][pos[1]][pos[2]] = volt[0][OpPos[0]][OpPos[1]][OpPos[2]]/delta[0];// /Op->gridDelta;
+						E_T[0][pos[0]][pos[1]][pos[2]] = Eng->GetVolt(0,OpPos[0],OpPos[1],OpPos[2])/delta[0];// /Op->gridDelta;
 					if (delta[1])
-						E_T[1][pos[0]][pos[1]][pos[2]] = volt[1][OpPos[0]][OpPos[1]][OpPos[2]]/delta[1];// /Op->gridDelta;
+						E_T[1][pos[0]][pos[1]][pos[2]] = Eng->GetVolt(1,OpPos[0],OpPos[1],OpPos[2])/delta[1];// /Op->gridDelta;
 					if (delta[2])
-						E_T[2][pos[0]][pos[1]][pos[2]] = volt[2][OpPos[0]][OpPos[1]][OpPos[2]]/delta[2];// /Op->gridDelta;
+						E_T[2][pos[0]][pos[1]][pos[2]] = Eng->GetVolt(2,OpPos[0],OpPos[1],OpPos[2])/delta[2];// /Op->gridDelta;
 				}
 			}
 		}
@@ -231,11 +225,11 @@ void ProcessFieldsTD::DumpNoInterpol(string filename)
 					delta[2]=Op->GetMeshDelta(2,OpPos,true);//fabs(Op->MainOp->GetIndexWidth(2,OpPos[2]));
 					//in x
 					if (delta[0])
-						H_T[0][pos[0]][pos[1]][pos[2]] = curr[0][OpPos[0]][OpPos[1]][OpPos[2]]/delta[0];// /Op->gridDelta;
+						H_T[0][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(0,OpPos[0],OpPos[1],OpPos[2])/delta[0];// /Op->gridDelta;
 					if (delta[1])
-						H_T[1][pos[0]][pos[1]][pos[2]] = curr[1][OpPos[0]][OpPos[1]][OpPos[2]]/delta[1];// /Op->gridDelta;
+						H_T[1][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(1,OpPos[0],OpPos[1],OpPos[2])/delta[1];// /Op->gridDelta;
 					if (delta[2])
-						H_T[2][pos[0]][pos[1]][pos[2]] = curr[2][OpPos[0]][OpPos[1]][OpPos[2]]/delta[2];// /Op->gridDelta;
+						H_T[2][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(2,OpPos[0],OpPos[1],OpPos[2])/delta[2];// /Op->gridDelta;
 				}
 			}
 		}
diff --git a/FDTD/processing.cpp b/FDTD/processing.cpp
index 8467da9..d557fd6 100644
--- a/FDTD/processing.cpp
+++ b/FDTD/processing.cpp
@@ -91,33 +91,53 @@ void Processing::DefineStartStopCoord(double* dstart, double* dstop)
 	if (Op->SnapToMesh(dstop,stop)==false) cerr << "Processing::DefineStartStopCoord: Warning: Snapped line outside field domain!!" << endl;
 }
 
-double Processing::CalcLineIntegral(unsigned int* start, unsigned int* stop, int field)
+double Processing::CalcLineIntegral(unsigned int* start, unsigned int* stop, int field) const
+{
+	switch (field) {
+	case 0:
+		return CalcLineIntegral_V( start, stop );
+	case 1:
+		return CalcLineIntegral_I( start, stop );
+	}
+	return 0;
+}
+
+double Processing::CalcLineIntegral_I(unsigned int* start, unsigned int* stop) const
 {
 	double result=0;
-	FDTD_FLOAT**** array;
-	if (field==0)
-		array=Eng->GetVoltages();
-	else if (field==1)
-		array=Eng->GetCurrents();
-	else return 0.0;
-
 	for (int n=0;n<3;++n)
 	{
 		if (start[n]<stop[n])
 		{
 			unsigned int pos[3]={start[0],start[1],start[2]};
 			for (;pos[n]<stop[n];++pos[n])
-			{
-				result+=array[n][pos[0]][pos[1]][pos[2]];
-			}
+				result += Eng->GetCurr(n,pos[0],pos[1],pos[2]);
 		}
 		else
 		{
 			unsigned int pos[3]={stop[0],stop[1],stop[2]};
 			for (;pos[n]<start[n];++pos[n])
-			{
-				result-=array[n][pos[0]][pos[1]][pos[2]];
-			}
+				result -= Eng->GetCurr(n,pos[0],pos[1],pos[2]);
+		}
+	}
+	return result;
+}
+double Processing::CalcLineIntegral_V(unsigned int* start, unsigned int* stop) const
+{
+	double result=0;
+	for (int n=0;n<3;++n)
+	{
+		if (start[n]<stop[n])
+		{
+			unsigned int pos[3]={start[0],start[1],start[2]};
+			for (;pos[n]<stop[n];++pos[n])
+				result += Eng->GetVolt(n,pos[0],pos[1],pos[2]);
+		}
+		else
+		{
+			unsigned int pos[3]={stop[0],stop[1],stop[2]};
+			for (;pos[n]<start[n];++pos[n])
+				result -= Eng->GetVolt(n,pos[0],pos[1],pos[2]);
 		}
 	}
 	return result;
diff --git a/FDTD/processing.h b/FDTD/processing.h
index 5224b51..95069a2 100644
--- a/FDTD/processing.h
+++ b/FDTD/processing.h
@@ -76,7 +76,9 @@ protected:
 	ofstream file;
 	string m_filename;
 
-	double CalcLineIntegral(unsigned int* start, unsigned int* stop, int field);
+	double CalcLineIntegral(unsigned int* start, unsigned int* stop, int field) const;
+	double CalcLineIntegral_I(unsigned int* start, unsigned int* stop) const;
+	double CalcLineIntegral_V(unsigned int* start, unsigned int* stop) const;
 };
 
 class ProcessingArray