initial copy from offline repository

Author: mgotz

CMakeLists.txt

@@ -0,0 +1,14 @@
+cmake_minimum_required(VERSION 3.0.2)
+
+project(IK_Simulation)
+
+set(Boost_USE_STATIC_LIBS ON)
+find_package(Boost 1.55.0 REQUIRED program_options)
+
+include_directories(include ${Boost_INCLUDE_DIRS})
+
+file(GLOB HEADERS "include/*.hpp")
+file(GLOB SOURCES "src/*.cpp")
+
+add_executable(IK_Simulation ${SOURCES} ${HEADERS})
+target_link_libraries(IK_Simulation ${Boost_LIBRARIES})

README.md

@@ -1,2 +1,26 @@
-# numerical-ks-calculator
+# IC_Simulation
 numerical solution of charge transport, creatiion and recombination inside a 1D ionization chamber model
+
+Building:
+run cmake from a directory you want to build in and give a path to this directory as argument
+cmake will create a makefile, which can be used to build by simply runn make
+
+Requirements: standard library and program_options from the boost library
+Tested on gcc 4.9.2
+
+
+Alternatively compile manually. Source files are in src, headers are in include directories.
+Link all files and don't forget -lboost_program_options
+
+Usage:
+The program's parameters can be given either on the command line or in a configuration file (or any combination of the two)
+Run IK_Simluation -h to see a list of the parameters.
+
+Beyond the straightfoward settings, the program needs a chamber file and an active medium file.
+The chamber file gives the dimensions and properties (mainly calibration factor to calculate liberated charge from a given dose) 
+of the chamber.
+See AdvMarkus_kammer.txt for an example.
+
+The active_medium file specifies the properties of the filling gas.
+Those can be simple constants, function definitions or point to additional files, if they are given in a spline format
+See active_medium_air_with_direct.txt for an example.

include/active_medium.hpp

@@ -0,0 +1,58 @@
+#ifndef ACTIVE_MEDIUM
+#define ACTIVE_MEDIUM
+/*
+ this class defines an active medium, with several methods to return its properties, such as ion mobility or electron attachment depending on the electri field strength
+ it is constructed by reading an ini-style file with boost::property_tree
+*/
+
+
+#include "active_medium_functiods.hpp"
+
+#include <fstream>
+#include <limits>
+#include <vector>
+#include <string>
+#include <stdexcept>
+
+
+using namespace std;
+
+
+class active_medium
+{
+	public:
+        active_medium(std::string); //constructor takes file path to the active medium ini file
+
+        //methods to get the different properties:
+        double mobility_e (double);
+        double mobility_pos (double);
+        double mobility_neg (double);
+        double attachment (double);
+        double volume_recombination (double);
+        double direct_recombination (double);
+
+        double relative_permittivity(double);
+
+        void write_log(ofstream&); //write table of function values in stream
+
+        double mMaxE;
+	private:
+        //functiods for the different properties
+        vector<evaluation_functiod*> mFunctiods;
+};
+
+class unknown_arg_error : public runtime_error
+{
+	public:
+		unknown_arg_error(const string &message) : runtime_error(message) {}
+};
+class file_open_error : public runtime_error
+{
+	public:
+		file_open_error(const string &message) : runtime_error(message) {}
+};
+
+
+
+
+#endif

include/active_medium_functiods.hpp

@@ -0,0 +1,137 @@
+#ifndef ACTIVE_MEDIUM_FUNCTIODS
+#define ACTIVE_MEDIUM_FUNCTIODS
+/*
+ There is one abstract base class the evaluation_functiod and several implementation
+ Each instance of an implementation is a function with some parameters that are defined at construction, such that the function can then be evaluated
+ at different value with the same parameters
+*/
+
+
+#include <fstream>
+#include <cmath>
+#include <limits>
+
+class evaluation_functiod{
+	public:
+        virtual double evaluate(double) = 0; //=0 makes this an abstract class, i.e. method is not defined and needs to be defined in the children
+		virtual ~evaluation_functiod() = 0;
+        double mMaxE;
+};
+
+inline evaluation_functiod::~evaluation_functiod() {}
+
+
+
+class constant : public evaluation_functiod{
+	public:
+        virtual double evaluate(double ){ return(mIntercept);}
+        constant(double pIntercept)
+        {
+            mIntercept=pIntercept;
+            mMaxE = std::numeric_limits<double>::max();
+        }
+	private:
+        double mIntercept;
+
+	};
+
+class linear : public evaluation_functiod{
+	public:
+        virtual double evaluate(double pEField){ return(mIntercept+mSlope*pEField);}
+        linear(double pIntercept, double pSlope)
+        {
+            mIntercept=pIntercept;
+            mSlope = pSlope;
+            mMaxE = std::numeric_limits<double>::max();
+        }
+	private:
+        double mIntercept;
+        double mSlope;
+
+	};
+
+class exponential : public evaluation_functiod{
+	public:
+        virtual double evaluate(double pEField){ return(m0+m1*exp(m2*pEField));}
+        exponential (double p0, double p1, double p2)
+        {
+            m0 = p0;
+            m1 = p1;
+            m2 = p2;
+            mMaxE = std::numeric_limits<double>::max();
+        }
+	private:
+        double m0;
+        double m1;
+        double m2;
+
+	};
+
+class pol4 : public evaluation_functiod{
+	public:
+        virtual double evaluate(double pEField){
+        return(m0+m1*pEField+m2*pEField*pEField+m3*pEField*pEField*pEField+m4*pEField*pEField*pEField*pEField);
+		}
+        pol4 (double p0, double p1, double p2, double p3, double p4)
+        {
+            m0 = p0;
+            m1 = p1;
+            m2 = p2;
+            m3 = p3;
+            m4 = p4;
+            mMaxE = std::numeric_limits<double>::max();
+        }
+	private:
+        double m0;
+        double m1;
+        double m2;
+        double m3;
+        double m4;
+};
+
+class pol5 : public evaluation_functiod{
+	public:
+        virtual double evaluate(double pEField){
+        return(m0+m1*pEField+m2*pEField*pEField+m3*pEField*pEField*pEField+m4*pEField*pEField*pEField*pEField+m5*pEField*pEField*pEField*pEField*pEField);
+		}
+        pol5 (double p0, double p1, double p2, double p3, double p4, double p5)
+        {
+            m0 = p0;
+            m1 = p1;
+            m2 = p2;
+            m3 = p3;
+            m4 = p4;
+            m5 = p5;
+            mMaxE = std::numeric_limits<double>::max();
+        }
+	private:
+        double m0;
+        double m1;
+        double m2;
+        double m3;
+        double m4;
+        double m5;
+
+};
+
+class spline_eval_only : public evaluation_functiod{
+	public:
+        virtual double evaluate(double pEField){
+            i = (size_t)((pEField-mX0)/mGridspacing);
+            return (mA[i]*pEField*pEField*pEField+mB[i]*pEField*pEField+mC[i]*pEField+mD[i]);
+			}
+		spline_eval_only(std::ifstream&);
+		~spline_eval_only();
+
+	private:
+        double* mA;
+        double* mB;
+        double* mC;
+        double* mD;
+        double mGridspacing;
+        double mX0;
+        size_t i;
+};
+
+
+#endif

include/beam_model.hpp

@@ -0,0 +1,42 @@
+#ifndef BEAM_MODEL
+#define BEAM_MODEL
+/*
+describe the beam, for now either a conitnuous rate over the entire pulse duration or a micropulse structure like at ELBE
+*/
+#include "settings.hpp"
+
+#include <string>
+
+class beam_model{
+public:
+    virtual double ionization(double, double) = 0; //=0 makes this an abstract class, i.e. method is not defined and needs to be defined in the children
+    virtual void print_characeristics() = 0;
+    virtual ~beam_model() = 0;
+protected:
+    double mPulseDuration;
+};
+
+inline beam_model::~beam_model() {}
+
+class continuous_pulse : public beam_model{
+public:
+    virtual double ionization(double, double);
+    virtual void print_characeristics();
+    continuous_pulse (double, double); //give the charge created per ns and bin and the pulse duration
+private:
+    double mIonizationRate;
+};
+
+
+class elbe_micropulses : public beam_model{
+public:
+    virtual double ionization(double , double );
+    virtual void print_characeristics();
+    elbe_micropulses (double, double); //give the charge created per micro pulse and bin and the pulse duration
+private:
+    double mChargePerPulse;
+
+
+};
+
+#endif // BEAM_MODEL

include/chamber.hpp

@@ -0,0 +1,28 @@
+#ifndef CHAMBER
+#define CHAMBER
+/*
+ Takes care of loading the chamber file and returning the chamber parameters pertinent to the calculation
+*/
+
+
+#include <fstream>
+
+const static double pi = 3.14159265358979323846;
+
+class chamber
+{
+	private:
+        double plateDistance;
+        double diameter;
+        double kQ;
+        double Nw;
+	public:
+		chamber(std::ifstream& );
+        double calc_bin_size(int numberOfBins) {return (plateDistance*1.0e6/numberOfBins);} //in nm
+        double calc_E0(double voltage) { return (voltage*10./plateDistance);} //in [Voltage]/cm
+        double calc_charge_liberated(double dose) {return (1.0e9*dose/(Nw*kQ));} // in pC for a dose given in mGy
+        double get_area() { return (pi/4.*diameter*diameter);} //in mm^2
+
+};
+
+#endif

include/functions.hpp

@@ -0,0 +1,40 @@
+#ifndef FUNCTIONS_H
+#define FUNCTIONS_H
+/*
+ assorted helper functions, so far only the reading of the input paramters is here
+*/
+
+
+#include "active_medium.hpp"
+#include "beam_model.hpp"
+#include <stdexcept>
+
+//struct to return all the user input variables from parse_cmd_parameters to main
+struct variable_package{
+    double dose;
+    double voltage;
+    double pulseduration;
+    int numberOfBins;
+    double binSize;
+    double rateLimit;
+    double speedLimit;
+    double electrodeArea;
+    double Q;
+    double E0;
+    streambuf* logBuffer;
+    string summaryFileName;
+    int summaryIndepVar;
+    string mediumName;
+    active_medium* fillgas;
+    beam_model* this_beam;
+    string settingsLogFileName;
+};
+
+//exceptions calls to signal an abort
+class aborted: public std::exception
+{};
+
+//return the parsed command line arguments as a struct of the variable_package form
+variable_package parse_cmd_parameters(int, char**, ofstream &pLogStream, string);
+
+#endif // FUNCTIONS_H

include/msvc_preheader.hpp

@@ -0,0 +1 @@
+#pragma execution_character_set("UTF-8")

include/settings.hpp

@@ -0,0 +1,33 @@
+#ifndef SETTINGS
+#define SETTINGS
+
+//use to define some global preprocessor switches
+
+//turn off the recalculation of the electric field, i.e. deactivate shield by charges
+//#define CONSTFIELD
+
+// logging levels are 0, 1 and 2
+// at 0 only the inital settings and the results are outputted
+// at 1 a status every 1000 steps is added
+// at 2 the calculation of the stepsize is logged at every iteration
+// in order to allow a definition of DEBUGLEVEL at the compilation command line this is only defined if it has not happened already
+#ifndef DEBUGLEVEL
+    #define DEBUGLEVEL 2
+#endif
+
+//output files containing E-Field, electron and ion distributions, every 1000 steps
+//#define CHARGETIMELINE
+
+//undefine saves computations compared to just setting electron recombination coefficient=0
+//#define ELECTRONCOMBO
+
+
+// define if using constant ion mobilities, constant ion recombination rate or constant electron recombination rate (no E-field dependence)
+// significantly speeds up the calculation
+#define CONSTIONMOB
+#define CONSTRECOMB
+#define CONSTELECTRONCOMBO
+
+
+
+#endif // SETTINGS

input/AdvMarkus_kammer.txt

@@ -0,0 +1,5 @@
+Dicke: 1.00 mm
+Durchmesser: 5.0 mm
+Nw: 1.385e+09 Gy/C
+kQ: 1.0 Strahlqualität*Tiefenkorrektur
+