diff --git a/EPR_script.m b/EPR_script.m
deleted file mode 100644
index d207223..0000000
--- a/EPR_script.m
+++ /dev/null
@@ -1,217 +0,0 @@
-% high-spin S=1 simulation
-% Inputs requested in command line at certain points
-
-clear variables
-close all
-
- % window positions (currently optimised for dual WQHD with main on right)
- % also uncomment all figure(gcf) when working with single monitor
-% Give desired position of figure window as number of pixels [pos_x pos_y size_x size_y]:
-position = [-1250,50,1200,800];
-% Give desired position of figure(2) window (will have two stacked subplots)
-% as number of pixels [pos_x pos_y size_x size_y]: 
-position2 = [-2000,50,700,800];
-% Give desired position of figure(3) window (will contain EPR spectrum and 
-% simulation) as number of pixels [pos_x pos_y size_x size_y]: 
-position3 = [-1250,50,1200,600];
-
-% specify dir for printing figures
-figdir = './';
-% specifiy excel-file for saving parameters
-table_path = 'example_results.xlsx';
-
-%% loading Data
-path = input('Path to dataset: ','s');
-load(path)
-
-whos % what variables have been loaded
-params % what information is contained in the structure called 'params'
-
-% get name of dataset
-dataname = string(extractBefore(extractAfter(path,asManyOfPattern(wildcardPattern + "/")),'.'));
-% using Parallel Computing toolbox to speed up (check with "gpuDevice" if you can use this)
-gpuData = gpuArray(Data);
-
-%% Baseline Correcting
-% plot the raw data & check the number of points before the signal (pre-trigger)
-plot(gpuData)
-title('raw data')
-set(gcf,'Position',position)
-
-% substract the pre-trigger
-pre_trigger = input('Number of pre-trigger points: ');
-signal_baseline_time = bsxfun(@minus, gpuData, mean(gpuData(1:pre_trigger,:)));
-plot(signal_baseline_time) % plot the corrected data set
-title('time corrected data')
-set(gcf,'Position',position)
-% figure(gcf) % bring figure to foreground
-
-ready = input('Proceed?');
-
- % plot the transpose and check the number of points to the lower and higher fields of the signal
-plot(signal_baseline_time.')
-title('transposed time corrected data')
-set(gcf,'Position',position)
-% figure(gcf) % bring figure to foreground
-
-
- % BASELINE correction
-baseline_points = input('Number of baseline points (use smaller value from left and right): ');
-l1 = mean(signal_baseline_time(:,1:baseline_points),2); % calculate the mean on the left along the time axis
-l2 = mean(signal_baseline_time(:,end-baseline_points:end),2); %calculate the mean on the right along the time axis
-baseline_time = (l1 +l2)/2; %take the average
-
-signal_baseline_time_field = bsxfun(@minus, signal_baseline_time, baseline_time); % subtract the background in the time-domain 
-
- % plot the corrected data set
-plot(signal_baseline_time_field.')
-title('transposed fully corrected data')
-set(gcf,'Position',position)
-% figure(gcf) % bring figure to foreground
-
-clear ready
-ready = input('Proceed?');
-
- % plot the transpose to find the region of maximum signal. Use this below
-plot(signal_baseline_time_field)
-title('fully corrected data')
-set(gcf,'Position',position)
-% figure(gcf) % bring figure to foreground
-
- % contour plot: The index gives the number of contours
-% contourf(signal_baseline_field_time,6)
-
- % NORMALISING
-max_region = input('Region of the maximum signal as [x1:x2]: ');
-% take the mean over the maxium region. You can decide how wide it is
-signal_baseline_time_field_mean = (mean(signal_baseline_time_field(max_region,:)));
-% normalise the amplitude to 1
-signal_baseline_time_field_mean_norm = signal_baseline_time_field_mean/max(signal_baseline_time_field_mean);
-
-%% Creating figure with two subplots
-figure(2)
-set(gcf,'PaperUnits','centimeters')
-set(gcf,'Position',position2)
-set(gcf,'InvertHardcopy','off','Color',[1 1 1])
-set(0,'DefaultAxesFontSize', 12,'DefaultAxesLineWidth',2)
-
-cont_or_surf = input('Should lower subplot be contour(1) or surface(2) plot? (1/2): ');
-subplot(2,1,2)
-if cont_or_surf == 1
-    % contour plot: add the time and field axes
-    contourf(0.1*params.Field_Vector, TimeBase*1e6 ,signal_baseline_time_field,'LineColor','none')
-elseif cont_or_surf == 2
-    % surface plot: add the time and field axes
-    surf(0.1*params.Field_Vector, TimeBase*1e6 ,signal_baseline_time_field)
-    colormap default
-    shading interp
-end
-
-xlabel('Magnetic Field / mT')
-ylabel('Time / \mus')
-
-subplot(2,1,1)
- % plot the spectrum
-plot(0.1*params.Field_Vector,signal_baseline_time_field_mean_norm,'LineWidth',2)
-
-xlabel('Magnetic Field / mT')
-axis('tight')
-box off
-
-%% Simulation section
-
-Exp.mwFreq = params.mwFreq; % GHz
-Exp.nPoints = length(params.Field_Vector);
-Exp.CenterSweep = 0.1*[params.Field_Center params.Field_Sweep]; % mT (converted from Gauss)
-Exp.Harmonic = 0; % zeroth harmonic
-
-init_proceed = 'n';
-while init_proceed == 'n'
-    % populations of the triplet sub-levels
-    % these need to be varied manually to get the right shape
-    Exp.Temperature = input('Input population of triplett sublevels as [T_x T_y T_z]: ');
-    % initial simulation settings
-    Sys.S = 1; % Total Spin
-    Sys.g = input('g value: '); % needs to be optimised
-    Sys.D = input('D and E value as [D E]: '); % mT, The D and E values need to be optimised
-    Sys.lw = input('Isotropic line broadening at FWHM as [Gaussian Lorentzian]: '); % mT, linewidth needs to be optimised
-    
-    [bfield,spec] = pepper(Sys,Exp); % perform a simulation with the parameters above
-    spec_norm = spec/max(spec); % normalize the simulation
-    
-    figure(3)
-    set (gcf,'PaperUnits','centimeters')
-    set (gcf,'Position',position3) % set the position, size and shape of the plot
-    set (gcf,'InvertHardcopy','off','Color',[1 1 1])
-    set(0,'DefaultAxesFontSize', 16,'DefaultAxesLineWidth',1.5)
-    plot(0.1*params.Field_Vector,signal_baseline_time_field_mean_norm,'r', bfield,spec_norm,'b','LineWidth',1);
-    axis('tight')
-    legend('experimental','simulation')
-    legend boxoff
-    xlabel('Magnetic Field / mT')
-    ylabel('EPR signal / A. U.')
-    set(gca,'Box','Off', 'XMinorTick','On', 'YMinorTick','On', 'TickDir','Out', 'YColor','k')
-    pause(2);
-    init_proceed = input('Spectrum shape manually fitted? [y/n]: ','s');
-end
-
- % variation settings for simulation
-Vary.g = 0.01; 
-Vary.D = [10 10];
-Vary.lw = [1 0];
- % further setup
-FitOpt.Method = 'simplex fcn';
-FitOpt.Scaling = 'lsq';
-
- % When you have got a good fit by eye, use esfit to optimise
-simu_proceed = 'n';
-while simu_proceed == 'n'
-    % fitting routine
-    [BestSys,BestSpc] = esfit('pepper',signal_baseline_time_field_mean_norm,Sys,Vary,Exp,[],FitOpt);
-    % plot best fit
-    figure(3)
-    plot(0.1*params.Field_Vector,signal_baseline_time_field_mean_norm,'r',...
-         0.1*params.Field_Vector,BestSpc,'b','LineWidth',1);
-    axis('tight')
-    legend('experimental','simulation')
-    legend boxoff
-    xlabel('Magnetic Field / mT')
-    ylabel('EPR signal / A. U.')
-    set(gca,'Box','Off', 'XMinorTick','On', 'YMinorTick','On', 'TickDir','Out', 'YColor','k')
-    
-    simu_proceed = input('Did the simulation converge? [y/n]: ','s');
-    if simu_proceed == 'n'
-        simu_val = input('Do you want to repeat the simulation with new best values? [y/n]: ','s');
-        if simu_val == 'y'
-            Sys.g = BestSys.g;
-            Sys.D = BestSys.D;
-            Sys.lw = BestSys.lw;
-        end
-    end
-end
-
-%% printing figures
-printing = input('Do you want to print figure(3)? [y/n]: ','s');
-if printing == 'y'
-    figure(3)
-    set(gcf,'Units','Inches');
-    pos = get(gcf,'Position');
-    set(gcf,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)]);
-    print(gcf,strcat(figdir,dataname),'-dpdf','-r0');
-end
-
-%% saving parameters
- % concatenate data to existing table
-table_old = readtable(table_path);
-table_old.Properties.VariableNames = {'filename', 'date', 'pre-trigger', ...
-    'baseline_points', 'max_area_left', 'max_area_right', 'T_x', 'T_y', 'T_z', ...
-    'sim. g-value', 'sim_D', 'sim_E', 'sim_lw_gauss', 'sim_lw_lorentz'};
-% new data as table
-table_new = table(dataname, string(datestr(clock)), pre_trigger, baseline_points, ...
-    max_region(1), max_region(end), Exp.Temperature(1), Exp.Temperature(2), Exp.Temperature(3), ...
-    BestSys.g, BestSys.D(1), BestSys.D(2), BestSys.lw(1), BestSys.lw(2), ...
-    'VariableNames', {'filename', 'date', 'pre-trigger', 'baseline_points', ...
-    'max_area_left', 'max_area_right', 'T_x', 'T_y', 'T_z', 'sim. g-value', ...
-    'sim_D', 'sim_E', 'sim_lw_gauss', 'sim_lw_lorentz'});
-table_conc = [table_old;table_new];
-writetable(table_conc,table_path)
diff --git a/double_2D_3D_plot.m b/double_2D_3D_plot.m
deleted file mode 100644
index 2d6d8d1..0000000
--- a/double_2D_3D_plot.m
+++ /dev/null
@@ -1,7 +0,0 @@
-function [outputArg1,outputArg2] = double_2D_3D_plot(inputArg1,inputArg2)
-%DOUBLE_2D_3D_PLOT Summary of this function goes here
-%   Detailed explanation goes here
-outputArg1 = inputArg1;
-outputArg2 = inputArg2;
-end
-
diff --git a/example_results.xlsx b/example_results.xlsx
deleted file mode 100644
index 7f6af96..0000000
Binary files a/example_results.xlsx and /dev/null differ