diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml
index c61f828e4..4f49efe66 100644
--- a/.github/workflows/tests.yml
+++ b/.github/workflows/tests.yml
@@ -1,5 +1,5 @@
 # This is a basic workflow to help you get started with Actions
-name: Testing
+name: Tests
 # Controls when the action will run. 
 on: [push, pull_request, workflow_dispatch]  
 jobs:
diff --git a/CITATION.cff b/CITATION.cff
new file mode 100644
index 000000000..4ff97e53f
--- /dev/null
+++ b/CITATION.cff
@@ -0,0 +1,123 @@
+cff-version: 1.2.0
+message: "If you use matRad, consider citing it and the corresponding research paper"
+authors:
+- family-names: "Ackermann"
+  given-names: "Benjamin"
+- family-names: "Bangert"
+  given-names: "Mark"
+- family-names: "Bennan"
+  given-names: "Amit Ben Antony"
+- family-names: "Burigo"
+  given-names: "Lucas"
+- family-names: "Cabal"
+  given-names: "Gonzalo"
+- family-names: "Cisternas"
+  given-names: "Eduardo"
+- family-names: "Charton"
+  given-names: "Louis"
+- family-names: "Christiansen"
+  given-names: "Eric"
+- family-names: "Ecker"
+  given-names: "Swantje"
+- family-names: "Ellerbrock"
+  given-names: "Malte"
+- family-names: "Gabryś"
+  given-names: "Hubert"
+- family-names: "Handrack"
+  given-names: "Josefine"
+- family-names: "Heath"
+  given-names: "Emily"
+- family-names: "Hermann"
+  given-names: "Cindy"
+- family-names: "Homolka"
+  given-names: "Noa"
+- family-names: "Jäkel"
+  given-names: "Oliver"
+- family-names: "Klinge"
+  given-names: "Thomas"
+- family-names: "Mairani"
+  given-names: "Andrea"
+- family-names: "Mescher"
+  given-names: "Hennig"
+- family-names: "Müller"
+  given-names: "Lucas-Raphael"
+- family-names: "Neishabouri"
+  given-names: "Ahmad"
+- family-names: "Palkowitsch"
+  given-names: "Martina"
+- family-names: "Parodi"
+  given-names: "Katia"
+- family-names: "Pezzano"
+  given-names: "Giuseppe"
+- family-names: "Ramirez"
+  given-names: "Daniel"
+- family-names: "Stadler"
+  given-names: "Alexander"
+- family-names: "Ulrich"
+  given-names: "Silke"
+- family-names: "Wahl"
+  given-names: "Niklas"
+  orcid: "https://orcid.org/0000-0002-1451-223X"
+- family-names: "Welsch"
+  given-names: "Jona"
+- family-names: "Wieser"
+  given-names: "Hans-Peter"
+- family-names: "Winter"
+  given-names: "Johanna"
+- family-names: "Xu"
+  given-names: "Tong" 
+
+title: "matRad"
+version: 2.10.1
+doi: 10.5281/zenodo.3879616
+date-released: 2020-06-05
+url: "http://www.matRad.org"
+location: "Heidelberg"
+organization: "Deutsches Krebsforschungszentrum"
+preferred-citation:
+  type: article
+  authors:
+  - family-names: "Wieser"
+    given-names: "Hans-Peter"
+  - family-names: "Cisternas"
+    given-names: "Eduardo"
+  - family-names: "Wahl"
+    given-names: "Niklas"
+    orcid: "https://orcid.org/0000-0002-1451-223X"
+  - family-names: "Ulrich"
+    given-names: "Silke"
+  - family-names: "Stadler"
+    given-names: "Alexander"
+  - family-names: "Mescher"
+    given-names: "Hennig"
+  - family-names: "Müller"
+    given-names: "Lucas-Raphael"
+  - family-names: "Klinge"
+    given-names: "Thomas"
+  - family-names: "Gabryś"
+    given-names: "Hubert"
+  - family-names: "Burigo"
+    given-names: "Lucas"
+  - family-names: "Mairani"
+    given-names: "Andrea"
+  - family-names: "Ecker"
+    given-names: "Swantje"
+  - family-names: "Ackermann"
+    given-names: "Benjamin"
+  - family-names: "Ellerbrock"
+    given-names: "Malte"
+  - family-names: "Parodi"
+    given-names: "Katia"
+  - family-names: "Jäkel"
+    given-names: "Oliver"
+  - family-names: "Bangert"
+    given-names: "Mark"
+  
+  doi: 10.1002/mp.12251
+  journal: "Medical Physics"
+  start: 2556 
+  end: 2568
+  title: "Development of the open-source dose calculation and optimization toolkit matRad"
+  issue: 6
+  volume: 46
+  year: 2017
diff --git a/MatRad_Config.m b/MatRad_Config.m
index 89fdaf1e0..564e5da50 100644
--- a/MatRad_Config.m
+++ b/MatRad_Config.m
@@ -179,8 +179,6 @@ function setDefaultProperties(obj)
             obj.propOpt.defaultRunDAO = 0;
             obj.propOpt.defaultRunSequencing = 0;
 
-            obj.propMC.defaultCarbonEnergySpread = 0; %[%]
-
             obj.propMC.ompMC_defaultHistories = 1e6;
             obj.propMC.ompMC_defaultOutputVariance = false;
 
diff --git a/MatRad_MCemittanceBaseData.m b/MatRad_MCemittanceBaseData.m
index 024ad3c3e..ac03e390b 100644
--- a/MatRad_MCemittanceBaseData.m
+++ b/MatRad_MCemittanceBaseData.m
@@ -23,16 +23,16 @@
     % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     
     properties
-        machine         %matRad base data machine struct
-        bdl_path = ''   %stores path to generated file
-        nozzleToIso     %Nozzle to Isocenter Distance
-        smx             %Scanning magnet X to isocenter Distance
-        smy             %Scanning magnet y to isocenter Distance
-        monteCarloData  %MC Phase space data struct
-        selectedFocus   %array containing selected focus indices per energy
-        energyspread    %custom energy spread
-        matRad_cfg      %matRad config
-        rangeShifters   %Stores range shifters
+        machine                             %matRad base data machine struct
+        bdl_path = ''                       %stores path to generated file
+        nozzleToIso                         %Nozzle to Isocenter Distance
+        smx                                 %Scanning magnet X to isocenter Distance
+        smy                                 %Scanning magnet y to isocenter Distance
+        monteCarloData                      %MC Phase space data struct
+        selectedFocus                       %array containing selected focus indices per energy
+        defaultRelativeEnergySpread = 0;    %default energy spread
+        matRad_cfg                          %matRad config
+        rangeShifters                       %Stores range shifters
     end
     
     properties (SetAccess = private)
@@ -58,7 +58,6 @@
             end
             
             obj.matRad_cfg = MatRad_Config.instance();
-            obj.energyspread = obj.matRad_cfg.propMC.defaultCarbonEnergySpread;
                         
             obj.machine = machine;
             obj.problemSigma = false;
@@ -258,7 +257,9 @@
                     % https://www.nist.gov/system/files/documents/2017/04/26/newstar.pdf
                     meanEnergy = @(x) 11.39 * x^0.628 + 11.24;
                     mcDataEnergy.MeanEnergy = meanEnergy(r80);
-                    mcDataEnergy.EnergySpread = obj.energyspread;
+                    % reading in a potential given energyspread could go here directly. How would you parse the energyspread
+                    % into the function? Through a field in the machine?
+                    mcDataEnergy.EnergySpread = obj.defaultRelativeEnergySpread;
                 case 'helium'
                     % Fit to Range-Energy relationship
                     % Data from "Update to ESTAR, PSTAR, and ASTAR Databases" - ICRU Report 90, 2014
@@ -266,7 +267,7 @@
                     % https://www.nist.gov/system/files/documents/2017/04/26/newstar.pdf
                     meanEnergy = @(x) 7.57* x.^0.5848 + 3.063;
                     mcDataEnergy.MeanEnergy = meanEnergy(r80);
-                    mcDataEnergy.EnergySpread = obj.energyspread;
+                    mcDataEnergy.EnergySpread = obj.defaultRelativeEnergySpread;
                 otherwise
                     error('not implemented')
             end
diff --git a/README.md b/README.md
index 7eb6d73cb..48b827919 100644
--- a/README.md
+++ b/README.md
@@ -3,8 +3,8 @@
 [![Contributors](https://img.shields.io/github/contributors/e0404/matRad)](https://github.com/e0404/matRad/graphs/contributors)
 [![GitPitch](https://gitpitch.com/assets/badge.svg)](https://gitpitch.com/e0404/matRad/master) 
 
-[![TravisCI Build Status](https://travis-ci.org/e0404/matRad.svg?branch=dev_varRBErobOpt)](https://travis-ci.org/e0404/matRad)
-[![Azure Pipelines Build Status](https://dev.azure.com/e0404/matRad/_apis/build/status/e0404.matRad?branchName=dev_varRBErobOpt)](https://dev.azure.com/e0404/matRad)
+[![GitHub Build Status](https://github.com/e0404/matRad/actions/workflows/tests.yml/badge.svg)](https://github.com/e0404/matRad/actions/workflows/tests.yml)
+[![Azure Pipelines Build Status](https://dev.azure.com/e0404/matRad/_apis/build/status/e0404.matRad)](https://dev.azure.com/e0404/matRad)
 
 DOIs:
 - General DOI: [![DOI](https://zenodo.org/badge/doi/10.5281/zenodo.3879615.svg)](https://doi.org/10.5281/zenodo.3879615)
@@ -111,7 +111,7 @@ Copyright 2020 the matRad development team.
 matrad@dkfz.de
 
 All the elements of the compilation of matRad and Ipopt are free software. You can redistribute and/or modify matRad's source code version provided as files with .m and .mat extension under the terms of the GNU GENERAL PUBLIC LICENSE Version 3 (GPL v3). You can also add to matRad the Ipopt functionality by using the precompiled mex files of the Ipopt optimizer in object code version which are licensed under the Eclipse Public License Version 1.0 (EPL v1.0), also made available for download via https://projects.coin-or.org/Ipopt.
-matRad also contains interfaces to an open-source photon Monte Carlo dose calculation engine developed by Edgardo Doerner hosted on GitHub (http://github.com/edoerner/ompMC) and to the open-source proton Monte Carlo project MCsquare (www.openmcsquare.org) from UCLouvain, Louvain-la-Neuve, Belgium. Both interfaces are integrated into matRad as submodules.
+matRad also contains interfaces to an open-source photon Monte Carlo dose calculation engine developed by Edgardo Dörner hosted on GitHub (http://github.com/edoerner/ompMC) and to the open-source proton Monte Carlo project MCsquare (www.openmcsquare.org) from UCLouvain, Louvain-la-Neuve, Belgium. Both interfaces are integrated into matRad as submodules.
 
 In addition, we provide a matlab standalone version of the compilation of matRad and Ipopt, where the files of matRad and Ipopt are licensed under GPL v3 and EPL v1.0 respectively. The matlab standalone version is meant to be used by students for learning and practicing scientific programming and does not yet contain the interfaces to the aforementioned Monte Carlo dose calculation engines.
 
diff --git a/dicom/matRad_importDicom.m b/dicom/matRad_importDicom.m
index 9a7301aaf..a8a3dac68 100644
--- a/dicom/matRad_importDicom.m
+++ b/dicom/matRad_importDicom.m
@@ -141,10 +141,8 @@
     if ~(cellfun(@isempty,files.rtdose(1,1:2))) 
         fprintf('loading Dose files \n', structures(i).structName);
         % parse plan in order to scale dose cubes to a fraction based dose
-        if exist('pln','var')
-            if isfield(pln,'numOfFractions')
-                resultGUI = matRad_importDicomRTDose(ct, files.rtdose, pln);
-            end
+        if exist('pln','var') && ~isempty(pln) && isfield(pln,'numOfFractions')
+            resultGUI = matRad_importDicomRTDose(ct, files.rtdose, pln);
         else
             resultGUI = matRad_importDicomRTDose(ct, files.rtdose);
         end
diff --git a/matRad_calcCubes.m b/matRad_calcCubes.m
index 8908d9890..466fcca17 100644
--- a/matRad_calcCubes.m
+++ b/matRad_calcCubes.m
@@ -146,24 +146,6 @@
 end
 
 %% Final processing
-% Add non-processed MC tallies
-% Note that the tallies are already computed per beam and altogether
-if isfield(dij,'MC_tallies')
-    for f = 1:numel(dij.MC_tallies)
-        tally = dij.MC_tallies{f};
-        % skip tallies processed above
-        if ~isfield(resultGUI,tally) && isempty(strfind(lower(tally),'std'))
-            tallyCut = strsplit(tally,'_');
-            if size(tallyCut,2) > 1
-                beamNum = str2num(cell2mat(regexp(tally,'\d','Match')));
-                resultGUI.(tally) = reshape(full(dij.(tallyCut{1}){scenNum}(:,beamNum)),dij.doseGrid.dimensions);
-            else
-                resultGUI.(tally) = reshape(full(dij.(tally){scenNum} * wBeam),dij.doseGrid.dimensions);
-            end
-        end
-    end
-end
-
 % Remove suffix for RBExD if there's only one available
 if any(cellfun(@(teststr) ~isempty(strfind(lower(teststr),'alpha')), fieldnames(dij))) && isfield(dij,'RBE_models') && length(dij.RBE_models) == 1
     % Get fieldnames that include the specified RBE model
diff --git a/matRad_calcDoseDirectMC.m b/matRad_calcDoseDirectMC.m
index 3068a085d..1d33b8ca1 100644
--- a/matRad_calcDoseDirectMC.m
+++ b/matRad_calcDoseDirectMC.m
@@ -87,7 +87,7 @@
 
 % calc resulting dose
 if ~isprop(pln.propMC,'externalCalculation') || ~pln.propMC.externalCalculation
-    if pln.multScen.totNumScen == 1
+    if pln.multScen.numOfCtScen == 1
         % calculate cubes; use uniform weights here, weighting with actual fluence
         % already performed in dij construction
         if size(dij.physicalDose{1},2) ~= dij.numOfBeams || size(dij.physicalDose{1},2) ~= numel(dij.beamNum)
@@ -114,9 +114,9 @@
 
     end
 
-    if pln.multScen.totNumScen ~= 1
+    if pln.multScen.numOfCtScen ~= 1
         resultGUI.accPhysicalDose = zeros(size(resultGUI.phaseDose{1}));
-        for i = 1:pln.multScen.totNumScen
+        for i = 1:pln.multScen.numOfCtScen
             resultGUI.accPhysicalDose = resultGUI.accPhysicalDose + resultGUI.phaseDose{i};
         end
     end
diff --git a/matRad_calcParticleDoseMCsquare.m b/matRad_calcParticleDoseMCsquare.m
index 89d7489c7..d7e0024ff 100644
--- a/matRad_calcParticleDoseMCsquare.m
+++ b/matRad_calcParticleDoseMCsquare.m
@@ -398,8 +398,6 @@
                             dij.doseGrid.dimensions, ...
                             dij.totalNumOfBixels, ...
                             mask);
-                        
-                        dij.MC_tallies{1} = 'LET';
                     end
                 else
                     %Read dose cube
@@ -416,20 +414,17 @@
                             sparse(VdoseGrid,ones(numel(VdoseGrid),1), ...
                                 cube(VdoseGrid), ...
                                 dij.doseGrid.numOfVoxels,1);
-                        
-                        dij.MC_tallies{1} = 'LET';
                     end
 
-                    % Postprocessing for dij: Remove fields since this is already the combined dose over all bixels
-                    fields = {'MCsquareCalcOrder','bixelNum','rayNum','totalNumOfRays','totalNumOfBixels','numOfRaysPerBeam','numOfBeams'};
-                    for f = 1:length(fields)
-                        if isfield(dij,fields{f})
-                            dij = rmfield(dij,fields{f});
-                        end
-                    end
+                    % Postprocessing for dij:
+                    % This is already the combined dose over all bixels, so all parameters are 1 in this case
+                    dij = rmfield(dij,'MCsquareCalcOrder');
+
                     dij.numOfBeams = 1;
                     dij.beamNum = 1;
                     dij.bixelNum = 1;
+                    dij.rayNum = 1;
+                    dij.totalNumOfBixels = 1;
                     dij.totalNumOfRays = 1;
                     dij.numOfRaysPerBeam = 1;
                 end
@@ -470,6 +465,9 @@
     end   
 end
 
+% Order fields for easier comparison between different dijs
+dij = orderfields(dij);
+
 %% cd back
 cd(currFolder);
 
diff --git a/matRad_calcParticleDoseMCtopas.m b/matRad_calcParticleDoseMCtopas.m
index 46c1147e0..876d8a723 100644
--- a/matRad_calcParticleDoseMCtopas.m
+++ b/matRad_calcParticleDoseMCtopas.m
@@ -1,4 +1,4 @@
-function topasCubes = matRad_calcParticleDoseMCtopas(ct,stf,pln,cst,calcDoseDirect)
+function dij = matRad_calcParticleDoseMCtopas(ct,stf,pln,cst,calcDoseDirect)
 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 % matRad TOPAS Monte Carlo proton dose calculation wrapper
 %   This calls a TOPAS installation (not included in matRad due to
@@ -91,7 +91,7 @@
 matRad_calcDoseInit;
 
 % for TOPAS we explicitly downsample the ct to the dose grid (might not be necessary in future versions with separated grids)
-[ctR,~,~] = pln.propMC.resampleGrid(ct,cst,pln,stf);
+[ctR,~,~] = matRad_resampleCTtoGrid(ct,cst,pln,stf);
 
 % overwrite CT grid in dij in case of modulation.
 if isfield(ctR,'ctGrid')
@@ -167,9 +167,11 @@
     % change director back to original directory
     cd(pln.propMC.workingDir);
 
-    % save dij and weights, they are needed for later reading the data back in
+    % Skip local calculation and data readout with this parameter. All necessary parameters to read the data back in
+    % later are stored in the MCparam file that is stored in the folder. The folder is generated in the working
+    % directory and the matRad_plan*.txt file can be manually called with TOPAS.
     if pln.propMC.externalCalculation
-        matRad_cfg.dispInfo('TOPAS simulation skipped for external calculation\n');
+        matRad_cfg.dispInfo(['TOPAS simulation skipped for external calculation\nFiles have been written to: "',replace(pln.propMC.workingDir,'\','\\'),'"']);
     else
         for ctScen = 1:ct.numOfCtScen
             for beamIx = 1:numel(stf)
@@ -233,11 +235,14 @@
 %% Simulation(s) finished - read out volume scorers from topas simulation
 % Skip readout if external files were generated
 if ~pln.propMC.externalCalculation
-    topasCubes = pln.propMC.readFiles(pln.propMC.workingDir);
+    dij = pln.propMC.readFiles(pln.propMC.workingDir);
 else
-    topasCubes = [];
+    dij = [];
 end
 
+% Order fields for easier comparison between different dijs
+dij = orderfields(dij);
+
 % manipulate isocenter back
 for k = 1:length(stf)
     stf(k).isoCenter = stf(k).isoCenter - pln.multScen.isoShift(shiftScen,:);
diff --git a/matRad_calcPhotonDose.m b/matRad_calcPhotonDose.m
index 0000a96e3..b60ffadaf 100644
--- a/matRad_calcPhotonDose.m
+++ b/matRad_calcPhotonDose.m
@@ -320,7 +320,7 @@
                         
                         % sample dose only for bixel based dose calculation
                         if ~isFieldBasedDoseCalc
-                            r0   = 25;   % [mm] sample beyond the inner core
+                            r0   = 20 + stf(i).bixelWidth;   % [mm] sample beyond the inner core
                             Type = 'radius';
                             [ix,bixelDose] = matRad_DijSampling(ix,bixelDose,manipulatedRadDepthCube,rad_distancesSq,Type,r0);
                         end
diff --git a/matRad_fluenceOptimization.m b/matRad_fluenceOptimization.m
index ac3e977de..1d6107715 100644
--- a/matRad_fluenceOptimization.m
+++ b/matRad_fluenceOptimization.m
@@ -108,16 +108,7 @@
     if ~isfield(dij,'RBE')
         dij.RBE = 1.1;
     end
-
-    if exist('wInit','var')
-        matRad_cfg.dispInfo('Initial weights found as input...\n');
-    else
-        doseTmp = dij.physicalDose{1}*wOnes;
-        bixelWeight =  (doseTarget)/(dij.RBE * mean(doseTmp(V)));
-        wInit       = wOnes * bixelWeight;
-        matRad_cfg.dispInfo('chosen uniform weight of %f!\n',bixelWeight);
-    end
-
+    
     doseTmp = dij.physicalDose{1}*wOnes;
     bixelWeight =  (doseTarget)/(dij.RBE * mean(doseTmp(V)));
     wInit       = wOnes * bixelWeight;
@@ -147,17 +138,13 @@
 
     if isequal(pln.bioParam.quantityOpt,'effect')
 
-        if exist('wInitIn','var')
-            matRad_cfg.dispInfo('Initial weights found as input...\n');
-        else
-            effectTarget = cst{ixTarget,5}.alphaX * doseTarget + cst{ixTarget,5}.betaX * doseTarget^2;
-            aTmp = dij.mAlphaDose{1}*wOnes;
-            bTmp = dij.mSqrtBetaDose{1} * wOnes;
-            p = sum(aTmp(V)) / sum(bTmp(V).^2);
-            q = -(effectTarget * length(V)) / sum(bTmp(V).^2);
+        effectTarget = cst{ixTarget,5}.alphaX * doseTarget + cst{ixTarget,5}.betaX * doseTarget^2;
+        aTmp = dij.mAlphaDose{1}*wOnes;
+        bTmp = dij.mSqrtBetaDose{1} * wOnes;
+        p = sum(aTmp(V)) / sum(bTmp(V).^2);
+        q = -(effectTarget * length(V)) / sum(bTmp(V).^2);
 
-            wInit        = -(p/2) + sqrt((p^2)/4 -q) * wOnes;
-        end
+        wInit        = -(p/2) + sqrt((p^2)/4 -q) * wOnes;
 
     elseif isequal(pln.bioParam.quantityOpt,'RBExD')
 
@@ -165,34 +152,26 @@
         dij.gamma             = zeros(dij.doseGrid.numOfVoxels,dij.numOfScenarios);
         dij.gamma(dij.ixDose) = dij.ax(dij.ixDose)./(2*dij.bx(dij.ixDose));
 
-        if exist('wInit','var')
-            matRad_cfg.dispInfo('Initial weights found as input...\n');
-        else
-            % calculate current effect in target
-            aTmp = dij.mAlphaDose{1}*wOnes;
-            bTmp = dij.mSqrtBetaDose{1} * wOnes;
-            doseTmp = dij.physicalDose{1}*wOnes;
-
-            CurrEffectTarget = aTmp(V) + bTmp(V).^2;
-            % ensure a underestimated biological effective dose
-            TolEstBio        = 1.2;
-            % calculate maximal RBE in target
-            maxCurrRBE = max(-cst{ixTarget,5}.alphaX + sqrt(cst{ixTarget,5}.alphaX^2 + ...
-                4*cst{ixTarget,5}.betaX.*CurrEffectTarget)./(2*cst{ixTarget,5}.betaX*doseTmp(V)));
-            wInit    =  ((doseTarget)/(TolEstBio*maxCurrRBE*max(doseTmp(V))))* wOnes;
-        end
+        % calculate current effect in target
+        aTmp = dij.mAlphaDose{1}*wOnes;
+        bTmp = dij.mSqrtBetaDose{1} * wOnes;
+        doseTmp = dij.physicalDose{1}*wOnes;
+
+        CurrEffectTarget = aTmp(V) + bTmp(V).^2;
+        % ensure a underestimated biological effective dose
+        TolEstBio        = 1.2;
+        % calculate maximal RBE in target
+        maxCurrRBE = max(-cst{ixTarget,5}.alphaX + sqrt(cst{ixTarget,5}.alphaX^2 + ...
+            4*cst{ixTarget,5}.betaX.*CurrEffectTarget)./(2*cst{ixTarget,5}.betaX*doseTmp(V)));
+        wInit    =  ((doseTarget)/(TolEstBio*maxCurrRBE*max(doseTmp(V))))* wOnes;
     end
 
     matRad_cfg.dispInfo('chosen weights adapted to biological dose calculation!\n');
 else
-    if exist('wInit','var')
-        matRad_cfg.dispInfo('Initial weights found as input...\n');
-    else
-        doseTmp = dij.physicalDose{1}*wOnes;
-        bixelWeight =  (doseTarget)/mean(doseTmp(V));
-        wInit       = wOnes * bixelWeight;
-        matRad_cfg.dispInfo('chosen uniform weight of %f!\n',bixelWeight);
-    end
+    doseTmp = dij.physicalDose{1}*wOnes;
+    bixelWeight =  (doseTarget)/mean(doseTmp(V));
+    wInit       = wOnes * bixelWeight;
+    matRad_cfg.dispInfo('chosen uniform weight of %f!\n',bixelWeight);
 end
 
 
diff --git a/tools/matRad_resampleCTtoGrid.m b/tools/matRad_resampleCTtoGrid.m
new file mode 100644
index 000000000..d6281ebe5
--- /dev/null
+++ b/tools/matRad_resampleCTtoGrid.m
@@ -0,0 +1,77 @@
+function [ctR,cst,stf] = matRad_resampleCTtoGrid(ct,cst,pln,stf)
+% function to resample the ct grid for example for faster MC computation
+%
+% call
+%   [ctR,cst,stf] = matRad_resampleGrid(ct,cst,stf)
+%
+% input
+%   ct:             Path to folder where TOPAS files are in (as string)
+%   cst:            matRad segmentation struct
+%   pln:            matRad plan struct
+%   stf:            matRad steering struct
+%
+% output
+%   ctR:            resampled CT
+%   cst:            updated ct struct (due to calcDoseInit)
+%   stf:            updated stf struct (due to calcDoseInit)
+%
+% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% Copyright 2022 the matRad development team.
+%
+% This file is part of the matRad project. It is subject to the license
+% terms in the LICENSE file found in the top-level directory of this
+% distribution and at https://github.com/e0404/matRad/LICENSES.txt. No part
+% of the matRad project, including this file, may be copied, modified,
+% propagated, or distributed except according to the terms contained in the
+% LICENSE file.
+%
+% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+matRad_cfg = MatRad_Config.instance(); %Instance of matRad configuration class
+
+% Load calcDoseInit as usual
+% Note of warning, even though the pln is marked as unused, it is needed for calcDoseInit!
+matRad_calcDoseInit;
+
+% Check if CT has already been resampled
+if ~isfield(ct,'resampled')
+    % Allpcate resampled cubes
+    cubeHUresampled = cell(1,ct.numOfCtScen);
+    cubeResampled = cell(1,ct.numOfCtScen);
+
+    % Perform resampling to dose grid
+    for s = 1:ct.numOfCtScen
+        cubeHUresampled{s} =  matRad_interp3(dij.ctGrid.x,  dij.ctGrid.y',  dij.ctGrid.z,ct.cubeHU{s}, ...
+            dij.doseGrid.x,dij.doseGrid.y',dij.doseGrid.z,'linear');
+        cubeResampled{s} =  matRad_interp3(dij.ctGrid.x,  dij.ctGrid.y',  dij.ctGrid.z,ct.cube{s}, ...
+            dij.doseGrid.x,dij.doseGrid.y',dij.doseGrid.z,'linear');
+    end
+
+    % Allocate temporary resampled CT
+    ctR = ct;
+    ctR.cube = cell(1);
+    ctR.cubeHU = cell(1);
+
+    % Set CT resolution to doseGrid resolution
+    ctR.resolution = dij.doseGrid.resolution;
+    ctR.cubeDim = dij.doseGrid.dimensions;
+    ctR.x = dij.doseGrid.x;
+    ctR.y = dij.doseGrid.y;
+    ctR.z = dij.doseGrid.z;
+
+    % Write resampled cubes
+    ctR.cubeHU = cubeHUresampled;
+    ctR.cube = cubeResampled;
+
+    % Set flag for complete resampling
+    ctR.resampled = 1;
+    ctR.ctGrid = dij.doseGrid;
+
+    % Save original grid
+    ctR.originalGrid = dij.ctGrid;
+else
+    ctR = ct;
+    matRad_cfg.dispWarning('CT already resampled.');
+end
+end
diff --git a/topas/MatRad_TopasConfig.m b/topas/MatRad_TopasConfig.m
index ca81408bc..aad30c192 100644
--- a/topas/MatRad_TopasConfig.m
+++ b/topas/MatRad_TopasConfig.m
@@ -304,89 +304,6 @@ function writeAllFiles(obj,ct,cst,pln,stf,machine,w)
             matRad_cfg.dispInfo('Successfully written TOPAS setup files!\n')
         end
 
-        function [ctR,cst,stf] = resampleGrid(~,ct,cst,pln,stf)
-            % function to resample the TOPAS grid for faster computation
-            %
-            % call
-            %   [ctR,cst,stf] = topasConfig.resampleGrid(ct,cst,pln,stf)
-            %   [ctR,cst,stf] = obj.resampleGrid(ct,cst,pln,stf)
-            %
-            % input
-            %   ct:             Path to folder where TOPAS files are in (as string)
-            %   cst:            matRad segmentation struct
-            %   pln:            matRad plan struct
-            %   stf:            matRad steering struct
-            %
-            % output
-            %   ctR:            resampled CT
-            %   cst:            updated ct struct (due to calcDoseInit)
-            %   stf:            updated stf struct (due to calcDoseInit)
-            %
-            % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-            %
-            % Copyright 2022 the matRad development team.
-            %
-            % This file is part of the matRad project. It is subject to the license
-            % terms in the LICENSE file found in the top-level directory of this
-            % distribution and at https://github.com/e0404/matRad/LICENSES.txt. No part
-            % of the matRad project, including this file, may be copied, modified,
-            % propagated, or distributed except according to the terms contained in the
-            % LICENSE file.
-            %
-            % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-            matRad_cfg = MatRad_Config.instance(); %Instance of matRad configuration class
-
-            % Set flag in case a modulated CT is detected
-            if isfield(ct,'modulated') && ct.modulated
-                pln.propDoseCalc.useGivenEqDensityCube = true;
-            end
-
-            % Load calcDoseInit as usual
-            matRad_calcDoseInit;
-
-            % Check if CT has already been resampled
-            if ~isfield(ct,'resampled')
-                % Allpcate resampled cubes
-                cubeHUresampled = cell(1,ct.numOfCtScen);
-                cubeResampled = cell(1,ct.numOfCtScen);
-
-                % Perform resampling to dose grid
-                for s = 1:ct.numOfCtScen
-                    cubeHUresampled{s} =  matRad_interp3(dij.ctGrid.x,  dij.ctGrid.y',  dij.ctGrid.z,ct.cubeHU{s}, ...
-                        dij.doseGrid.x,dij.doseGrid.y',dij.doseGrid.z,'linear');
-                    cubeResampled{s} =  matRad_interp3(dij.ctGrid.x,  dij.ctGrid.y',  dij.ctGrid.z,ct.cube{s}, ...
-                        dij.doseGrid.x,dij.doseGrid.y',dij.doseGrid.z,'linear');
-                end
-
-                % Allocate temporary resampled CT
-                ctR = ct;
-                ctR.cube = cell(1);
-                ctR.cubeHU = cell(1);
-
-                % Set CT resolution to doseGrid resolution
-                ctR.resolution = dij.doseGrid.resolution;
-                ctR.cubeDim = dij.doseGrid.dimensions;
-                ctR.x = dij.doseGrid.x;
-                ctR.y = dij.doseGrid.y;
-                ctR.z = dij.doseGrid.z;
-
-                % Write resampled cubes
-                ctR.cubeHU = cubeHUresampled;
-                ctR.cube = cubeResampled;
-
-                % Set flag for complete resampling
-                ctR.resampled = 1;
-                ctR.ctGrid = dij.doseGrid;
-
-                % Save original grid
-                ctR.originalGrid = dij.ctGrid;
-            else
-                ctR = ct;
-                matRad_cfg.dispWarning('CT already resampled.');
-            end
-        end
-
         function dij = readFiles(obj,folder)
             % function to read out TOPAS data
             %
@@ -900,7 +817,7 @@ function writeAllFiles(obj,ct,cst,pln,stf,machine,w)
                         end
                     end
                 end
-
+                
                 % Remove processed physDoseFields from total tallies
                 topasCubesTallies = topasCubesTallies(~physDoseFields);
 
@@ -918,11 +835,7 @@ function writeAllFiles(obj,ct,cst,pln,stf,machine,w)
                                         dij.([topasCubesTallies{j} processedQuantities{p}]){ctScen,1}(:,d) = sum(w)*reshape(topasCubes.([topasCubesTallies{j} '_ray' num2str(dij.rayNum(d)) '_bixel' num2str(dij.bixelNum(d)) processedQuantities{p} '_beam' num2str(dij.beamNum(d))]){ctScen},[],1);
                                     end
                                 end
-                                % Handle RBE-related quantities (not multiplied by sum(w)!)
                             elseif ~isempty(strfind(lower(topasCubesTallies{j}),'alpha'))
-                                modelName = strsplit(topasCubesTallies{j},'_');
-                                modelName = modelName{end};
-                                if isfield(topasCubes,[topasCubesTallies{j} '_ray' num2str(dij.rayNum(d)) '_bixel' num2str(dij.bixelNum(d)) '_beam' num2str(dij.beamNum(d))]) ...
                                         && iscell(topasCubes.([topasCubesTallies{j} '_ray' num2str(dij.rayNum(d)) '_bixel' num2str(dij.bixelNum(d)) '_beam' num2str(dij.beamNum(d))]))
                                     dij.(['mAlphaDose_' modelName]){ctScen,1}(:,d) = reshape(topasCubes.([topasCubesTallies{j} '_ray' num2str(dij.rayNum(d)) '_bixel' num2str(dij.bixelNum(d)) '_beam' num2str(dij.beamNum(d))]){ctScen},[],1) .* dij.physicalDose{ctScen,1}(:,d);
                                 end
diff --git a/topas/linuxScripts/install/README.txt b/topas/linuxScripts/install/README.txt
deleted file mode 100644
index 70e1c40b6..000000000
--- a/topas/linuxScripts/install/README.txt
+++ /dev/null
@@ -1,11 +0,0 @@
-## installTopas.sh
-1. Place all files in the directory you want TOPAS to be installed (both file parts)
-2. You probably have to change execute permissions to the installTopas.sh file (chmod +x installTopas.sh)
-3. call "sudo ./installTopas.sh"
-
-## installTelegram.sh - installs python Linux to communicate back from the server with updates about the simulation
-USAGE: ./installTelegram.sh
-
-## WARNING: There might be an error occuring that is attributed to line endings, so make sure to double check that the line endings are LF and not CRLF.
-
-If you have any questions don't hesitate to write me an email under n.homolka@dkfz-heidelberg.de
\ No newline at end of file
diff --git a/topas/linuxScripts/install/installTelegram.sh b/topas/linuxScripts/install/installTelegram.sh
deleted file mode 100644
index 4934f60d1..000000000
--- a/topas/linuxScripts/install/installTelegram.sh
+++ /dev/null
@@ -1,5 +0,0 @@
-#!/bin/bash
-
-sudo apt install python3-pip
-sudo pip3 install telegram-send
-pip3 install telegram-send
\ No newline at end of file
diff --git a/topas/linuxScripts/install/installTopas.sh b/topas/linuxScripts/install/installTopas.sh
deleted file mode 100644
index 97658173e..000000000
--- a/topas/linuxScripts/install/installTopas.sh
+++ /dev/null
@@ -1,71 +0,0 @@
-#!/bin/bash
-# Install-script for TOPAS 3.7
-# use as follows:
-# put install script in folder together with the topas tar package (only home directory works apparently)
-# run "sudo ./installTopas.sh"
-
-install_path="$( cd "$( dirname "$0" )" && pwd )"
-cat topas_*_debian9.tar.gz.part_* > topas_debian9.tar.gz
-
-# Create folder structure
-mkdir $install_path/topas $install_path/topas_extensions
-
-# Install prerequisites
-apt install -y libexpat1-dev
-apt install -y libgl1-mesa-dev
-apt install -y libglu1-mesa-dev
-apt install -y libxt-dev
-apt install -y xorg-dev
-apt install -y build-essential
-apt install -y libharfbuzz-dev
-
-# Install TOPAS
-apt install -y unzip
-tar -zxvf $install_path/topas_debian9.tar.gz -C $install_path
-
-# Install extensions
-wget -O $install_path/topas_extensions/extensions.zip https://github.com/topasmc/extensions/archive/master.zip
-unzip $install_path/topas_extensions/extensions.zip -d $install_path/topas_extensions/
-mv $install_path/topas_extensions/extensions-master/* $install_path/topas_extensions/
-rm -r $install_path/topas_extensions/extensions-master
-rm $install_path/topas_extensions/extensions.zip
-
-# Install Geant4 Data
-unzip $install_path/topas/Geant4Header*.zip -d $install_path/topas
-mkdir $install_path/G4Data
-cd $install_path/G4Data
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4NDL.4.6.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4EMLOW.7.9.1.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4PhotonEvaporation.5.5.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4RadioactiveDecay.5.4.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4PARTICLEXS.2.1.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4SAIDDATA.2.0.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4ABLA.3.1.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4INCL.1.0.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4PII.1.3.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4ENSDFSTATE.2.2.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4RealSurface.2.1.1.tar.gz
-wget -nc -4 http://geant4-data.web.cern.ch/geant4-data/datasets/G4TENDL.1.3.2.tar.gz
-
-tar -zxf G4NDL.4.6.tar.gz
-tar -zxf G4EMLOW.7.9.1.tar.gz
-tar -zxf G4PhotonEvaporation.5.5.tar.gz
-tar -zxf G4RadioactiveDecay.5.4.tar.gz
-tar -zxf G4PARTICLEXS.2.1.tar.gz
-tar -zxf G4SAIDDATA.2.0.tar.gz
-tar -zxf G4ABLA.3.1.tar.gz
-tar -xzf G4INCL.1.0.tar.gz
-tar -zxf G4PII.1.3.tar.gz
-tar -zxf G4ENSDFSTATE.2.2.tar.gz
-tar -zxf G4RealSurface.2.1.1.tar.gz
-tar -zxf G4TENDL.1.3.2.tar.gz
-
-export TOPAS_G4_GATA_DIR= $install_path/G4Data
-
-# Build
-cd $install_path/topas
-apt install -y cmake
-cmake -DTOPAS_EXTENSIONS_DIR=$install_path/topas_extensions
-make
-
-echo "Installer finished."
diff --git a/topas/linuxScripts/run/README.txt b/topas/linuxScripts/run/README.txt
deleted file mode 100644
index a3db547ab..000000000
--- a/topas/linuxScripts/run/README.txt
+++ /dev/null
@@ -1,8 +0,0 @@
-## calcAll.sh - wrapper script to run calcTopas.sh over multiple folders
-Usage:  ./calcAll.sh path-to-folder1 path-to-folder2 ...
-
-## calcTopas.sh - runs Topas for all subfolders in given folder and sends updates to telegram (see doc
-## for "telegram-send") with the estimated time remaining.
-Usage:  ./calcTopas.sh path-to-folder
-
-If you have any questions don't hesitate to write me an email under n.homolka@dkfz-heidelberg.de
\ No newline at end of file
diff --git a/topas/linuxScripts/run/calcAll.sh b/topas/linuxScripts/run/calcAll.sh
deleted file mode 100644
index 9046422c2..000000000
--- a/topas/linuxScripts/run/calcAll.sh
+++ /dev/null
@@ -1,9 +0,0 @@
-#!/bin/bash
-
-for var in "$@"
-do
-	./calcTopas.sh "$var"
-done
-
-telegram-send --format markdown "*** ${HOSTNAME:11:7}: All queued simulations finished ***"
-echo "All done!"
diff --git a/topas/linuxScripts/run/calcTopas.sh b/topas/linuxScripts/run/calcTopas.sh
deleted file mode 100644
index 35885b554..000000000
--- a/topas/linuxScripts/run/calcTopas.sh
+++ /dev/null
@@ -1,87 +0,0 @@
-#!/bin/bash
-# heterogeneity
-
-export TOPAS_G4_DATA_DIR=~/G4Data
-
-counter=0
-for f in $(find $PWD/$1 -maxdepth 3 -type d | sort -n)
-do
-	res=($f/*run*.txt) #check if folder contains txt file to run
-	if [ -f "$res" ]; then
-                directories[counter]=$f
-		let counter=counter+1
-	fi
-done
-
-echo "Submitted calculation of ${#directories[@]} folders."
-counter=1
-
-for directory in "${directories[@]}"
-do	
-	cd $directory
-	seconds1=$(date '+%s')
-	fileTest=(*run5.txt)
-	if [ -f "$fileTest" ]; then
-		subFiles=(matRad_plan_field*.txt)
-		counterSub=1
-	fi
-
-	for s in $directory/matRad_plan_*.txt
-	do
-		# Start second timer if multiple runs are present (long simulation with substeps)
-		if [ -f "$fileTest" ]; then
-			secondsSub1=$(date '+%s')
-		fi
-		
-		# Start TOPAS calculation of specific txt file
-		echo "Starting calculation for file $(basename "$s") in folder $directory"
-		~/topas/bin/topas ${s%.*}.txt > ${s%.*}.out > ${s%.*}.err
-
-		# check if produced bin file contains values or is empty
-		file=${s##*/}
-		if [ -s $directory/score_${file%.*}_physicalDose.bin ]; then
-        		echo "Calculation successful.\n"
-		else
-			seconds2=$(date '+%s')
-			let difference[counter]=seconds2-seconds1
-        		printf "#######################################\nERROR:\n${HOSTNAME:11:7}: No results in file ${s%.*} of folder '${1:0:-1}'! Finished in ${difference[counter]} seconds.\n#######################################" | telegram-send --stdin
-			continue 2
-		fi
-
-		# output for multiple subfolders/runs to give updates on the length
-		if [ -f "$fileTest" ]; then
-			secondsSub2=$(date '+%s')
-			let differenceSub[counterSub]=secondsSub2-secondsSub1
-			sumSub=$(IFS=+; echo "$((${differenceSub[*]}))")
-			let estimatedTimeSub=(${#subFiles[@]}-counterSub)
-			let estimatedTimeSub=estimatedTimeSub*sumSub/counterSub
-
-			if [ "$estimatedTimeSub" -ne "0" ]; then
-				finishDateSub=$(date '+%a, %d.%m., %H:%Mh' --date="@$((secondsSub2+ estimatedTimeSub))")
-				printf "${HOSTNAME:11:7}: Simulation of subFile $file ($counterSub/${#subFiles[@]}) in '$(basename $directory)' finished in ${differenceSub[counterSub]} seconds.\nEstimated date finished: $finishDateSub." | telegram-send --silent --stdin
-			else
-				printf "${HOSTNAME:11:7}: Simulation of subFile $file ($counterSub/${#subFiles[@]}) in '$(basename $directory)' finished in ${differenceSub[counterSub]} seconds." | telegram-send --silent --stdin	
-			fi
-			let counterSub=counterSub+1
-		fi
-
-	done
-
-	seconds2=$(date '+%s')
-	let difference[counter]=seconds2-seconds1
-	sum=$(IFS=+; echo "$((${difference[*]}))")
-	let estimatedTime=(${#directories[@]}-counter)
-	let estimatedTime=estimatedTime*sum/counter
-
-	if [ "$estimatedTime" -ne "0" ]; then
-		finishDate=$(date '+%a, %d.%m., %H:%Mh' --date="@$((seconds2 + estimatedTime))")
-		printf "${HOSTNAME:11:7}: Simulation of folder $counter/${#directories[@]} in '${1:0:-1}' finished in $(date -d@${difference[counter]} -u +%Hh%Mm).\nEstimated date finished: $finishDate." | telegram-send --silent --stdin
-	else
-		printf "${HOSTNAME:11:7}: Simulation of folder $counter/${#directories[@]} in '${1:0:-1}' finished in $(date -d@${difference[counter]} -u +%Hh%Mm)." | telegram-send --silent --stdin	
-	fi
-	let counter=counter+1
-done
-
-cd $PWD
-printf "***************************************\n${HOSTNAME:11:7}: Simulation of folder '${1:0:-1}' completed. Total time: $(date -d@$sum -u +%Hh%Mm).\n***************************************" | telegram-send --stdin
-echo "All done!"
\ No newline at end of file