ausa/eloss/ allows you to calculate energy losses for practically any ion-target combination using several different tabulations of stopping powers and ranges. The basic structure of the program is shown the figure below. Tabulations of stopping power and range may be generated for any ion-target combination using one of three sources: SRIM, GEANT or ICRU. The tabulations are fed to the `StoppingPowerCalculator` and `RangeCalculator` which calculate the stopping power and range for an arbitrary energy, typically, by cubic spline interpolation. Finally, the `EnergLossCalculator` calculates the energy loss for an arbitrary energy (initial or final) and thickness, either by integrating the stopping powers or by inverting the range function. For complete documentation of the code see http://docs.kern.phys.au.dk/.
First you must tell AUSAlib where to look for the tabulated data. This is most easily accomplished by adding the following two lines to the .profile file (or .bashrc or similar) in your home directory
The range calculator, on the other hand, has two different implementations (cubic spline interpolation and B-spline least-squares fit), which may be initialized in the following way:
In case you are confused - what the above chunks of code do is to create three calculators and provide us with (unique) pointers to them called `spc`, `rc1`, and `rc2`.
## 4. Energy-loss calculators
Next, we create three energy-loss calculators:
```cpp
#include<ausa/eloss/GslIntegrator.h>
#include<ausa/eloss/EnergyLossRangeInverter.h>
GslIntegratorDEsp(move(spc));
EnergyLossRangeInverterDEr1(move(rc1));
EnergyLossRangeInverterDEr2(move(rc2));
```
The calculators provide identical functionalities, but rely on different methods. The first calculator (`DEsp`) integrates the stopping powers numerically, whereas the second and the third calculator (`DEr1` and `DEr2`) use the inverse-range method, though they differ in the way they determine the range function.
If you do not have [GSL](https://git.kern.phys.au.dk/ausa/ausalib/wikis/GSL) installed on your computer, you can still use our home-made integrator:
```cpp
#include<ausa/eloss/EnergyLossIntegrator.h>
EnergyLossIntegratorDEsp(move(spc));
```
which, however, is not as fast as the GSL integrator.
## 5. Calculate!
We are finally ready to do some calculations. As an example, let us calculate the energy loss of 3 MeV alpha particles in 100 nm silicon by integrating the stopping powers numerically:
```cpp
doublede=DEsp.getTotalEnergyLoss(3000.,100E-6);
```
The result should be close to 100 keV, the exact value depending on the tabulation you are using.
More examples can be found [here](uploads/66bc7fb78ed597a789249a0e2870cc03/snippet7.cpp).
For a complete list of available get methods see [EnergyLossCalculator.h](http://docs.kern.phys.au.dk/classAUSA_1_1EnergyLoss_1_1EnergyLossCalculator.html).
## 6. Default calculator
If you are in a rush and do not care too much about how the energy loss is being determined, you can skip steps 2-4 and instead use the [default energy-loss calculator](https://git.kern.phys.au.dk/ausa/ausalib/blob/master/include/ausa/eloss/Default.h) which determines energy loss by numerical integration of GEANT stopping powers. To determine, say, the energy loss of 1 MeV protons in 100 nm aluminium, you would write
You can chose between three different tabulations: SRIM, GEANT, and ICRU. The quantities provided by each tabulation are listed in the table below. Figure 3 below gives an idea of the level of agreement that can be expected.
SRIM provides both electronic and nuclear stopping powers and is the only one that also provides ranges. Each ion-target combination requires its own unique SRIM file. AUSAlib only offers a very limited selection of SRIM files, which may be found in `res/eloss/SRIM08/`. If you need a SRIM file that is not in AUSAlib you must generate it yourself using the [SRIM program](http://www.srim.org/). Note that `EnergyLossCalculator` has a method called `setScalingFactor` that allows you to scale all stopping powers and/or ranges by a common factor. This is particularly convenient if you are using SRIM because it allows you to change the density of the stopping material without having to generate a new SRIM file.
## GEANT
GEANT provides both electronic and nuclear stopping powers but not ranges. Stopping powers are given in units of MeV*cm<sup>2</sup>/mg. They are assumed to be independent of the ion mass and are tabulated as a function of the energy per nucleon. AUSAlib contains files for all ion-target combinations from hydrogen to uranium (92 x 92 = 8464 files). The files may be found in `res/eloss/GEANT4/`. They have been generated using a slightly modified version of [TestEm0](http://geant4.web.cern.ch/geant4/UserDocumentation/Doxygen/examples_doc/html/ExampleTestEm0.html) and assuming [standard atomic weights](http://www.nist.gov/pml/data/comp.cfm).
Within an isotopic chain the stopping power is assumed to scale inversely proportional to the atomic weight.
For composite target materials the stopping power is calculated as the sum of the individual stopping powers weighted by the [mass fraction](http://en.wikipedia.org/wiki/Mass_fraction_%28chemistry%29).
## ICRU
ICRU only provides electronic stopping powers. AUSAlib contains files for most but not all ion-target combinations. For H and He the stopping powers are the ones used by the programs [http://physics.nist.gov/PhysRefData/Star/Text/PSTAR.html](PSTAR) and [ASTAR](http://physics.nist.gov/PhysRefData/Star/Text/ASTAR.html) and are (as far as I understand) the ones recommended in ICRU report 49.
The files have been obtained from the [libdedx](http://sourceforge.net/projects/libdedx/) repository and can be found inside AUSAlib at `res/eloss/ICRU/`. The stopping powers of Li and heavier ions have been generated using the program PASS, based on the binary theory of Sigmund and Schinner (NIMB 195 (2002) 64 and EPJD 12 (2000) 425), and are the ones recommended in ICRU report 73. The files have been obtained from the GEANT4 database on low-energy electromagnetic processes (G4EMLOW6.32). Like GEANT stopping powers, ICRU stopping powers are independent of ion mass, scale inversely proportional to the atomic weight, and are additive.
# Comparison of calculators
Currently, you can choose between two calculators, the __integrator__ and the __range inverter__. The integrator exists in a GSL version and a home-made version. The range function may be determined using either cubic spline interpolation or B-spline least-squares fitting. So we actually have four different methods to calculate energy losses. The strengths and weaknesses of these methods are discussed below.
Comparison of CPU times (1e6 repetitions):
| Process | Integrator (GSL) | Integrator (home-made) | Range inverter (interpolation) | Range inverter (least-squares fit) |
| 3 MeV alphas in 500 nm Si | 0.34 | 0.58 | 0.06 | 0.38 |
| 6 MeV alphas in 10 um Si | 0.37 | 0.69 | 0.04 | 0.38 |
See also [this comparison](https://git.kern.phys.au.dk/ausa/ausalib/issues/30) of the two integrators.
## Integrator
The main advantages of the integrator are that it is numerically reliable and it yields results for all three tabulations, including GEANT which is readily available for any ion-target combinations. Another advantage is that it can separate the energy loss into its electronic and nuclear components. The main drawbacks are that it is slow (the thicker the material, the longer it takes) and that it systematically underestimates the energy loss of light ions at low energies because it does not take into account angular straggling (see Figure 1 below). It is possible to speed up the integrator by adjusting various parameters such as the minimum step size. The complete list of adjustable parameters can be found in [GslIntegrator.h](https://git.kern.phys.au.dk/ausa/ausalib/blob/gsl-integrator/include/ausa/eloss/GslIntegrator.h)(GSL integrator) and [EnergyLossIntegrator.h](https://git.kern.phys.au.dk/ausa/ausalib/blob/gsl-integrator/include/ausa/eloss/EnergyLossIntegrator.h)(home-made integrator).
## Range inverter
The main advantages of the integrator are that it is fast and does not suffer from the above mentioned systematic error. The main drawbacks are that it only yields results for SRIM and that it has some numerical problems around the maximum of the stopping power curve (see Figures 1 and 2 below). Also it does not allow separation of the energy loss into electronic and nuclear components. The numerical problems can be improved but not fully eliminated, at the cost of computational speed, by using B-spline least-squares fitting instead of cubic spline interpolation to determine the range function. The implementation of the B-spline least-square fitting has been taken from [Geometric Tools Engine](http://www.geometrictools.com/Source/CurvesSurfacesVolumes.html)(see also [license agreement](http://www.geometrictools.com/License/Boost/LICENSE_1_0.txt)) while the implementation of the cubic spline interpolation has been taken from [ROOT](http://project-mathlibs.web.cern.ch/project-mathlibs/sw/5_18_00/html/group__Interpolation.html).
