Commit ffda2630 authored by Ruben Undheim's avatar Ruben Undheim

Import GIT HEAD of openEMS sub-project (with Python interface)

parent 086965be
......@@ -122,15 +122,7 @@ find_package(Boost 1.46 COMPONENTS
)
# vtk
if (WIN32)
find_package(VTK 6.1 REQUIRED)
else()
# prefer >=6.1, fallback to >=5.4
find_package(VTK 6.1 COMPONENTS vtkIOXML vtkIOGeometry vtkIOLegacy vtkIOPLY NO_MODULE)
IF (NOT ${VTK_FOUND})
find_package(VTK REQUIRED)
endif()
endif()
find_package(VTK COMPONENTS vtkIOXML vtkIOGeometry vtkIOLegacy vtkIOPLY NO_MODULE REQUIRED)
message(STATUS "Found package VTK. Using version " ${VTK_VERSION})
if("${VTK_MAJOR_VERSION}" GREATER 5)
......
......@@ -56,6 +56,10 @@ public:
virtual double* GetJField(const unsigned int* pos, double* out) const =0;
//! Get the total current density field by rot(H) at \p pos. \sa SetInterpolationType
virtual double* GetRotHField(const unsigned int* pos, double* out) const =0;
//! Get the (interpolated) electric flux density field at \p pos. \sa SetInterpolationType
virtual double* GetDField(const unsigned int* pos, double* out) const =0;
//! Get the (interpolated) magnetic flux density field at \p pos. \sa SetInterpolationType
virtual double* GetBField(const unsigned int* pos, double* out) const =0;
//! Calculate the electric field integral along a given line
virtual double CalcVoltageIntegral(const unsigned int* start, const unsigned int* stop) const =0;
......
......@@ -71,6 +71,10 @@ string ProcessFields::GetFieldNameByType(DumpType type)
return "J-Field";
case ROTH_FIELD_DUMP:
return "RotH-Field";
case D_FIELD_DUMP:
return "D-Field";
case B_FIELD_DUMP:
return "B-Field";
case SAR_LOCAL_DUMP:
return "SAR-local";
case SAR_1G_DUMP:
......@@ -95,6 +99,30 @@ bool ProcessFields::NeedConductivity() const
return false;
}
bool ProcessFields::NeedPermittivity() const
{
switch (m_DumpType)
{
case D_FIELD_DUMP:
return true;
default:
return false;
}
return false;
}
bool ProcessFields::NeedPermeability() const
{
switch (m_DumpType)
{
case B_FIELD_DUMP:
return true;
default:
return false;
}
return false;
}
void ProcessFields::InitProcess()
{
if (Enabled==false) return;
......@@ -333,6 +361,44 @@ FDTD_FLOAT**** ProcessFields::CalcField()
}
}
return field;
case D_FIELD_DUMP:
for (unsigned int i=0; i<numLines[0]; ++i)
{
pos[0]=posLines[0][i];
for (unsigned int j=0; j<numLines[1]; ++j)
{
pos[1]=posLines[1][j];
for (unsigned int k=0; k<numLines[2]; ++k)
{
pos[2]=posLines[2][k];
m_Eng_Interface->GetDField(pos,out);
field[0][i][j][k] = out[0];
field[1][i][j][k] = out[1];
field[2][i][j][k] = out[2];
}
}
}
return field;
case B_FIELD_DUMP:
for (unsigned int i=0; i<numLines[0]; ++i)
{
pos[0]=posLines[0][i];
for (unsigned int j=0; j<numLines[1]; ++j)
{
pos[1]=posLines[1][j];
for (unsigned int k=0; k<numLines[2]; ++k)
{
pos[2]=posLines[2][k];
m_Eng_Interface->GetBField(pos,out);
field[0][i][j][k] = out[0];
field[1][i][j][k] = out[1];
field[2][i][j][k] = out[2];
}
}
}
return field;
default:
cerr << "ProcessFields::CalcField(): Error, unknown dump type..." << endl;
return field;
......
......@@ -40,7 +40,7 @@ public:
Current dump types are electric field (E_FIELD_DUMP), magnetic field (H_FIELD_DUMP),
(conduction) electric current density (kappa*E) (J_FIELD_DUMP) and total current density (rotH)
*/
enum DumpType { E_FIELD_DUMP=0, H_FIELD_DUMP=1, J_FIELD_DUMP=2, ROTH_FIELD_DUMP=3, SAR_LOCAL_DUMP=20, SAR_1G_DUMP=21, SAR_10G_DUMP=22, SAR_RAW_DATA=29};
enum DumpType { E_FIELD_DUMP=0, H_FIELD_DUMP=1, J_FIELD_DUMP=2, ROTH_FIELD_DUMP=3, D_FIELD_DUMP=4, B_FIELD_DUMP=5, SAR_LOCAL_DUMP=20, SAR_1G_DUMP=21, SAR_10G_DUMP=22, SAR_RAW_DATA=29};
virtual std::string GetProcessingName() const {return "common field processing";}
......@@ -75,6 +75,8 @@ public:
static std::string GetFieldNameByType(DumpType type);
virtual bool NeedConductivity() const;
virtual bool NeedPermittivity() const;
virtual bool NeedPermeability() const;
protected:
DumpType m_DumpType;
......
......@@ -174,7 +174,7 @@ void ProcessModeMatch::InitProcess()
for (int n=0; n<2; ++n)
{
m_ModeDist[n][posP][posPP] = m_ModeParser[n]->Eval(var); //calc mode template
if ((isnan(m_ModeDist[n][posP][posPP])) || (isinf(m_ModeDist[n][posP][posPP])))
if ((std::isnan(m_ModeDist[n][posP][posPP])) || (std::isinf(m_ModeDist[n][posP][posPP])))
m_ModeDist[n][posP][posPP] = 0.0;
norm += pow(m_ModeDist[n][posP][posPP],2) * area;
}
......
......@@ -47,6 +47,11 @@ double* Engine_Interface_FDTD::GetJField(const unsigned int* pos, double* out) c
return GetRawInterpolatedField(pos, out, 1);
}
double* Engine_Interface_FDTD::GetDField(const unsigned int* pos, double* out) const
{
return GetRawInterpolatedField(pos, out, 3);
}
double* Engine_Interface_FDTD::GetRotHField(const unsigned int* pos, double* out) const
{
return GetRawInterpolatedField(pos, out, 2);
......@@ -116,6 +121,27 @@ double* Engine_Interface_FDTD::GetRawInterpolatedField(const unsigned int* pos,
}
double* Engine_Interface_FDTD::GetHField(const unsigned int* pos, double* out) const
{
return GetRawInterpolatedDualField(pos, out, 0);
}
double* Engine_Interface_FDTD::GetBField(const unsigned int* pos, double* out) const
{
return GetRawInterpolatedDualField(pos, out, 1);
}
double Engine_Interface_FDTD::GetRawDualField(unsigned int n, const unsigned int* pos, int type) const
{
double value = m_Eng->GetCurr(n,pos[0],pos[1],pos[2]);
double delta = m_Op->GetEdgeLength(n,pos,true);
if ((type==0) && (delta))
return value/delta;
if ((type==1) && (m_Op->m_mueR) && (delta))
return value*m_Op->m_mueR[n][pos[0]][pos[1]][pos[2]]/delta;
return 0.0;
}
double* Engine_Interface_FDTD::GetRawInterpolatedDualField(const unsigned int* pos, double* out, int type) const
{
unsigned int iPos[] = {pos[0],pos[1],pos[2]};
int nP,nPP;
......@@ -124,9 +150,9 @@ double* Engine_Interface_FDTD::GetHField(const unsigned int* pos, double* out) c
{
default:
case NO_INTERPOLATION:
out[0] = m_Eng->GetCurr(0,pos) / m_Op->GetEdgeLength(0,pos,true);
out[1] = m_Eng->GetCurr(1,pos) / m_Op->GetEdgeLength(1,pos,true);
out[2] = m_Eng->GetCurr(2,pos) / m_Op->GetEdgeLength(2,pos,true);
out[0] = GetRawDualField(0, pos, type);
out[1] = GetRawDualField(1, pos, type);
out[2] = GetRawDualField(2, pos, type);
break;
case NODE_INTERPOLATE:
for (int n=0; n<3; ++n)
......@@ -138,13 +164,13 @@ double* Engine_Interface_FDTD::GetHField(const unsigned int* pos, double* out) c
out[n] = 0;
continue;
}
out[n]=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
out[n] = GetRawDualField(n, iPos, type);
--iPos[nP];
out[n]+=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
out[n]+= GetRawDualField(n, iPos, type);
--iPos[nPP];
out[n]+=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
out[n]+= GetRawDualField(n, iPos, type);
++iPos[nP];
out[n]+=m_Eng->GetCurr(n,iPos)/m_Op->GetEdgeLength(n,iPos,true);
out[n]+= GetRawDualField(n, iPos, type);
++iPos[nPP];
out[n]/=4;
}
......@@ -153,7 +179,7 @@ double* Engine_Interface_FDTD::GetHField(const unsigned int* pos, double* out) c
for (int n=0; n<3; ++n)
{
delta = m_Op->GetEdgeLength(n,iPos,true);
out[n] = m_Eng->GetCurr(n,iPos);
out[n] = GetRawDualField(n, iPos, type);
if ((pos[n]>=m_Op->GetNumberOfLines(n,true)-1))
{
out[n] = 0; //magnetic field on the outer boundaries is always zero
......@@ -162,7 +188,7 @@ double* Engine_Interface_FDTD::GetHField(const unsigned int* pos, double* out) c
++iPos[n];
double deltaUp = m_Op->GetEdgeLength(n,iPos,true);
double deltaRel = delta / (delta+deltaUp);
out[n] = out[n]*(1.0-deltaRel)/delta + (double)m_Eng->GetCurr(n,iPos)/deltaUp*deltaRel;
out[n] = out[n]*(1.0-deltaRel) + (double)GetRawDualField(n, iPos, type)*deltaRel;
--iPos[n];
}
break;
......@@ -201,6 +227,8 @@ double Engine_Interface_FDTD::GetRawField(unsigned int n, const unsigned int* po
return value/delta;
if ((type==1) && (m_Op->m_kappa) && (delta))
return value*m_Op->m_kappa[n][pos[0]][pos[1]][pos[2]]/delta;
if ((type==3) && (m_Op->m_epsR) && (delta))
return value*m_Op->m_epsR[n][pos[0]][pos[1]][pos[2]]/delta;
if (type==2) //calc rot(H)
{
int nP = (n+1)%3;
......
......@@ -44,6 +44,8 @@ public:
virtual double* GetHField(const unsigned int* pos, double* out) const;
virtual double* GetJField(const unsigned int* pos, double* out) const;
virtual double* GetRotHField(const unsigned int* pos, double* out) const;
virtual double* GetDField(const unsigned int* pos, double* out) const;
virtual double* GetBField(const unsigned int* pos, double* out) const;
virtual double CalcVoltageIntegral(const unsigned int* start, const unsigned int* stop) const;
......@@ -56,10 +58,15 @@ protected:
Operator* m_Op;
Engine* m_Eng;
//! Internal method to get an interpolated field of a given type. (0: E, 1: J, 2: rotH)
//! Internal method to get an interpolated field of a given type. (0: E, 1: J, 2: rotH, 3: D)
virtual double* GetRawInterpolatedField(const unsigned int* pos, double* out, int type) const;
//! Internal method to get a raw field of a given type. (0: E, 1: J, 2: rotH)
//! Internal method to get a raw field of a given type. (0: E, 1: J, 2: rotH, 3: D)
virtual double GetRawField(unsigned int n, const unsigned int* pos, int type) const;
//! Internal method to get an interpolated dual field of a given type. (0: H, 1: B)
virtual double* GetRawInterpolatedDualField(const unsigned int* pos, double* out, int type) const;
//! Internal method to get a raw dual field of a given type. (0: H, 1: B)
virtual double GetRawDualField(unsigned int n, const unsigned int* pos, int type) const;
};
#endif // ENGINE_INTERFACE_FDTD_H
......@@ -32,7 +32,10 @@
#include "boost/date_time/posix_time/posix_time.hpp"
#include "boost/date_time/gregorian/gregorian.hpp"
#include <iomanip>
#ifndef SSE_CORRECT_DENORMALS
#include <xmmintrin.h>
#endif
//! \brief construct an Engine_Multithread instance
//! it's the responsibility of the caller to free the returned pointer
......
......@@ -15,7 +15,10 @@
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef SSE_CORRECT_DENORMALS
#include <xmmintrin.h>
#endif
#include "engine_sse.h"
//! \brief construct an Engine_sse instance
......
......@@ -53,7 +53,7 @@ void Excitation::Reset( double timestep )
m_foi = 0;
}
bool Excitation::SetupGaussianPulse(double f0, double fc)
void Excitation::SetupGaussianPulse(double f0, double fc)
{
m_Excit_Type = Excitation::GaissianPulse;
m_f0 = f0;
......@@ -62,7 +62,7 @@ bool Excitation::SetupGaussianPulse(double f0, double fc)
m_SignalPeriod = 0;
}
bool Excitation::SetupSinusoidal(double f0)
void Excitation::SetupSinusoidal(double f0)
{
m_Excit_Type = Excitation::Sinusoidal;
m_f0 = f0;
......@@ -70,21 +70,21 @@ bool Excitation::SetupSinusoidal(double f0)
m_SignalPeriod = 1/f0;
}
bool Excitation::SetupDiracPulse(double fmax)
void Excitation::SetupDiracPulse(double fmax)
{
m_Excit_Type = Excitation::DiracPulse;
m_SignalPeriod = 0;
m_f_max = fmax;
}
bool Excitation::SetupStepExcite(double fmax)
void Excitation::SetupStepExcite(double fmax)
{
m_Excit_Type = Excitation::Step;
m_SignalPeriod = 0;
m_f_max = fmax;
}
bool Excitation::SetupCustomExcite(string str, double f0, double fmax)
void Excitation::SetupCustomExcite(string str, double f0, double fmax)
{
m_Excit_Type = Excitation::CustomExcite;
m_CustomExc_Str = str;
......
......@@ -31,11 +31,11 @@ public:
virtual void Reset( double timestep );
bool SetupGaussianPulse(double f0, double fc);
bool SetupSinusoidal(double f0);
bool SetupDiracPulse(double fmax);
bool SetupStepExcite(double fmax);
bool SetupCustomExcite(std::string str, double f0, double fmax);
void SetupGaussianPulse(double f0, double fc);
void SetupSinusoidal(double f0);
void SetupDiracPulse(double fmax);
void SetupStepExcite(double fmax);
void SetupCustomExcite(std::string str, double f0, double fmax);
double GetCenterFreq() {return m_f0;}
double GetCutOffFreq() {return m_fc;}
......
......@@ -225,6 +225,9 @@ bool Operator_Ext_LorentzMaterial::BuildExtension()
// CSProperties* prop = m_Op->GetGeometryCSX()->GetPropertyByCoordPriority(coord,(CSProperties::PropertyType)(CSProperties::METAL | CSProperties::MATERIAL), true);
CSProperties* prop = m_Op->GetGeometryCSX()->GetPropertyByCoordPriority(coord, vPrims, true);
if (prop==NULL) continue;
if ((mat = prop->ToLorentzMaterial()))
{
w_plasma = mat->GetEpsPlasmaFreqWeighted(order,n,coord) * 2 * PI;
......@@ -277,6 +280,9 @@ bool Operator_Ext_LorentzMaterial::BuildExtension()
// CSProperties* prop = m_Op->GetGeometryCSX()->GetPropertyByCoordPriority(coord,(CSProperties::PropertyType)(CSProperties::METAL | CSProperties::MATERIAL), true);
CSProperties* prop = m_Op->GetGeometryCSX()->GetPropertyByCoordPriority(coord, vPrims, true);
if (prop==NULL) continue;
if ((mat = prop->ToLorentzMaterial()))
{
w_plasma = mat->GetMuePlasmaFreqWeighted(order,n,coord) * 2 * PI;
......
......@@ -159,7 +159,7 @@ bool Operator_Ext_TFSF::BuildExtension()
else
m_PhVel=m_Op->CalcNumericPhaseVelocity(m_Start,m_Stop,m_PropDir,m_Frequency);
if ((m_PhVel<0) || (m_PhVel>__C0__/ref_index) || isnan(m_PhVel))
if ((m_PhVel<0) || (m_PhVel>__C0__/ref_index) || std::isnan(m_PhVel))
{
cerr << "Operator_Ext_TFSF::BuildExtension: Warning, invalid phase velocity found, resetting to c0! " << endl;
m_PhVel = __C0__/ref_index;
......
......@@ -1395,7 +1395,7 @@ bool Operator::AverageMatCellCenter(int ny, const unsigned int* pos, double* Eff
if (EffMat[3]) EffMat[3]=length / EffMat[3];
for (int n=0; n<4; ++n)
if (isnan(EffMat[n]) || isinf(EffMat[n]))
if (std::isnan(EffMat[n]) || std::isinf(EffMat[n]))
{
cerr << "Operator::" << __func__ << ": Error, an effective material parameter is not a valid result, this should NOT have happend... exit..." << endl;
cerr << ny << "@" << n << " : " << pos[0] << "," << pos[1] << "," << pos[2] << endl;
......@@ -1508,7 +1508,7 @@ bool Operator::AverageMatQuarterCell(int ny, const unsigned int* pos, double* Ef
if (EffMat[3]) EffMat[3]=length / EffMat[3];
for (int n=0; n<4; ++n)
if (isnan(EffMat[n]) || isinf(EffMat[n]))
if (std::isnan(EffMat[n]) || std::isinf(EffMat[n]))
{
cerr << "Operator::" << __func__ << ": Error, An effective material parameter is not a valid result, this should NOT have happend... exit..." << endl;
cerr << ny << "@" << n << " : " << pos[0] << "," << pos[1] << "," << pos[2] << endl;
......@@ -1555,7 +1555,7 @@ bool Operator::Calc_LumpedElements()
if (R<0)
R = NAN;
if ((isnan(R)) && (isnan(C)))
if ((std::isnan(R)) && (std::isnan(C)))
{
cerr << "Operator::Calc_LumpedElements(): Warning: Lumped Element R or C not specified! skipping. "
<< " ID: " << prims.at(bn)->GetID() << " @ Property: " << PLE->GetName() << endl;
......
......@@ -41,11 +41,11 @@ Operator_sse::~Operator_sse()
Delete();
}
Engine* Operator_sse::CreateEngine() const
Engine* Operator_sse::CreateEngine()
{
//! create a special sse-engine
Engine_sse* eng = Engine_sse::New(this);
return eng;
m_Engine = Engine_sse::New(this);
return m_Engine;
}
void Operator_sse::Init()
......
......@@ -29,7 +29,7 @@ public:
static Operator_sse* New();
virtual ~Operator_sse();
virtual Engine* CreateEngine() const;
virtual Engine* CreateEngine();
inline virtual FDTD_FLOAT GetVV( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return f4_vv[n][x][y][z%numVectors].f[z/numVectors]; }
inline virtual FDTD_FLOAT GetVI( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { return f4_vi[n][x][y][z%numVectors].f[z/numVectors]; }
......
......@@ -42,16 +42,14 @@ Operator_SSE_Compressed::~Operator_SSE_Compressed()
Delete();
}
Engine* Operator_SSE_Compressed::CreateEngine() const
Engine* Operator_SSE_Compressed::CreateEngine()
{
if (!m_Use_Compression)
{
//! create a default sse-engine
Engine_sse* eng = Engine_sse::New(this);
return eng;
}
Engine_SSE_Compressed* eng = Engine_SSE_Compressed::New(this);
return eng;
m_Engine = Engine_sse::New(this);
else
m_Engine = Engine_SSE_Compressed::New(this);
return m_Engine;
}
int Operator_SSE_Compressed::CalcECOperator( DebugFlags debugFlags )
......
......@@ -43,7 +43,7 @@ public:
static Operator_SSE_Compressed* New();
virtual ~Operator_SSE_Compressed();
virtual Engine* CreateEngine() const;
virtual Engine* CreateEngine();
inline virtual FDTD_FLOAT GetVV( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { if (m_Use_Compression) return f4_vv_Compressed[n][m_Op_index[x][y][z%numVectors]].f[z/numVectors]; else return Operator_sse::GetVV(n,x,y,z);}
inline virtual FDTD_FLOAT GetVI( unsigned int n, unsigned int x, unsigned int y, unsigned int z ) const { if (m_Use_Compression) return f4_vi_Compressed[n][m_Op_index[x][y][z%numVectors]].f[z/numVectors]; else return Operator_sse::GetVI(n,x,y,z);}
......
function [delay, fidelity, nf2ff_out] = DelayFidelity(nf2ff, port, path, weight_theta, weight_phi, theta, phi, f_0, f_c, varargin)
% [delay, fidelity] = DelayFidelity(nf2ff, port, path, theta, phi, f_lo, f_hi, varargin)
% [delay, fidelity] = DelayFidelity(nf2ff, port, path, weight_theta, weight_phi, theta, phi, f_lo, f_hi, varargin)
%
%
% This function calculates the time delay from the source port to the phase center of the antenna and the fidelity.
......@@ -11,8 +11,8 @@ function [delay, fidelity, nf2ff_out] = DelayFidelity(nf2ff, port, path, weight_
% nf2ff: return value of CreateNF2FFBox.
% port: return value of AddLumpedPort
% path: path of the simulation results.
% weight_theta: weight if the E_theta component
% weight_phi: eight of the E_phi component
% weight_theta: weight of the E_theta component
% weight_phi: weight of the E_phi component
% -> with both (possibly complex) parameters any polarization can be examined
% theta: theta values to be simulated
% phi: phi values to be simulated
......@@ -27,7 +27,8 @@ function [delay, fidelity, nf2ff_out] = DelayFidelity(nf2ff, port, path, weight_
% example:
% theta = [-180:10:180] * pi / 180;
% phi = [0, 90] * pi / 180;
% [delay, fidelity] = DelayFidelity2(nf2ff, port, Sim_Path, sin(tilt), cos(tilt), theta, phi, f_0, f_c, 'Mode', 1);
% % use circular right handed polarization
% [delay, fidelity] = DelayFidelity2(nf2ff, port, Sim_Path, -1i, 1, theta, phi, f_0, f_c, 'Mode', 1);
% figure
% polar(theta.', delay(:,1) * 3e11); % delay in mm
% figure
......
......@@ -170,11 +170,11 @@ if (enable_ssh)
warning('openEMS:RunOpenEMS','remote cleanup failed!');
end
else
args = [Sim_File ' ' opts];
args = ['"' Sim_File '" ' opts]
if isfield(Settings,'LogFile') && isfield(Settings,'Silent')
invoke_openEMS(args,Settings.LogFile,Settings.Silent);
invoke_openEMS(args,['"' Settings.LogFile '"'],Settings.Silent);
elseif isfield(Settings,'LogFile')
invoke_openEMS(args,Settings.LogFile);
invoke_openEMS(args,['"' Settings.LogFile '"']);
elseif isfield(Settings,'Silent')
invoke_openEMS(args,[],Settings.Silent);
else
......
......@@ -18,9 +18,11 @@ function h = plotRefl(port, varargin)
% for the frequency (always in MHz), default is 0
% 'threshold': - Threshold value (in dB) for the upper and lower
% cutoff frequency, default is -3
% 'fmarkers': - set lower and upper frequency marker in Hz manually,
% like so: [4e9, 6.5e9]
% example:
% myport = calcPort(myport, Sim_Path, linspace(f_0-f_c, f_0+f_c, 200));
% plotRefl(myport);
% plotRefl(myport, 'fmarkers', [4e9, 6.5e9]);
%
% See also calcPort
%
......@@ -31,14 +33,16 @@ function h = plotRefl(port, varargin)
%defaults
precision = 0;
threshold = -3;
fmarkers = [];
for n=1:2:numel(varargin)
if (strcmp(varargin{n},'precision')==1);
precision = varargin{n+1};
elseif (strcmp(varargin{n},'threshold')==1);
threshold = varargin{n+1};
else
elseif (strcmp(varargin{n},'fmarkers')==1);
fmarkers = varargin{n+1};
else
warning('openEMS:polarFF',['unknown argument key: ''' varargin{n} '''']);
end
end
......@@ -111,8 +115,14 @@ text (real (rho), imag (rho), '0');
s11dB = 20*log10(abs(s11));
if(isempty(fmarkers))
upperind = s11dB(1:end-1) < threshold & s11dB(2:end) > threshold;
lowerind = s11dB(1:end-1) > threshold & s11dB(2:end) < threshold;
else
upperind = [nthargout(2, @min, abs(fmarkers(2)-port.f))];
lowerind = [nthargout(2, @min, abs(fmarkers(1)-port.f))];
end
minind = nthargout(2, @min, s11dB);
handle1 = plot(s11(lowerind),['<','b']);
handle2 = plot(s11(upperind),['>','b']);
......@@ -120,6 +130,7 @@ handle3 = plot(s11(minind),['*', 'b']);
llegend = num2str(port.f(lowerind)(1)/1e6, ffmt);
ulegend = num2str(port.f(upperind)(1)/1e6, ffmt);
if nnz(lowerind) > 1
for i= 2:nnz(lowerind)
llegend = strjoin({llegend, num2str(port.f(lowerind)(i)/1e6, ffmt)}, ', ');
......@@ -142,4 +153,4 @@ if (nargout == 0)
clear h;
end
end
\ No newline at end of file
end
......@@ -19,14 +19,14 @@ if isOctave()
fflush(stdout);
if isunix
[res, fn_so] = unix('find /usr/lib -name libhdf5.so');
[res, fn_h] = unix('find /usr/include -name hdf5.h');
[res, fn_h] = unix('find /usr/include -name hdf5.h | head -1');
if length(fn_so)>0 && length(fn_h)>0
[hdf5lib_dir, hdf5lib_fn] = fileparts(fn_so);
[hdf5lib_dir, hdf5lib_fn, ext] = fileparts(fn_so);
disp(["HDF5 library path found at: " hdf5lib_dir])
[hdf5inc_dir, hdf5inc_fn] = fileparts(fn_h);
[hdf5inc_dir, hdf5inc_fn, ext] = fileparts(fn_h);
disp(["HDF5 include path found at: " hdf5inc_dir])
mkoctfile(["-L" hdf5lib_dir " -I" hdf5inc_dir],"-lhdf5", "h5readatt_octave.cc")
mkoctfile("h5readatt_octave.cc", ["-L" hdf5lib_dir], ["-I" hdf5inc_dir], "-L hdf5")
else
mkoctfile -lhdf5 h5readatt_octave.cc
end
......
......@@ -313,6 +313,10 @@ bool nf2ff_calc::AddPlane(float **lines, unsigned int* numLines, complex<float>*
if ((m_MirrorType[n]!=MIRROR_OFF) && (m_MirrorType[nP]==MIRROR_OFF) && (m_MirrorType[nPP]==MIRROR_OFF))
{
this->AddMirrorPlane(n, lines, numLines, E_field, H_field, MeshType);
for (unsigned int i=0;i<numLines[n];++i)
lines[n][i] = 2.0*m_MirrorPos[n] - lines[n][i];
break;
}
//check if two planes are on
......@@ -321,6 +325,10 @@ bool nf2ff_calc::AddPlane(float **lines, unsigned int* numLines, complex<float>*
this->AddMirrorPlane(nP, lines, numLines, E_field, H_field, MeshType);
this->AddMirrorPlane(nPP, lines, numLines, E_field, H_field, MeshType);
this->AddMirrorPlane(nP, lines, numLines, E_field, H_field, MeshType);
for (unsigned int i=0;i<numLines[nPP];++i)
lines[nPP][i] = 2.0*m_MirrorPos[nPP] - lines[nPP][i];
break;
}
}
......@@ -334,6 +342,10 @@ bool nf2ff_calc::AddPlane(float **lines, unsigned int* numLines, complex<float>*
this->AddMirrorPlane(0, lines, numLines, E_field, H_field, MeshType);
this->AddMirrorPlane(1, lines, numLines, E_field, H_field, MeshType);
this->AddMirrorPlane(0, lines, numLines, E_field, H_field, MeshType);
for (unsigned int i=0;i<numLines[2];++i)
lines[2][i] = 2.0*m_MirrorPos[2] - lines[2][i];
}
//cleanup E- & H-Fields
......@@ -496,7 +508,7 @@ bool nf2ff_calc::AddSinglePlane(float **lines, unsigned int* numLines, complex<f
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*fZ0) * abs((m_E_theta[tn][pn]*conj(m_E_theta[tn][pn])+m_E_phi[tn][pn]*conj(m_E_phi[tn][pn])));
m_P_rad[tn][pn] = abs((m_E_theta[tn][pn]*conj(m_E_theta[tn][pn])+m_E_phi[tn][pn]*conj(m_E_phi[tn][pn])))/(2*fZ0);
if (m_P_rad[tn][pn]>P_max)
P_max = m_P_rad[tn][pn];
}
......
......@@ -182,13 +182,13 @@ bool openEMS::parseCommandLineArgument( const char *argv )
else if (strcmp(argv,"--debug-PEC")==0)
{
cout << "openEMS - dumping PEC info to 'PEC_dump.vtk'" << endl;
m_debugPEC = true;
DebugPEC();
return true;
}
else if (strcmp(argv,"--debug-CSX")==0)
{
cout << "openEMS - dumping CSX geometry to 'debugCSX.xml'" << endl;
m_debugCSX = true;
DebugCSX();
return true;
}
else if (strcmp(argv,"--engine=basic")==0)
......@@ -217,7 +217,7 @@ bool openEMS::parseCommandLineArgument( const char *argv )
}
else if (strncmp(argv,"--numThreads=",13)==0)
{
m_engine_numThreads = atoi(argv+13);
this->SetNumberOfThreads(atoi(argv+13));
cout << "openEMS - fixed number of threads: " << m_engine_numThreads << endl;
return true;
}
......@@ -284,7 +284,7 @@ void openEMS::WelcomeScreen()
cout << " ---------------------------------------------------------------------- " << endl;
cout << " | openEMS " << bits << " -- version " GIT_VERSION << endl;
cout << " | (C) 2010-2016 Thorsten Liebig <thorsten.liebig@gmx.de> GPL license" << endl;
cout << " | (C) 2010-2018 Thorsten Liebig <thorsten.liebig@gmx.de> GPL license" << endl;
cout << " ---------------------------------------------------------------------- " << endl;
cout << openEMS::GetExtLibsInfo("\t") << endl;
}
......@@ -472,9 +472,9 @@ bool openEMS::SetupProcessing()