initial files

This commit is contained in:
bartels-rausch committed 2023-10-31 11:16:30 +01:00
1 parent 62181940c5
commit b23894e06c
88 files changed
+13263

No files matched your search

Binary file not shown.
Binary file not shown.
Binary file not shown.
+82
View File
@@ -0,0 +1,82 @@
function LopapTimeT = fitLopapBkgrd(LopapTimeT,lopapFitParameter)
%FITLOPAPBKGRD fits and subtracts the background from Lopap data. Input is
%a timetable with the processed data for both channels.
%
% LopapTimeT = fitLopapBkgrd(LopapTimeT) fits a linear background through 2
% periods a 15 minutes at the beginning and at the end of the data set for
% both channels.
%
% LopapTimeT = fitLopapBkgrd(LopapTimeT,'BkgrdPeriod',[datenum('14-Sep-2020 14:50:00')
% datenum('14-Sep-2020 15:15:00'); datenum('14-Sep-2020 16:50:00')
% datenum('14-Sep-2020 17:15:00')]) fits a linear background through 2
% periods as specified by the time matrix for both channels: This matrix
% holds the start time and end time in col1 and col2 of each line for
% individual regions that are used to fit. The number of lines defines the
% number of regions to be used to fit. All data outside these regions are
% excluded from the fits.
%
% LopapTimeT = fitLopapBkgrd(LopapTimeT,'ChannelNumber',2) Processing is
% restricted data of channel 2. Possible input are 1, 2, [1 2]. Default is [1 2]
%
% The optional input can also be combined:
%
% localFitBkgrd(LopapTimeT,'ChannelNumber', 1, 'FitType','exp1') fits data
% in channel 1 with an exponential function. Standard FitType is linear
% ('poly1').
% Check the input: Is channel reasonable, is end-date after start-date?
arguments
LopapTimeT timetable
lopapFitParameter.ChannelNumber double {mustBeMember(lopapFitParameter.ChannelNumber,[1,2])} = [1 2]
lopapFitParameter.BkgrdPeriod double {mustBeLopapLaterDate(lopapFitParameter.BkgrdPeriod)} = ...
[datenum(LopapTimeT.Time(1)+minutes(10)) datenum(LopapTimeT.Time(1)+minutes(25)); ...
datenum(LopapTimeT.Time(end)-minutes(25)) datenum(LopapTimeT.Time(end)-minutes(10))]
lopapFitParameter.FitType char ...
{mustBeMember(lopapFitParameter.FitType,{'poly1','poly2','poly3','poly4','poly5','poly6','poly7','poly8','poly9','exp1','exp2'})} = "poly1"
lopapFitParameter.TimeRangeToPlot
end
% gets index (row number) of data for given time periods
for ii = 1:size(lopapFitParameter.BkgrdPeriod,1)
index(ii,1) = find(datenum(LopapTimeT.Time) >= lopapFitParameter.BkgrdPeriod(ii,1),1);
index(ii,2) = find(datenum(LopapTimeT.Time) >= lopapFitParameter.BkgrdPeriod(ii,2),1);
end
% gets index (row number) of data for given time periods
for ii = 1:size(lopapFitParameter.TimeRangeToPlot,1)
indexPlot(ii,1) = find(datenum(LopapTimeT.Time) >= lopapFitParameter.TimeRangeToPlot(ii,1),1);
indexPlot(ii,2) = find(datenum(LopapTimeT.Time) >= lopapFitParameter.TimeRangeToPlot(ii,2),1);
end
% define data
for ii = 1:size(lopapFitParameter.ChannelNumber,2)
[xData, yData] = prepareCurveData( [], LopapTimeT.(['Ch #' num2str(lopapFitParameter.ChannelNumber(ii))]));
clearvars excluded*
%opts = rmfield(opts,'Exclude');
% Set up fittype and options.
ft = fittype( lopapFitParameter.FitType );
% exclude all data outside the given time periods (based on the index):
for ij = 1:size(lopapFitParameter.BkgrdPeriod,1)
excludedRegion(:,ij) = excludedata( xData, yData, 'domain',...
[index(ij) index(ij+size(lopapFitParameter.BkgrdPeriod,1))]);
end
excludedPoints = all(excludedRegion,2);
opts = fitoptions( 'Method', 'LinearLeastSquares' );
opts.Exclude = excludedPoints;
% Fit model to data.
[fitresult, ~] = fit( xData, yData, ft, opts );
% Plot fit with data.
figure( 'Name', ['Fit background Ch #' num2str(lopapFitParameter.ChannelNumber(ii))]);
h = plot( fitresult, xData, yData, excludedPoints );
legend( h, 'y', 'Excluded y', 'untitled fit 1', 'Location', 'NorthEast', 'Interpreter', 'none' );
% Label axes
ylabel( 'y', 'Interpreter', 'none' );
grid on
xlim([indexPlot(1,1) indexPlot(1,2)])
% Subtract Background
LopapTimeT.(['Ch #' num2str(lopapFitParameter.ChannelNumber(ii)) ' Bkgrd']) = ...
LopapTimeT.(['Ch #' num2str(lopapFitParameter.ChannelNumber(ii))]) -fitresult(1:height(LopapTimeT));
end
end
+32
View File
@@ -0,0 +1,32 @@
function [LopapTimeT, fitresult] = fitLopapCalibration(LopapTimeT,calib)
%CREATEFIT(X,Y)
% Create a fit.
%
% Data for 'untitled fit 1' fit:
% X Input : x
% Y Output: y
% Output:
% fitresult : a fit object representing the fit.
% gof : structure with goodness-of fit info.
%
% See also FIT, CFIT, SFIT.
% Auto-generated by MATLAB on 21-Sep-2020 11:28:06
%% Fit: 'untitled fit 1'.
[xData, yData] = prepareCurveData([0 calib.mean(3)],[0 calib.mixingRatio(1)]);
% Set up fittype and options.
ft = fittype( 'poly1' );
% Fit model to data.
[fitresult.calibChannel1, ~ ] = fit( xData, yData, ft );
LopapTimeT.('Ch #1 Calib') = fitresult.calibChannel1(LopapTimeT.('Ch #1 Bkgrd'));
[xData, yData] = prepareCurveData([0 calib.mean(4)],[0 calib.mixingRatio(2)]);
% Fit model to data.
[fitresult.calibChannel2, ~] = fit( xData, yData, ft );
% calibrate
LopapTimeT.('Ch #2 Calib') = fitresult.calibChannel2(LopapTimeT.('Ch #2 Bkgrd'));
+37
View File
@@ -0,0 +1,37 @@
function [lopapReferenceWavelength, LopapDataCh1,LopapDataCh2,LopapSpecCh1,LopapSpecCh2]...
= importLopap(lopapDateNum,fileLocation)
%IMPORTLOPAP wraps import functions to read data from all Lopap files for one experiments
% Check input to function, is the path correct?
arguments
lopapDateNum double
fileLocation struct {mustBeLopapFolder(fileLocation, lopapDateNum)}
end
% generate full path name
lopapPath = fullfile(fileLocation.server,fileLocation.dataPath,datestr(lopapDateNum,'yyyymmdd_HHMM'));
% Import reference wavelength
%
% This reads the reference wavelength for both channels as recorded in the
% *.dat file
lopapReferenceWavelength = importLopapReferenceWavelength(fullfile(lopapPath,[datestr(lopapDateNum,'yyyymmdd_HHMM') '.dat']), [2, 3]);
% Import data
%
% Here, we read the data from the Lopap .dat file into timetables, one for
% each channel. This data holds the intensity at selected wavelength
% divided by the reference wavelength. The information is thus the same as
% shown in the "measurement plot" window of the Lopap software.
[LopapDataCh1,LopapDataCh2] = importLopapDatFile(fullfile(lopapPath,[datestr(lopapDateNum,'yyyymmdd_HHMM') '.dat']), [6, Inf]);
% Import spectra
%
% Here, we read the data from the Lopap .spec file into timetables, one for
% each channel. The data hold the all spectra acquired during the
% experiment. The information is thus the same as in the "spectrum plot"
% window of the Lopap software.
[LopapSpecCh1,LopapSpecCh2] = importLopapSpecFile(fullfile(lopapPath,[datestr(lopapDateNum,'yyyymmdd_HHMM') '.spec']), [1, Inf]);
end
+98
View File
@@ -0,0 +1,98 @@
function [LopapDataCh1, LopapDataCh2] = importLopapDatFile(filename, dataLines)
%IMPORTLOPAPDATFILE Import data from a Lopap .dat file
% [LopapDataCh1, LopapDataCh2] = importLopapDatFile(FILENAME) reads data
% from Lopap .dat text file FILENAME for the default selection. Returns
% the data as timetables for each Lopap channel with time in rows and
% wavelength in columns.
%
% [LOPAPDATACH1, LOPAPDATACH2] = importLopapDatFile(FILE, DATALINES) reads
% data for the specified row interval(s) of the Lopap .dat text file
% FILENAME. Specify DATALINES as a positive scalar integer or a N-by-2
% array of positive scalar integers for dis-contiguous row intervals.
%
% Example:
% [LopapDataCh1, LopapDataCh2] = importLopapDatFile("20200619_1320.dat", [5, Inf]);
%
% See also READTABLE.
%
% Auto-generated by MATLAB on 14-Sep-2020 13:27:21
%% Lopap Measurement Data Files (.dat)
% The measurement data files contain the processed measurement values,
% which were written while the running measurement. This file is constantly
% updated while the measurement is running, thus even if the application is
% terminated accidentally, this contains the measurement values until the
% time of termination.
%
% *File Format* The file format is very straightforward.
% A Measurement Data File consists of a header in which the file format and
% the reference wavelengths for every spectrometer channel is recorded.
% After this header, the data section starts. The data is arranged in a
% chronological semicolon separated value list, with one line per sample
% time. The first two columns give the date and the time of the data
% sample. The following columns give the measurement values for every
% measurement wavelength of the spectrometer channels. The spectrometer
% channels are separated with the string ChX where "X" stands for the
% spectrometer channel index. The very first line of the data section gives
% the meaning of the columns in a sort form.
%% Read header to string
opts = delimitedTextImportOptions("NumVariables", 26);
% Specify range and delimiter
opts.DataLines = [5, 5];
opts.Delimiter = ";";
% Specify column names and types
opts.VariableNames = ["Date", "Time", "Ch1", "VarName4", "VarName5", "VarName6", "VarName7", "VarName8", "VarName9", "VarName10", "VarName11", "VarName12", "VarName13", "VarName14", "Ch2", "VarName16", "VarName17", "VarName18", "VarName19", "VarName20", "VarName21", "VarName22", "VarName23", "VarName24", "VarName25", "VarName26"];
opts.VariableTypes = ["string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string"];
% Specify file level properties
opts.ExtraColumnsRule = "ignore";
opts.EmptyLineRule = "read";
% Specify variable properties
opts = setvaropts(opts, ["Date", "Time", "Ch1", "VarName4", "VarName5", "VarName6", "VarName7", "VarName8", "VarName9", "VarName10", "VarName11", "VarName12", "VarName13", "VarName14", "Ch2", "VarName16", "VarName17", "VarName18", "VarName19", "VarName20", "VarName21", "VarName22", "VarName23", "VarName24", "VarName25", "VarName26"], "WhitespaceRule", "preserve");
opts = setvaropts(opts, ["Date", "Time", "Ch1", "VarName4", "VarName5", "VarName6", "VarName7", "VarName8", "VarName9", "VarName10", "VarName11", "VarName12", "VarName13", "VarName14", "Ch2", "VarName16", "VarName17", "VarName18", "VarName19", "VarName20", "VarName21", "VarName22", "VarName23", "VarName24", "VarName25", "VarName26"], "EmptyFieldRule", "auto");
% Import the data
variableNames = readmatrix(filename, opts);
%% Input handling
% If dataLines is not specified, define defaults
if nargin < 2
dataLines = [6, Inf];
end
%% Setup the Import Options and import the data
clearvars opts
opts = delimitedTextImportOptions("NumVariables", 26);
% Specify range and delimiter
opts.DataLines = dataLines;
opts.Delimiter = ";";
% Specify column names and types
opts.VariableNames = variableNames(:);
opts.VariableTypes = ["datetime", "datetime", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
% Specify file level properties
opts.ExtraColumnsRule = "ignore";
opts.EmptyLineRule = "read";
opts.PreserveVariableNames = 1;
% Specify variable properties
opts = setvaropts(opts, "Date", "InputFormat", "dd.MM.yyyy");
opts = setvaropts(opts, "Time", "InputFormat", "HH:mm:ss");
opts = setvaropts(opts, ["Ch1", "Ch2"], "TrimNonNumeric", true);
opts = setvaropts(opts, ["Ch1", "Ch2"], "ThousandsSeparator", ",");
% Import the data
T = readtable(filename, opts);
T.Date.Format = 'dd.MM.uuuu HH:mm:ss';
T.Time.Format = 'dd.MM.uuuu HH:mm:ss';
columnCh2 = find(strcmp(opts.VariableNames,'Ch2') >=1,1);
LopapDataCh1 = table2timetable(T(:,4:columnCh2-1),'RowTimes',T.Date + timeofday(T.Time));
LopapDataCh2 = table2timetable(T(:,columnCh2+1:end),'RowTimes',T.Date + timeofday(T.Time));
end
+51
View File
@@ -0,0 +1,51 @@
function VarName2 = importLopapPrjFile(filename, dataLines)
%IMPORTLOPAPPRJFILE Import data from a Lopap .prj text file
% VARNAME2 = importLopapPrjFile(FILENAME) reads data from Lopap .prj text
% file FILENAME for the default selection. Returns the data as column
% vectors.
%
% VARNAME2 = importLopapPrjFile(FILE, DATALINES) reads data for the specified
% row interval(s) of Lopap .prj text file FILENAME. Specify DATALINES as a
% positive scalar integer or a N-by-2 array of positive scalar integers
% for dis-contiguous row intervals.
%
% Example:
% wavelengthCh1 = importLopapPrjFile("20200619_1320.prj", [10, 10]);
%
% See also READTABLE.
%
% Auto-generated by MATLAB on 15-Sep-2020 08:38:14
%% Input handling
% If dataLines is not specified, define defaults
if nargin < 2
dataLines = [10, Inf];
end
%% Setup the Import Options and import the data
opts = delimitedTextImportOptions("NumVariables", 2);
% Specify range and delimiter
opts.DataLines = dataLines;
opts.Delimiter = "=";
% Specify column names and types
opts.VariableNames = ["Var1", "VarName2"];
opts.SelectedVariableNames = "VarName2";
opts.VariableTypes = ["string", "string"];
% Specify file level properties
opts.ExtraColumnsRule = "ignore";
opts.EmptyLineRule = "read";
% Specify variable properties
opts = setvaropts(opts, ["Var1", "VarName2"], "WhitespaceRule", "preserve");
opts = setvaropts(opts, ["Var1", "VarName2"], "EmptyFieldRule", "auto");
% Import the data
tbl = readtable(filename, opts);
%% Convert to output type
VarName2 = tbl.VarName2;
end
@@ -0,0 +1,50 @@
function referenceWavelength = importLopapReferenceWavelength(filename, dataLines)
%IMPORTLOPAPREFERENCEWAVELENGTH IMPORT DATA FROM A TEXT FILE
% REFERENCEWAVELENGTH = importReferenceWavelength(FILENAME) reads data from Lopap
% .dat text file FILENAME for the default selection. Returns the data as
% a string array with text string for channel 1 in 1st row and for channel
% 2 in second row.
%
% REFERENCEWAVELENGTH = importReferenceWavelength(FILE, DATALINES) reads data for the specified
% row interval(s) of text file FILENAME. Specify DATALINES as a
% positive scalar integer or a N-by-2 array of positive scalar integers
% for dis-contiguous row intervals.
%
% Example:
% referenceWavelength = importReferenceWavelength("20200619_1320.dat", [2, 3]);
%
% See also READTABLE.
%
% Auto-generated by MATLAB on 14-Sep-2020 13:25:26
%% Input handling
% If dataLines is not specified, define defaults
if nargin < 2
dataLines = [2, Inf];
end
%% Setup the Import Options and import the data
opts = delimitedTextImportOptions("NumVariables", 26);
% Specify range and delimiter
opts.DataLines = dataLines;
opts.Delimiter = ";";
% Specify column names and types
opts.VariableNames = ["LOPAPlicmeasurementdatav01", "Var2", "Var3", "Var4", "Var5", "Var6", "Var7", "Var8", "Var9", "Var10", "Var11", "Var12", "Var13", "Var14", "Var15", "Var16", "Var17", "Var18", "Var19", "Var20", "Var21", "Var22", "Var23", "Var24", "Var25", "Var26"];
opts.SelectedVariableNames = "LOPAPlicmeasurementdatav01";
opts.VariableTypes = ["string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string", "string"];
% Specify file level properties
opts.ExtraColumnsRule = "ignore";
opts.EmptyLineRule = "read";
% Specify variable properties
opts = setvaropts(opts, ["LOPAPlicmeasurementdatav01", "Var2", "Var3", "Var4", "Var5", "Var6", "Var7", "Var8", "Var9", "Var10", "Var11", "Var12", "Var13", "Var14", "Var15", "Var16", "Var17", "Var18", "Var19", "Var20", "Var21", "Var22", "Var23", "Var24", "Var25", "Var26"], "WhitespaceRule", "preserve");
opts = setvaropts(opts, ["LOPAPlicmeasurementdatav01", "Var2", "Var3", "Var4", "Var5", "Var6", "Var7", "Var8", "Var9", "Var10", "Var11", "Var12", "Var13", "Var14", "Var15", "Var16", "Var17", "Var18", "Var19", "Var20", "Var21", "Var22", "Var23", "Var24", "Var25", "Var26"], "EmptyFieldRule", "auto");
% Import the data
referenceWavelength = readmatrix(filename, opts);
end
+86
View File
@@ -0,0 +1,86 @@
function [LopapSpecCh1, LopapSpecCh2] = importLopapSpecFile(filename, dataLines)
%IMPORTLOPAPSPECFILE imports data from a Lopap .spec file
% [LopapSpecCh1, LopapSpecCh2] = importLopapSpecFile(FILENAME) reads data
% from Lopap .spec text file FILENAME for the default selection. Returns
% the data as a timetable for each channel with time in rows and
% wavelength in nm in columns.
%
% [LOPAPSPECCH1, LOPAPSPECCH2] = importLopapSpecFile(FILE, DATALINES) reads data for the specified
% row interval(s) of text file FILENAME. Specify DATALINES as a
% positive scalar integer or a N-by-2 array of positive scalar integers
% for dis-contiguous row intervals.
%
% Example:
% [LopapSpecCh1, LopapSpecCh2] = importLopapSpecFile("20200619_1320.spec", [1, Inf]);
%
% See also READTABLE.
%
% Auto-generated by MATLAB on 14-Sep-2020 16:40:15
%% Input handling
% If dataLines is not specified, define defaults
if nargin < 2
dataLines = [1, Inf];
end
%% Setup the Import Options and import the data
opts = delimitedTextImportOptions("NumVariables", 2051);
% Specify range and delimiter
opts.DataLines = dataLines;
opts.Delimiter = ";";
% Specify column names and types
opts.VariableNames = ["VarName1", "VarName2", "VarName3", "VarName4", "VarName5", "VarName6", "VarName7", "VarName8", "VarName9", "VarName10", "VarName11", "VarName12", "VarName13", "VarName14", "VarName15", "VarName16", "VarName17", "VarName18", "VarName19", "VarName20", "VarName21", "VarName22", "VarName23", "VarName24", "VarName25", "VarName26", "VarName27", "VarName28", "VarName29", "VarName30", "VarName31", "VarName32", "VarName33", "VarName34", "VarName35", "VarName36", "VarName37", "VarName38", "VarName39", "VarName40", "VarName41", "VarName42", "VarName43", "VarName44", "VarName45", "VarName46", "VarName47", "VarName48", "VarName49", "VarName50", "VarName51", "VarName52", "VarName53", "VarName54", "VarName55", "VarName56", "VarName57", "VarName58", "VarName59", "VarName60", "VarName61", "VarName62", "VarName63", "VarName64", "VarName65", "VarName66", "VarName67", "VarName68", "VarName69", "VarName70", "VarName71", "VarName72", "VarName73", "VarName74", "VarName75", "VarName76", "VarName77", "VarName78", "VarName79", "VarName80", "VarName81", "VarName82", "VarName83", "VarName84", "VarName85", "VarName86", "VarName87", "VarName88", "VarName89", "VarName90", "VarName91", "VarName92", "VarName93", "VarName94", "VarName95", "VarName96", "VarName97", "VarName98", "VarName99", "VarName100", "VarName101", "VarName102", "VarName103", "VarName104", "VarName105", "VarName106", "VarName107", "VarName108", "VarName109", "VarName110", "VarName111", "VarName112", "VarName113", "VarName114", "VarName115", "VarName116", "VarName117", "VarName118", "VarName119", "VarName120", "VarName121", "VarName122", "VarName123", "VarName124", "VarName125", "VarName126", "VarName127", "VarName128", "VarName129", "VarName130", "VarName131", "VarName132", "VarName133", "VarName134", "VarName135", "VarName136", "VarName137", "VarName138", "VarName139", "VarName140", "VarName141", "VarName142", "VarName143", "VarName144", "VarName145", "VarName146", "VarName147", "VarName148", "VarName149", "VarName150", "VarName151", "VarName152", "VarName153", "VarName154", "VarName155", "VarName156", "VarName157", "VarName158", "VarName159", "VarName160", "VarName161", "VarName162", "VarName163", "VarName164", "VarName165", "VarName166", "VarName167", "VarName168", "VarName169", "VarName170", "VarName171", "VarName172", "VarName173", "VarName174", "VarName175", "VarName176", "VarName177", "VarName178", "VarName179", "VarName180", "VarName181", "VarName182", "VarName183", "VarName184", "VarName185", "VarName186", "VarName187", "VarName188", "VarName189", "VarName190", "VarName191", "VarName192", "VarName193", "VarName194", "VarName195", "VarName196", "VarName197", "VarName198", "VarName199", "VarName200", "VarName201", "VarName202", "VarName203", "VarName204", "VarName205", "VarName206", "VarName207", "VarName208", "VarName209", "VarName210", "VarName211", "VarName212", "VarName213", "VarName214", "VarName215", "VarName216", "VarName217", "VarName218", "VarName219", "VarName220", "VarName221", "VarName222", "VarName223", "VarName224", "VarName225", "VarName226", "VarName227", "VarName228", "VarName229", "VarName230", "VarName231", "VarName232", "VarName233", "VarName234", "VarName235", "VarName236", "VarName237", "VarName238", "VarName239", "VarName240", "VarName241", "VarName242", "VarName243", "VarName244", "VarName245", "VarName246", "VarName247", "VarName248", "VarName249", "VarName250", "VarName251", "VarName252", "VarName253", "VarName254", "VarName255", "VarName256", "VarName257", "VarName258", "VarName259", "VarName260", "VarName261", "VarName262", "VarName263", "VarName264", "VarName265", "VarName266", "VarName267", "VarName268", "VarName269", "VarName270", "VarName271", "VarName272", "VarName273", "VarName274", "VarName275", "VarName276", "VarName277", "VarName278", "VarName279", "VarName280", "VarName281", "VarName282", "VarName283", "VarName284", "VarName285", "VarName286", "VarName287", "VarName288", "VarName289", "VarName290", "VarName291", "VarName292", "VarName293", "VarName294", "VarName295", "VarName296", "VarName297", "VarName298", "VarName299", "VarName300", "VarName301", "VarName302", "VarName303", "VarName304", "VarName305", "VarName306", "VarName307", "VarName308", "VarName309", "VarName310", "VarName311", "VarName312", "VarName313", "VarName314", "VarName315", "VarName316", "VarName317", "VarName318", "VarName319", "VarName320", "VarName321", "VarName322", "VarName323", "VarName324", "VarName325", "VarName326", "VarName327", "VarName328", "VarName329", "VarName330", "VarName331", "VarName332", "VarName333", "VarName334", "VarName335", "VarName336", "VarName337", "VarName338", "VarName339", "VarName340", "VarName341", "VarName342", "VarName343", "VarName344", "VarName345", "VarName346", "VarName347", "VarName348", "VarName349", "VarName350", "VarName351", "VarName352", "VarName353", "VarName354", "VarName355", "VarName356", "VarName357", "VarName358", "VarName359", "VarName360", "VarName361", "VarName362", "VarName363", "VaLine truncated
opts.VariableTypes = ["double", "datetime", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "doubLine truncated
% Specify file level properties
opts.ExtraColumnsRule = "ignore";
opts.EmptyLineRule = "read";
% Specify variable properties
opts = setvaropts(opts, "VarName2", "InputFormat", "HH:mm:ss dd.MM.yyyy");
opts = setvaropts(opts, ["VarName1", "VarName2051"], "TrimNonNumeric", true);
opts = setvaropts(opts, ["VarName1", "VarName2051"], "ThousandsSeparator", ",");
% Import the data
T = readtable(filename, opts);
%% Import wavelength of spectrometer channels
[pfad,name,~] = fileparts(filename);
wavelengthCh1 = importLopapPrjFile([pfad filesep name '.prj'], [10, 10]);
wavelengthCh1 = split(wavelengthCh1,",")';
idx = cellfun(@isempty,wavelengthCh1);
wavelengthCh1(idx) = "2000";
wavelengthCh2 = importLopapPrjFile([pfad filesep name '.prj'], [13, 13]);
wavelengthCh2 = split(wavelengthCh2,",")';
idx = cellfun(@isempty,wavelengthCh2);
wavelengthCh2(idx) = "2000";
rows = T.VarName1 == 0;
T.VarName2.Format = 'HH:mm:ss dd.MM.uuuu';
LopapSpecCh1 = table2timetable(T(rows,3:end),'RowTimes',T.VarName2(rows));% LopapSpecCh1 = renamevars(LopapSpecCh1,1:width(LopapSpecCh1),wavelengthCh1)
LopapSpecCh1.Properties.VariableNames = wavelengthCh1;
rows = T.VarName1 == 1;
T.VarName2.Format = 'HH:mm:ss dd.MM.uuuu';
LopapSpecCh2 = table2timetable(T(rows,3:end),'RowTimes',T.VarName2(rows));
LopapSpecCh2.Properties.VariableNames = wavelengthCh2;
% Subtract spectrometer offset
% Wenn Ihr direkt die
% Spektren benutzen wollt (log(Iref/I), dann müßt Ihr erst den Dunkelstrom
% des Spektrometers abziehen. Den kannst Du gut sehen, wenn die Lampe aus
% ist. Der ist je nach Spektrometer aber meist unabhängig von der
% Wellenlänge. Daher nehmen wir bei den alten Spektrometern (S2000 bei
% Euch?) die Pixel 11-61 und ziehen diesen Mittelwert von allen Signalen
% ab. Erst dann wird log(Iref/I0) gebildet:
LopapSpecCh1.Variables = LopapSpecCh1.Variables - mean(LopapSpecCh1{:,10:59},2);
LopapSpecCh2.Variables = LopapSpecCh2.Variables - mean(LopapSpecCh2{:,10:59},2);
end
+64
View File
@@ -0,0 +1,64 @@
function TimeRange = modifyLopapTimeRanges(TimeRange, LopapParameters, LopapDataCh1, LopapDataCh2)
%
if LopapParameters.Interactive == 1
close all
clc
list = {'Data Range','Background Fit',...
'Compare Spectrum #1','Compare Spectrum #2'};
indx = listdlg('ListString',list,'PromptString',{'Choose parameter to modify,','then select region in figure'}...
,'SelectionMode','single');
if indx == 2
prompt = {'Enter number of regions to fit the background:'};
dlgtitle = 'Background Fit';
dims = [1 40];
definput = {'2'};
answer = inputdlg(prompt,dlgtitle,dims,definput);
end
if indx >= 1
clc
close all
figure( 'Name', 'LopapChannel 550','NumberTitle','off','Position',[10 100 550 800]);
plot(datenum(LopapDataCh1.Time),LopapDataCh1.("485.5"));
hold on
datetick('x','dd.mm hh:MM','keeplimits','keepticks');
ax = gca;
yyaxis right
plot(datenum(LopapDataCh2.Time),LopapDataCh2.("550.5"));
datetick('x','dd.mm hh:MM','keeplimits','keepticks');
grid ON
% Add title and axis labels
xlabel('Lab Time')
legend({'MeasData Ch#1','MeasData Ch#2'})
% select regions in plot
if indx == 1
Roi.calibrationRange = drawline(ax,'Color','r');
TimeRange.dataRange = [Roi.calibrationRange.Position(1,1) Roi.calibrationRange.Position(2,1)];
pause(3)
close all
elseif indx == 2
for ii = 1:str2double(answer)
Roi.backgroundFit = drawline(ax,'Color','r');
TimeRange.backgroundFit(ii,:) = [Roi.backgroundFit.Position(1,1) Roi.backgroundFit.Position(2,1)];
end
pause(3)
close all
elseif indx == 3
Roi.specWithReagents = drawpoint(ax,'Color','r');
TimeRange.specCompareOne = Roi.specWithReagents.Position(1,1);
pause(3)
close all
elseif indx == 4
Roi.specWithNitrite = drawpoint(ax,'Color','r');
TimeRange.specCompareTwo = Roi.specWithNitrite.Position(1,1);
pause(3)
close all
end
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'],'TimeRange');
end
else
disp('Variable LopapParameters.Interactive not set to 1: Modify data ranges manually has been skiped');
end
end
+117
View File
@@ -0,0 +1,117 @@
function TimeRange = modifyLopapTimeRangesCalib(TimeRange, LopapParameters, LopapDataCh1, LopapDataCh2)
%
if LopapParameters.Interactive == 1
close all
clc
list = {'Calibration Range','Background Fit','Mean ZeroAir','Mean Nitrit',...
'Spectrum One','Spectrum Two'};
indx = listdlg('ListString',list,'PromptString',{'Choose parameter to modify,','then select region in figure'}...
,'SelectionMode','single');
if indx == 2
prompt = {'Enter number of regions to fit the background:'};
dlgtitle = 'Background Fit';
dims = [1 40];
definput = {'2'};
answer = inputdlg(prompt,dlgtitle,dims,definput);
end
if indx >= 1
clc
close all
figure( 'Name', 'LopapChannel 550','NumberTitle','off','Position',[10 100 550 800]);
plot(datenum(LopapDataCh1.Time),LopapDataCh1.("550.7"));
hold on
datetick('x','dd.mm hh:MM','keeplimits','keepticks');
ax = gca;
yyaxis right
plot(datenum(LopapDataCh2.Time),LopapDataCh2.("550.5"));
datetick('x','dd.mm hh:MM','keeplimits','keepticks');
grid ON
% Add title and axis labels
xlabel('Lab Time')
legend({'MeasData Ch#1','MeasData Ch#2'})
% select regions in plot
if indx == 1
Roi.calibrationRange = drawline(ax,'Color','r');
TimeRange.calibrationRange = [Roi.calibrationRange.Position(1,1) Roi.calibrationRange.Position(2,1)];
pause(3)
close all
elseif indx == 2
for ii = 1:str2double(answer)
Roi.backgroundFit = drawline(ax,'Color','r');
TimeRange.backgroundFit(ii,:) = [Roi.backgroundFit.Position(1,1) Roi.backgroundFit.Position(2,1)];
end
pause(3)
close all
elseif indx == 3
Roi.meanZeroAir = drawline(ax,'Color','r');
TimeRange.MeanZeroAir = [Roi.meanZeroAir.Position(1,1) Roi.meanZeroAir.Position(2,1)];
pause(3)
close all
elseif indx == 4
Roi.meanNitriteStandard = drawline(ax,'Color','r');
TimeRange.MeanNitriteStandard = [Roi.meanNitriteStandard.Position(1,1) Roi.meanNitriteStandard.Position(2,1)];
pause(3)
close all
elseif indx == 5
Roi.specCompareOne = drawpoint(ax,'Color','r');
TimeRange.specCompareOne = Roi.specCompareOne.Position(1,1);
pause(3)
close all
elseif indx == 6
Roi.specCompareTwo = drawpoint(ax,'Color','r');
TimeRange.specCompareTwo = Roi.specCompareTwo.Position(1,1);
pause(3)
close all
end
disp("***************************************************************")
disp("* *")
disp("* If modified, COPY THIS CODE TO YOUR SCRIPT *")
disp("* *")
disp("***************************************************************")
disp(" % State the start and end time of the calibration. Usually")
disp(" % this is included zero air with R1&2, zero air with nitrite standard,")
disp(" % zero air with R1&2")
fprintf(" TimeRange.calibrationRange = [datenum('%s') datenum('%s')]; \n",datestr(TimeRange.calibrationRange(1)), datestr(TimeRange.calibrationRange(2)));
disp(" % <Adapt>: State the time for the first spectrum that you like to")
disp(" % compare (to confirm nitrite)")
fprintf(" TimeRange.specCompareOne = datenum('%s'); \n", datestr(TimeRange.specCompareOne(1)));
disp(" % <Adapt>: State the time for one spectrum of the second spectrum that")
disp(" % you like to compare (to confirm nitrite)")
fprintf(" TimeRange.specCompareTwo = datenum('%s'); \n", datestr(TimeRange.specCompareTwo(1)));
disp(" % <Adapt>: Choose a start and end time for one as many intervals that can")
disp(" % be used to fit the background. Background is the signal when R1&R2 was")
disp(" % attached to the Lopap with zero air entering the sampling coil.")
if size(TimeRange.backgroundFit,1) ==1
fprintf(" TimeRange.backgroundFit = [datenum('%s') datenum('%s')]; \n",datestr(TimeRange.backgroundFit(1,1)), datestr(TimeRange.backgroundFit(1,2)))
end
if size(TimeRange.backgroundFit,1) ==2
fprintf(" TimeRange.backgroundFit = [datenum('%s') datenum('%s'); datenum('%s') datenum('%s')]; \n",datestr(TimeRange.backgroundFit(1,1)), datestr(TimeRange.backgroundFit(1,2)),datestr(TimeRange.backgroundFit(2,1)), datestr(TimeRange.backgroundFit(2,2)))
end
if size(TimeRange.backgroundFit,1) ==3
fprintf(" TimeRange.backgroundFit = [datenum('%s') datenum('%s'); datenum('%s') datenum('%s'); datenum('%s') datenum('%s')]; \n",datestr(TimeRange.backgroundFit(1,1)), datestr(TimeRange.backgroundFit(1,2)),datestr(TimeRange.backgroundFit(2,1)), datestr(TimeRange.backgroundFit(2,2)),datestr(TimeRange.backgroundFit(3,1)), datestr(TimeRange.backgroundFit(3,2)))
end
disp(" % <Adapt>: State the start and end time of data that can be used to")
disp(" % calculate the mean of the zero air signal without nitrite")
fprintf(" TimeRange.MeanZeroAir = [datenum('%s') datenum('%s')]; \n",datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)));
disp(" % <Adapt>: State the start and end time of data that can be used to")
disp(" % calculate the mean of the zero air signal with nitrite")
fprintf(" TimeRange.MeanNitriteStandard = [datenum('%s') datenum('%s')]; \n",datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)));
disp("***************************************************************")
if LopapParameters.Interactive == 1
distFig
fprintf('press any key after you have copied above text and when you want to proceed\n')
pause
%%
%
close all
end
end
else
disp('Variable LopapParameters.Interactive not set to 1: Modify data ranges manually has been skiped');
end
end
+13
View File
@@ -0,0 +1,13 @@
function mustBeLopapFolder(inputPath,inputDate)
% MUSTBEFOLDER checks two things: First that the drive hosting the data
% folder is accessible, i.e. that the windows server are mounted and that
% there is no typo in the path and secondly that the folder name generated
% from the experimental date exists.
if ~isfolder(fullfile(inputPath.server,inputPath.dataPath))
error(['The folder ' fullfile(inputPath.server,inputPath.dataPath) ' appears not to be mounted. Check your network connection.'])
end
if ~isfolder(fullfile(inputPath.server,inputPath.dataPath,datestr(inputDate,'yyyymmdd_HHMM')))
error(['The folder ' datestr(inputDate,'yyyymmdd_HHMM') ' appears not to exist. The windows shared drive seems to be mounted. Check the date.'])
end
end
+7
View File
@@ -0,0 +1,7 @@
function mustBeLopapLaterDate(input)
% MUSTBELATERDATE checks that the start date comes before end date in
% functions that require date ranges as input.
if input(1,1) >= input(1,2)
error(['First date (' datestr(input(1,1)) ') appears to be later than second one (',datestr(input(1,2)),').'])
end
end
+14
View File
@@ -0,0 +1,14 @@
function mustBeLopapReasonableDate(dateInput)
% MUSTBEREASONABLEATE checks that a date is reaonalbe, that is after 1990
% and before 2025
for ij = 1:length(dateInput)
if datenum(dateInput(1,ij)) <= datenum('01-Jan-1999 01:00:00')
error(['Date (' datestr(dateInput(1,ij)) ') appears to be earlier than 1999. Please check.'])
end
if datenum(dateInput(1,ij)) >= datenum('01-Jan-2025 01:00:00')
error(['Date (' datestr(dateInput(1,ij)) ') appears to be later than 2025. Please check.'])
end
end
end
File diff suppressed because it is too large. Load diff
Binary file not shown.

After

Width:  |  Height:  |  Size: 2.1 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 20 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 55 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 56 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 22 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 21 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 163 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 21 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 26 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 20 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 16 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 14 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 18 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 16 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 730 B

Binary file not shown.

After

Width:  |  Height:  |  Size: 3.8 KiB

BIN
View File
Binary file not shown.
Binary file not shown.
Binary file not shown.

After

Width:  |  Height:  |  Size: 143 KiB

Binary file not shown.
File diff suppressed because it is too large. Load diff
Binary file not shown.

After

Width:  |  Height:  |  Size: 2.3 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 22 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 55 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 54 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 23 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 22 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 145 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 24 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 30 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 22 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 16 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 14 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 19 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 17 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 443 B

Binary file not shown.

After

Width:  |  Height:  |  Size: 2.3 KiB

@@ -0,0 +1,584 @@
%% Lopap Template
% Uses data from 16.9.2020 13h26
%
%% Info: How to use this
%
% # Make a copy of this file for each Lopap measurement. Name it
% "Lopap"_date_time.m and store it in the folder with the Lopap data. %
% Remove all the general text and comments in your file and add notes and
% observations about the experiment.
%
% # You can run this script, publish it, or evaluate it section wise. For
% new data run it section wise to carefully update the settings and to
% observe the results.
% # Once the data is loaded, you can also run individual sections only.
% # Use the LabBook :: Chrono Notes section as lab book.
% # First, you need to give the date and time of the Lopap measurement to
% load and process the data.
% # All sections that can be modified are marked with the word "<Adapt>:".
% # You also need to give time regions such as for the background fit, as
% well as to calculate means. There are several options how to do this
%
% * While running the experiments and observing the signal in the Lopap
% software, you can enter the date and times directly in the LabBook :
% Select data ranges section.
% * You can manually select the regions in a plot. A window will pop up to
% select the region to modify. To modify several regions, re-run this
% section of the script. This overwrites the current settings.
% * The times of the regions are stored in the "TimeRanges.mat" file. If it
% exists, this is read; the LabBook : Select data ranges section is then
% ignored.
%% Info: Output & Variables
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
%
% * *TimeRanges.mat*: Time ranges after they have been manually modified. If
% this file exists, the time range definitions in this script are ignored.
% * *LopapMatlabDataDATE_TIME.mat*: Raw data and experimental settings.
% Background subtraction and calibration are not saved in this file. When
% this file exists, data will not be read from the windows server.
% * *LopapDATE_TIME.mat*: Processed data - absorption with time (column 1
% and 2: 'Ch #1' and 'Ch #2') after background (columns 3 and 4: 'Ch #1
% Bkgrd' and 'Ch #2 Bkgrd') subtraction and calibration (columns 5 and 6:
% 'Ch #1 Calib' and 'Ch #2 Calib') for one wavelength - that can be used
% for further comparison with other experiments and discussion.
%
%% Matlab: Initialise
%
% This starts Matlab and defines some common settings such as the path to
% the Lopap data. Any modifications can be done at the end of this file in
% where the local function _localInitializeMatlab_ is stored.
% clear variables in workspace, clear text in command window
clearvars
clc
% set the path to the data and define some custom colours for plotting
[fileLocation, myColorMap] = localInitializeMatlab;
% <Adapt>: Date when Lopap recording was started. Defines names of Lopap
% files and folders. Seconds are irrelevant, but format is important.
LopapParameters.MeasurementDate = datenum('16-Sep-2020 13:26:00');
%%
% Print info to command window
disp('The data is read from:');
disp(fullfile(fileLocation.server,fileLocation.dataPath,datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM')));
disp(' ')
disp('Date of experiment:')
disp(datestr(LopapParameters.MeasurementDate));
disp(' ')
disp('Date of analyis:')
disp(datestr(today));
%% Process: Import data
%
% Reads and imports data from Lopap *.dat, *.spec, and *.prj files into:
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
if isfile(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'])
load(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'])
else
[LopapParameters.instrumentReferenceWavelength, LopapDataCh1,LopapDataCh2,LopapSpecCh1,LopapSpecCh2] ...
= importLopap(LopapParameters.MeasurementDate,fileLocation);
end
%% Matlab: Run mode
%
% <Adapt>:
%
% * runs without user interaction and without ability to modify data
% ranges, when set to 0; use this to publish.
% * if set to 1: stops for user interaction,
% modify data ranges is enabled.Use this to run.
LopapParameters.Interactive = 0;
%% LabBook : ChronoNotes
%
% <Adapt>: Your notes go here
%
% * Nitritstandardlösung: 1g/l
% * Nitritlösung: 0.01mg/l
% * Volumenfluss 1: 26.08min für 5ml -> 0.192ml/min
% * Volumenfluss 2: 26.33min für 5ml -> 0.190ml/min
% * Pumpgeschwindigkeit: 10U/min
% * Gasfluss: 2000ml/min
% * *16.9.20 13h30*: Start instrument with water
% * time: 2s/2s
% * time: 80ms/4ms
% * *16.9.20 13h00*: water -> R1 & R2, ABS closed
%% LabBook : Measurement Settings
%
% <Adapt>: list the Lopap settings here:
LopapParameters.concStandard = 0.01; % mg/l concentration of nitrite in standard used for calibration
LopapParameters.flowReagents = 5./[26.08 26.33]; % ml/min, enter time for 5 ml
LopapParameters.flowGas = 2000; % ml/min, take average of 10.
LopapParameters.temperature = 298; % temp of gas flow in K
LopapParameters.pressure = 1013; % pressure of gas flow in mbar
% <Adapt>: Choose a reference wavelength. Default is 650 nm.
LopapParameters.wavelengthR = 650; %nm, reference intensity
% <Adapt>: Choose a signal wavelength to calculate the absorption with time
% from the spectra
LopapParameters.wavelengthS = 534; %nm, signal intensity
%% LabBook : Define data ranges
%
% These settings are ignored if the TimesRanges.mat file exist. If you do
% not want to use the data stored in the TimesRanges.mat file, you need to
% delete it before running this section.
%
% These settings can be modified by selecting data ranges in a plot later.
% So, you do not necessarily need to type the dates and times here.
if isfile(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'])
load(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'])
elseif ~exist('TimeRange','var')
% <Adapt>: State the start and end time of the calibration. Usually
% this is included zero air with R1&2, zero air with nitrite standard,
% zero air with R1&2
TimeRange.calibrationRange = [datenum('16-Sep-2020 13:55:00') datenum('16-Sep-2020 15:15:00')];
% <Adapt>: State the time for one spectrum with R1&R2 without nitrite
TimeRange.specWithReagents = datenum('2020-09-16 14:31:06');
% <Adapt>: State the time for one spectrum with R1&R2 and nitrite
TimeRange.specWithNitrite = datenum('2020-09-16 15:10:14');
% <Adapt>: Choose a start and end time for one as many intervals that can
% be used to fit the background. Background is the signal when R1&R2 was
% attached to the Lopap with zero air entering the sampling coil.
TimeRange.backgroundFit = [datenum('16-Sep-2020 14:29:28') datenum('16-Sep-2020 14:40:04');...
datenum('16-Sep-2020 15:40:50') datenum('16-Sep-2020 15:51:06')];
% <Adapt>: State the start and end time of data that can be used to
% calculate the mean of the zero air signal without nitrate
TimeRange.MeanZeroAir = [datenum('16-Sep-2020 14:28:42') datenum('16-Sep-2020 14:38:36')];
% <Adapt>: State the start and end time of data that can be used to
% calculate the mean of the zero air signal without nitrate
TimeRange.MeanNitriteStandard = [datenum('16-Sep-2020 15:09:20') datenum('16-Sep-2020 15:18:30')];
save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'],'TimeRange');
else
disp('Nothing done: TimeRange in Workspace used')
end
%% Process: Modify data ranges
%
% Select the range to modify in the pop-up window and the number of fit
% ranges for the background fit in the 2nd pop-up window. Then, a figure
% will show, select the data ranges in the figure:
%
% * a line for calibration range, mean zero air and mean
% nitrite standard. Select the range by keeping the mouse pushed while
% moving from start to end time.
% * a line for each background region, default is 2. Select the range by
% keeping the mouse pushed while moving from start to end time.
% * a point for spectrum with reagents and spectrum with nitrite. Click
% once to select a specific time.
%
% * LopapParameters.Interactive needs to be 1 for this to run.
%
% # Run this section several times to modify several time regions
% # Press cancel in the pop-up window if you do not want to modify any
% region.
TimeRange = modifyLopapTimeRanges(TimeRange, LopapParameters, LopapDataCh1, LopapDataCh2);
%% Process: Find closest match of selected dates in data set
%
% The data is sampled at specific times with a certain interval. The
% selected or chosen dates might not necessarily exactly match the data.
% Here, closest match to the given dates are selected.
tmp.d1 = TimeRange.specWithReagents;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specWithReagents = tmp.d2(ind1,:);
tmp.d1 = TimeRange.specWithNitrite ;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specWithNitrite = tmp.d2(ind1,:);
% save current state of
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'])
save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'],'TimeRange');
%% LabBook: Data ranges used in this work
%
% This plots the
fprintf('Calibration range: %s to %s \n',datestr(TimeRange.calibrationRange(1)), datestr(TimeRange.calibrationRange(2)))
for ii=1:size(TimeRange.backgroundFit,2)
fprintf('Background range: %s to %s \n',datestr(TimeRange.backgroundFit(ii,1)), datestr(TimeRange.backgroundFit(ii,2)))
end
fprintf('Zero mean: %s to %s \n',datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)))
fprintf('Nitrite mean: %s to %s \n',datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)))
fprintf('Zero spectrum: %s \n',datestr(TimeRange.specWithReagents(1)))
fprintf('Nitrite spectrum: %s \n',datestr(TimeRange.specWithNitrite(1)))
%% Plot all data
%
% Here, we plot the data as shown in the "Measurement window" of the Lopap
% software. That is, the processed $I^{*}_{i}$ data for all wavelength as
% defined in the settings of the Lopap software.
%
%
% For every measurement wavelength I_i the program calculates the
% measurement (absorption) value I_i^* from the spectral intensity on the
% measurement wavelength I_i and a reference wavelength I_r with the
% following operation:
%
%
% $$I^{*}_{i} = log_{10} (\frac{I_r}{I_i})$$
%
%
%
% The plots show time on the x-axis and the ABS intensity on the y axis for
% each wavelength in separate graphs. One Figure for each Channel.
figure( 'Name', 'LopapChannel1','NumberTitle','off','Position',[10 100 550 700]);
ch1s=stackedplot(LopapDataCh1);
title(LopapParameters.instrumentReferenceWavelength(1))
for ii = 1:width(LopapDataCh1)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch1s.AxesProperties(ii).YLimits = [-0.4 -0.1];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2','NumberTitle','off','Position',[600 100 550 700]);
ch2s=stackedplot(LopapDataCh2);
title(LopapParameters.instrumentReferenceWavelength(2))
for ii = 1:width(LopapDataCh2)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties(ii).YLimits = [-0.2 0.3];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Labbook: Notes to Figures
%
% * Kanal 2 ist stärker am rauschen.
% * Ev. wegen Leck am T-Stück in der Absorptionskammer.
%% Plot selected data
%
% Here we plot the $I^{*}_{i}$ data for one wavelength only.
figure( 'Name', 'LopapChannel1 550','NumberTitle','off','Position',[10 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure.
ch1s=stackedplot(LopapDataCh1,{'550.7'});
title(LopapParameters.instrumentReferenceWavelength(1))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
%ch1s.AxesProperties.YLimits = [-0.4 -0.1];
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2 550','NumberTitle','off','Position',[600 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure
ch2s=stackedplot(LopapDataCh2,{'550.5'});
title(LopapParameters.instrumentReferenceWavelength(2))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties.YLimits = [-0.2 0.3];
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap lopap* TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot several spectra
%
% * A spectrum is obtained from the spectrometer.
% * This spectrum is averaged from multiple raw spectra with a given
% integration time. Optionally the raw spectra can be boxcar averaged with
% an adjustable boxcar window width. The recommended boxcar window width is
% 10 Pixel.
% * Then a background correction is performed: The mean of the values of
% the background pixels, which contains no measurement signal, is
% subtracted from the whole spectrum. The background pixels are freely
% configurable.
clc
close all
S = timerange(datestr(TimeRange.calibrationRange(1)), datestr(TimeRange.calibrationRange(2)));
figure( 'Name', 'LopapSpectraWithTime','NumberTitle','off')
T = LopapSpecCh1(S,:);
plot3(repmat(str2double(T.Properties.VariableNames),[height(T),1]),repmat(T.Time,[1,width(T)]),T.Variables)
% Set the viewing angle and the axis limits
view(-25, 1)
xlim([400 700])
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('Time')
zlabel('Intensity')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot selected spectra
%
% Here, we plot a spectrum without ('R1&R2') and with nitrite
% ('R1+nitrite&R2'') for each channel.
clc
close all
figure('Name', 'LopapSelectedSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specWithReagents),:},...
'DisplayName','#1 R1&R2') % a time where only R1&R2 was connected to the Lopap
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specWithReagents),:},...
'DisplayName','#2 R1&R2','LineStyle','--') % a time where only R1&R2 was connected to the Lopap
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specWithNitrite),:},...
'DisplayName','#1 R1+nitrite&R2') % a time where nitrite in R1 was connected to the Lopap
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specWithNitrite),:},...
'DisplayName','#2 R1+nitrite&R2','LineStyle','--') % a time where nitrite in R1 was connected to the Lopap
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('Intensity')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process: Absorption spectrum
%
% Here, we reprocess the data to derive an absorption spectrum. This should
% match the spectra shown for nitrite in the Lopap manual.
%
%
%
% <html>
% <img src="../Bildschirmfoto2020-09-16um14.30.03.png" width="550">
% </html>
%
clc
close all
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
figure('Name', 'LopapAbsorbtionSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),log10(LopapSpecCh1{datestr(TimeRange.specWithNitrite),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specWithNitrite),:})-log10(LopapSpecCh1{datestr(TimeRange.specWithReagents),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specWithReagents),:}),'DisplayName','#1 nitrite')
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),log10(LopapSpecCh2{datestr(TimeRange.specWithNitrite),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specWithNitrite),:})-log10(LopapSpecCh2{datestr(TimeRange.specWithReagents),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specWithReagents),:}),'DisplayName','#2 nitrite','LineStyle','--')
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('ABS (log I_{r}/I_{i}')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Absorption with time
%
% Here, we reprocess the data to get the intensity with time for 1 chosen
% signal wavelength and 1 reference wavelength. The data is identical to that
% shown in the "measurement window" of the Lopap software.
%
% The data is written to the _LopapTimeT_ time table. Rows hold the time,
% columns hold the processed signal for channel 1 and channel 2.
%
% It is important that the same signal wavelength and reference wavelength
% are taken for the calibration and for the experimental data.
clc
close all
indexS(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexS(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
LopapTimeT = timetable(LopapSpecCh1.Time,log10(LopapSpecCh1{:,indexR(1)}./LopapSpecCh1{:,indexS(1)}),log10(LopapSpecCh2{:,indexR(2)}./LopapSpecCh1{:,indexS(2)}));
LopapTimeT.Properties.VariableNames ={'Ch #1','Ch #2'};
figure('Name', 'Lopap TimeTable')
stackedplot(LopapTimeT,{{'Ch #1','Ch #2'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Background
%
clc
close all
LopapTimeT = fitLopapBkgrd(LopapTimeT,'BkgrdPeriod',TimeRange.backgroundFit);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
%% Plot bkgrd corrected data
%
clc
close all
% Plot
figure('Name', 'Lopap TimeTable Bkgrd')
s=stackedplot(LopapTimeT,{{'Ch #1 Bkgrd','Ch #2 Bkgrd'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
s.AxesProperties.YLimits = [-0.01 0.3];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
%% Process :: Calibrate
%
clc
close all
S = timerange(datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)));
TZ = LopapTimeT(S,:);
S = timerange(datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)));
TN = LopapTimeT(S,:);
calib.mean = mean(TN{:,:},1) - mean(TZ{:,:},1);
calib.std = std(TN{:,:},1);
calib.mixingRatio = LopapParameters.concStandard./1000.*LopapParameters.flowReagents./1000./46.*6.022E23./...
(246E17.*LopapParameters.flowGas).*1E12;
[LopapTimeT, calib.fitresult] = fitLopapCalibration(LopapTimeT,calib);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
fprintf('Detection limit channel #1 is %.0f ppt \n',calib.fitresult.calibChannel1(3.*calib.std(3)))
fprintf('Detection limit channel #2 is %.0f ppt \n',calib.fitresult.calibChannel2(3.*calib.std(4)))
%% Plot calibrated data
%
clc
close all
% Plot
figure('Name', 'Lopap TimeTable Calibrated')
s=stackedplot(LopapTimeT,{{'Ch #1 Calib','Ch #2 Calib'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
s.AxesProperties.YLimits = [-100 2500];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Save the data
%
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'])
save(['Lopap' datestr(LopapParameters.MeasurementDate,'yyyyMMdd_hhmm') '.mat'],'LopapParameters','LopapTimeT','TimeRange')
%% Local functions
function [fileLocation, myColorMap] = localInitializeMatlab
%INITIALIZEMATLAB Sets Lopap specific settings
% Adds processing functions to path and defines file locations for mac, windows,
% and Matlab online (which runs on unix)
% change working directory
path = matlab.desktop.editor.getActiveFilename;
[filepath,~,~] = fileparts(path);
cd(filepath);
% Data Server Address
if ismac
fileLocation.server = '/Volumes/';
fileLocation.dataPath = 'data_2411/Lopap/Marc/';
% add Lopap functions and general add-ons to path
addpath('/Users/bartels-rausch/Documents/SwitchDrive/PSI_Matlab/Add-Ons');
addpath('/Users/bartels-rausch/Documents/SwitchDrive/PSI_Matlab/Lopap/1.LOPAP_Restliche Daten');
elseif ispc
fileLocation.server = '\\fs00\'; %UNC name for file server that hosts Napp shared drive
fileLocation.dataPath = 'data_2411\Lopap\Marc\';
% add Lopap functions and general add-ons to path
addpath('C:\Users\iezzi_l\switchdrive\Add-Ons');
addpath('\\fs02\5601_surface_chemistry\People\Marc\Switchdrive\PSI_Matlab\Lopap\1.LOPAP_Restliche Daten');
elseif isunix
fileLocation.server = '/';
fileLocation.dataPath = 'Lopap/';
end
% Color maps
myColorMap.red = colormap(autumn(7));
myColorMap.blue = colormap(winter(7));
close all
end
Binary file not shown.
+99
View File
@@ -0,0 +1,99 @@
%% Compile LOPAP data
% from different data set
%% Choose data
%
% * Choose the mat files with the Lopap timetables to be plotted
% * For each data set, define the 0-time. This could be the time when you
% turned on the UV, or when the flow was directed via the flow tube. This
% depends on what you want to compare and plot
%
% The data is then imported and then transferred to relative time
% clear variables in workspace, clear text in command window
clearvars
clc
% set the path to the data and define some custom colours for plotting
[fileLocation, myColorMap] = localInitializeMatlab;
% <Adapt>: data files to be loaded
tmpFileList = {'Lopap20211201_1033.mat';...
'Lopap20211115_1431.mat'};
% <Adapt>: define 0-time for each data set to be loaded. This needs one
% date for each data set.
tmpRefTime = [datetime('01-Dec-2021 15:24:23');...
datetime('15-Nov-2021 16:49:23')];
% <Adapt>: define time window to be plotted, in seconds
tmpPlotTime = -30*60:2:120*60; % 30 minutes prior, 2 s interval, 2 hours after event
%% Load data and process
%
for ii = 1: length(tmpFileList)
% load data from mat file
load(tmpFileList{ii},'LopapTimeT')
% get time
tmpX = duration(LopapTimeT.Time-tmpRefTime(ii));
% interpolate HONO data (to have same time axis)
tmpY = interp1(seconds(tmpX),LopapTimeT.HONO,...
tmpPlotTime);
if ~exist('LopapData','var') % make new timetable if timetable does not exist
% arrange in table
LopapData = timetable(minutes(tmpPlotTime./60)',tmpY');
LopapData = renamevars(LopapData,'Var1',tmpFileList{ii});
else % or add new variable to existing timetable
LopapData = addvars(LopapData,tmpY','NewVariableNames',tmpFileList{ii});
end
end
% tidy up
clearvars ii tmp* LopapTimeT
%% Plot data
figure
stackedplot(LopapData)
title('HONO in ppt')
%% Local functions
function [fileLocation, myColorMap] = localInitializeMatlab
%INITIALIZEMATLAB Sets Lopap specific settings
% Adds processing functions to path and defines file locations for mac, windows,
% and Matlab online (which runs on unix)
% change working directory
path = matlab.desktop.editor.getActiveFilename;
[filepath,~,~] = fileparts(path);
cd(filepath);
% Data Server Address
if ismac
fileLocation.server = '/Volumes/';
fileLocation.dataPath = 'data_2411/Lopap/';
% add Lopap functions and general add-ons to path
addpath('/Users/bartels-rausch/Documents/SwitchDrive/PSI_Matlab/Add-Ons');
addpath('/Users/bartels-rausch/Documents/SwitchDrive/PSI_Matlab/Lopap/1.LOPAP_Restliche Daten');
elseif ispc
fileLocation.server = '\\fs00\'; %UNC name for file server that hosts Napp shared drive
fileLocation.dataPath = 'data_2411\Lopap\';
% add Lopap functions and general add-ons to path
addpath('C:\Users\iezzi_l\switchdrive\Add-Ons');
addpath('C:\Users\bletry_c\switchdrive\Add-Ons');
addpath('C:\Users\iezzi_l\switchdrive\Lopap\1.LOPAP_Restliche Daten');
addpath('C:\Users\bletry_c\switchdrive\Lopap\1.LOPAP_Restliche Daten');
% addpath('\\fs02\5601_surface_chemistry\People\Marc\Switchdrive\PSI_Matlab\Add-Ons');
% addpath('\\fs02\5601_surface_chemistry\People\Marc\Switchdrive\PSI_Matlab\Lopap\1.LOPAP_Restliche Daten');
elseif isunix
fileLocation.server = '/';
fileLocation.dataPath = 'Lopap/';
end
% Color maps
myColorMap.red = colormap(autumn(7));
myColorMap.blue = colormap(winter(7));
close all
end
File diff suppressed because it is too large. Load diff
Binary file not shown.
+656
View File
@@ -0,0 +1,656 @@
%% Lopap Template
% Uses data from 11.03.2021 09h51
%
%% Info: How to use this
%
% # Make a copy of this file for each Lopap measurement. Name it
% "Lopap"_date_time.m and store it in the folder with the Lopap data. %
% Remove all the general text and comments in your file and add notes and
% observations about the experiment.
%
% # You can run this script, publish it, or evaluate it section wise. For
% new data run it section wise to carefully update the settings and to
% observe the results.
% # Once the data is loaded, you can also run individual sections only.
% # Use the LabBook :: Chrono Notes section as lab book.
% # First, you need to give the date and time of the Lopap measurement to
% load and process the data.
% # All sections that can be modified are marked with the word "<Adapt>:".
% # You also need to give time regions such as for the background fit, as
% well as to calculate means. There are several options how to do this
%
% * While running the experiments and observing the signal in the Lopap
% software, you can enter the date and times directly in the LabBook :
% Select data ranges section.
% * You can manually select the regions in a plot. A window will pop up to
% select the region to modify. To modify several regions, re-run this
% section of the script. This overwrites the current settings.
% * The times of the regions are stored in the "TimeRanges.mat" file. If it
% exists, this is read; the LabBook : Select data ranges section is then
% ignored.
%% ToDO
%
% * Analyse time delay of Lopap and and to exported variables. We need this
% when looking at flow tube data.
%% Info: Output & Variables
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
% * *LopapTimeT* This is the main results table:
%
% * time
% * raw signal channel 1
% * raw signal channel 2
% * background subtracted channel 1
% * background subtracted channel 2
% * calibrated signal channel 1 in ppt
% * calibrated signal channel 2 in ppt
% * HONO in ppt (channel 2 -channel 1)
%
%
% * *TimeRangesCalib.mat*: Time ranges after they have been manually modified. If
% this file exists, the time range definitions in this script are ignored.
% * *LopapMatlabDataDATE_TIMECalib.mat*: Raw data and experimental settings.
% Background subtraction and calibration are not saved in this file. When
% this file exists, data will not be read from the windows server.
% * *LopapDATE_TIMECalib.mat*: Processed data - absorption with time (column 1
% and 2: 'Ch #1' and 'Ch #2') after background (columns 3 and 4: 'Ch #1
% Bkgrd' and 'Ch #2 Bkgrd') subtraction and calibration (columns 5 and 6:
% 'Ch #1 Calib' and 'Ch #2 Calib') for one wavelength - that can be used
% for further comparison with other experiments and discussion.
%% Matlab: Initialise
%
% This starts Matlab and defines some common settings such as the path to
% the Lopap data. Any modifications can be done at the end of this file in
% where the local function _localInitializeMatlab_ is stored.
% clear variables in workspace, clear text in command window
clearvars
clc
% set the path to the data and define some custom colours for plotting
[fileLocation, myColorMap] = localInitializeMatlab;
% <Adapt>: Date when Lopap recording was started. Defines names of Lopap
% files and folders. Seconds are irrelevant, but format is important.
LopapParameters.MeasurementDate = datenum('11-March-2021 09:51:00');
%%
% Print info to command window
disp('The data is read from:');
disp(fullfile(fileLocation.server,fileLocation.dataPath,datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM')));
disp(' ')
disp('Date of experiment:')
disp(datestr(LopapParameters.MeasurementDate));
disp(' ')
disp('Date of analyis:')
disp(datestr(today));
%% Process: Import data
%
% Reads and imports data from Lopap *.dat, *.spec, and *.prj files into:
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
if isfile(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'])
load(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'])
else
[LopapParameters.instrumentReferenceWavelength, LopapDataCh1,LopapDataCh2,LopapSpecCh1,LopapSpecCh2] ...
= importLopap(LopapParameters.MeasurementDate,fileLocation);
end
%% Matlab: Run mode
%
% <Adapt>:
%
% * runs without user interaction and without ability to modify data
% ranges, when set to 0; use this to publish.
% * if set to 1: stops for user interaction,
% modify data ranges is enabled.Use this to run.
LopapParameters.Interactive = 1;
%% LabBook : Measurement Settings % from 20.12.2020 template
%
% <Adapt>: list the Lopap settings here:
LopapParameters.concStandard = 0.01; % mg/l concentration of nitrite in standard used for calibration
LopapParameters.flowReagents = 5./[22+0/60 25+26/60]; % ml/min, enter time for 5 ml in minutes,
% channel 1 and channel 2, first number is minutes, second is seconds.
LopapParameters.flowGas = 396; % ml/min, at room temp and pressure, take average of 10.
LopapParameters.temperature = 298; % temp of gas flow in K
LopapParameters.pressure = 1013; % pressure of gas flow in mbar
LopapParameters.flowGas = LopapParameters.flowGas.*LopapParameters.pressure ./LopapParameters.temperature...
.*298.15./1013.25; % STP ml/min (as used in original LOPAP documents)
% <Adapt>: Choose a reference wavelength. Default is 650 nm.
LopapParameters.wavelengthR = 650.3; %nm, reference intensity
% <Adapt>: Choose a signal wavelength to calculate the absorption with time
% from the spectra
LopapParameters.wavelengthS = 530; %nm, signal intensity
%% LabBook : Define data ranges
%
% These settings are ignored if the TimesRanges.mat file exist. If you do
% not want to use the data stored in the TimesRanges.mat file, you need to
% delete it before running this section.
%
% These settings can be modified by selecting data ranges in a plot later.
% So, you do not necessarily need to type the dates and times here.
if LopapParameters.Interactive == 1
answer = questdlg('How would you like to set the data ranges of the experiment?', ...
'Data range selection', ...
'Read from script code','Read from mat file','Use from workspace','Read from script code');
% Handle response
switch answer
case 'Read from script code'
disp([answer 'OK. Data range will set as coded in the script.'])
tmpDataRangeSource = 1;
case 'Read from mat file'
disp([answer 'OK. Data range will loaded from mat file in current directory.'])
tmpDataRangeSource = 2;
case 'Use from workspace'
disp([answer 'OK. Data range in workspace variable will be used.'])
tmpDataRangeSource = 3;
case 'Cancel'
disp('OK. Do nothing')
tmpDataRangeSource = 0;
end
else
if isfile(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'])
tmpDataRangeSource = 2;
disp([answer ' OK. Data range will loaded from mat file in current directory.'])
elseif ~exist('TimeRange','var')
tmpDataRangeSource = 1;
disp([answer 'OK. Data range will set as coded in the script.'])
elseif exist('TimeRange','var')
tmpDataRangeSource = 3;
disp(' OK. Data range in workspace variable will be used.')
end
end
if tmpDataRangeSource == 2
load(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'])
disp('TimeRange in TimeRange*.mat file in current directory used')
disp('You need to delete the TimeRange variable in the workspace AND the TimeRange*.mat file in the current directory')
disp('to use the dates entered in the script')
disp('You keep to delete the TimeRange variable in the workspace AND delete TimeRange*.mat file in the current directory')
disp('to use the dates saved in the workspace')
elseif tmpDataRangeSource == 1
% <Adapt> (Optional - comment if no needed) state the start and end
% time of the time that you want to look at (to exclude very noisy
% signal in the beginning and at the end of the data)
TimeRange.calibrationRange = [datenum('15-Nov-2021 16:30:00') datenum('15-Nov-2021 17:20:00')];
% <Adapt>: State the time for the first spectrum that you like to
% compare (to confirm nitrite)
TimeRange.specCompareOne = datenum('15-Nov-2021 16:49:30 ');
% <Adapt>: State the time for one spectrum of the second spectrum that
% you like to compare (to confirm nitrite)
TimeRange.specCompareTwo = datenum('15-Nov-2021 17:02:00 ');
% <Adapt>: Choose a start and end time for one as many intervals that can
% be used to fit the background. Background is the signal when R1&R2 was
% attached to the Lopap with zero air entering the sampling coil.
TimeRange.backgroundFit = [datenum('15-Nov-2021 15:10:00') datenum('15-Nov-2021 15:50:00');...
datenum('0-Dec-2021 16:26:00') datenum('0-Dec-2021 16:26:00')];
% <Adapt>: State the start and end time of data that can be used to
% calculate the mean of the zero air signal without nitrite
TimeRange.MeanZeroAir = [datenum('04-Feb-2022 11:49:07') datenum('04-Feb-2022 11:59:14 ')];
% <Adapt>: State the start and end time of data that can be used to
% calculate the mean of the zero air signal with nitrite
TimeRange.MeanNitriteStandard = [datenum('04-Feb-2022 12:26:28') datenum('04-Feb-2022 12:39:46')];
save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'TimeRange');
disp('TimeRange as coded in script used')
elseif tmpDataRangeSource == 3
if exist('TimeRange','var')
disp('TimeRange in Workspace used')
disp('You need to delete the TimeRange variable in the workspace AND the TimeRange*.mat file in the current directory')
disp('to use the dates entered in the script')
elseif ~exist('TimeRange','var')
warning('TimeRange in the Workspace was selected, but there is no TimeRange variable')
end
end
% confirm that all dates are reasonable
tmp = fieldnames(TimeRange);
for ii = 1:length(tmp)
mustBeLopapReasonableDate(TimeRange.(tmp{ii}))
end
clearvars tmp*
%% Process: Modify data ranges
%
% Select the range to modify in the pop-up window and the number of fit
% ranges for the background fit in the 2nd pop-up window. Then, a figure
% will show, select the data ranges in the figure:
%
% * a line for calibration range, mean zero air and mean
% nitrite standard. Select the range by keeping the mouse pushed while
% moving from start to end time.
% * a line for each background region, default is 2. Select the range by
% keeping the mouse pushed while moving from start to end time.
% * a point for spectrum with reagents and spectrum with nitrite. Click
% once to select a specific time.
%
% * LopapParameters.Interactive needs to be 1 for this to run.
%
% # Run this section several times to modify several time regions
% # Press cancel in the pop-up window if you do not want to modify any
% region.
TimeRange = modifyLopapTimeRangesCalib(TimeRange, LopapParameters, LopapDataCh1, LopapDataCh2);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'],'TimeRange');
%% Process: Find closest match of selected dates in data set
%
% The data is sampled at specific times with a certain interval. The
% selected or chosen dates might not necessarily exactly match the data.
% Here, closest match to the given dates are selected.
tmp.d1 = TimeRange.specCompareOne;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specCompareOne = tmp.d2(ind1,:);
tmp.d1 = TimeRange.specCompareTwo ;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specCompareTwo = tmp.d2(ind1,:);
% save current state of
clearvars -except Lopap* fileLocation myColorMap TimeRange*
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'LopapDataCh1','LopapDataCh2','LopapParameters','LopapSpecCh1','LopapSpecCh2')
save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'TimeRange');
%% LabBook: Data ranges used in this work
%
clc
fprintf('Calibration range: %s to %s \n',datestr(TimeRange.calibrationRange(1)), datestr(TimeRange.calibrationRange(2)))
for ii=1:size(TimeRange.backgroundFit,1)
fprintf('Background range: %s to %s \n',datestr(TimeRange.backgroundFit(ii,1)), datestr(TimeRange.backgroundFit(ii,2)))
end
fprintf('Zero mean: %s to %s \n',datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)))
fprintf('Nitrite mean: %s to %s \n',datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)))
fprintf('Comparison Spectrum #1: %s \n',datestr(TimeRange.specCompareOne(1)))
fprintf('Comparison Spectrum #2: %s \n',datestr(TimeRange.specCompareTwo(1)))
%% Plot all data
%
% Here, we plot the data as shown in the "Measurement window" of the Lopap
% software. That is, the processed $I^{*}_{i}$ data for all wavelength as
% defined in the settings of the Lopap software.
%
%
% For every measurement wavelength I_i the program calculates the
% measurement (absorption) value I_i^* from the spectral intensity on the
% measurement wavelength I_i and a reference wavelength I_r with the
% following operation:
%
%
% $$I^{*}_{i} = log_{10} (\frac{I_r}{I_i})$$
%
%
%
% The plots show time on the x-axis and the ABS intensity on the y axis for
% each wavelength in separate graphs. One Figure for each Channel.
figure( 'Name', 'LopapChannel1','NumberTitle','off','Position',[10 100 550 700]);
ch1s=stackedplot(LopapDataCh1);
title(LopapParameters.instrumentReferenceWavelength(1))
for ii = 1:width(LopapDataCh1)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch1s.AxesProperties(ii).YLimits = [-0.4 -0.1];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2','NumberTitle','off','Position',[600 100 550 700]);
ch2s=stackedplot(LopapDataCh2);
title(LopapParameters.instrumentReferenceWavelength(2))
for ii = 1:width(LopapDataCh2)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties(ii).YLimits = [-0.2 0.3];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot all data :: selected time range
%
TR = timerange(datestr(TimeRange.calibrationRange(1)),datestr(TimeRange.calibrationRange(2)));
figure( 'Name', 'LopapChannel1','NumberTitle','off','Position',[10 100 550 700]);
ch1s=stackedplot(LopapDataCh1(TR,:));
title(LopapParameters.instrumentReferenceWavelength(1))
for ii = 1:width(LopapDataCh1)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch1s.AxesProperties(ii).YLimits = [-0.4 -0.1];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2','NumberTitle','off','Position',[600 100 550 700]);
ch2s=stackedplot(LopapDataCh2(TR,:));
title(LopapParameters.instrumentReferenceWavelength(2))
for ii = 1:width(LopapDataCh2)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties(ii).YLimits = [-0.2 0.3];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot selected data in selected time range
%
% Here we plot the $I^{*}_{i}$ data for one wavelength only.
TR = timerange(datestr(TimeRange.calibrationRange(1)),datestr(TimeRange.calibrationRange(2)));
figure( 'Name', 'LopapChannel1 550','NumberTitle','off','Position',[10 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure.
ch1s=stackedplot(LopapDataCh1(TR,:),{'550.7'});
title(LopapParameters.instrumentReferenceWavelength(1))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
%ch1s.AxesProperties.YLimits = [-0.4 -0.1];
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2 550','NumberTitle','off','Position',[600 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure
ch2s=stackedplot(LopapDataCh2(TR,:),{'550.5'});
title(LopapParameters.instrumentReferenceWavelength(2))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties.YLimits = [-0.2 0.3];
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap lopap* TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot selected spectra
%
% Here, we plot a spectrum without ('R1&R2') and with nitrite
% ('R1+nitrite&R2'') for each channel or which ever combination of spectra
% that we have selected.
clc
close all
figure('Name', 'LopapSelectedSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specCompareOne),:},...
'DisplayName',['#1 ' datestr(TimeRange.specCompareOne) ]) % a time where only R1&R2 was connected to the Lopap
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specCompareOne),:},...
'DisplayName',['#2 ' datestr(TimeRange.specCompareOne)],'LineStyle','--') % a time where only R1&R2 was connected to the Lopap
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specCompareTwo),:},...
'DisplayName',['#1 ' datestr(TimeRange.specCompareTwo)]) % a time where nitrite in R1 was connected to the Lopap
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specCompareTwo),:},...
'DisplayName',['#2 ' datestr(TimeRange.specCompareTwo)],'LineStyle','--') % a time where nitrite in R1 was connected to the Lopap
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('Intensity')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process: Absorption spectrum
%
% Here, we reprocess the data to derive an absorption spectrum. This should
% match the spectra shown for nitrite in the Lopap manual.
clc
close all
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
figure('Name', 'LopapAbsorbtionSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),log10(LopapSpecCh1{datestr(TimeRange.specCompareTwo),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specCompareTwo),:})-log10(LopapSpecCh1{datestr(TimeRange.specCompareOne),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specCompareOne),:}),'DisplayName','#1')
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),log10(LopapSpecCh2{datestr(TimeRange.specCompareTwo),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specCompareTwo),:})-log10(LopapSpecCh2{datestr(TimeRange.specCompareOne),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specCompareOne),:}),'DisplayName','#2','LineStyle','--')
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('ABS (log I_{r}/I_{i}')
title([datestr(TimeRange.specCompareOne) ' & ' datestr(TimeRange.specCompareTwo)])
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Absorption with time
%
% Here, we reprocess the data to get the intensity with time for 1 chosen
% signal wavelength and 1 reference wavelength. The data is identical to that
% shown in the "measurement window" of the Lopap software.
%
% The data is written to the _LopapTimeT_ time table. Rows hold the time,
% columns hold the processed signal for channel 1 and channel 2.
%
% It is important that the same signal wavelength and reference wavelength
% are taken for the calibration and for the experimental data.
clc
close all
indexS(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexS(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
LopapTimeT = timetable(LopapSpecCh1.Time,log10(LopapSpecCh1{:,indexR(1)}./LopapSpecCh1{:,indexS(1)}),...
log10(LopapSpecCh2{:,indexR(2)}./LopapSpecCh2{:,indexS(2)}));
LopapTimeT.Properties.VariableNames ={'Ch #1','Ch #2'};
figure('Name', 'Lopap TimeTable')
stackedplot(LopapTimeT,{{'Ch #1','Ch #2'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Background
%
% Currently a linear fit forced through 0.
clc
close all
LopapTimeT = fitLopapBkgrd(LopapTimeT,'BkgrdPeriod',TimeRange.backgroundFit);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%%
%
clc
close all
% Plot
figure('Name', 'Lopap TimeTable Bkgrd')
s=stackedplot(LopapTimeT,{{'Ch #1 Bkgrd','Ch #2 Bkgrd'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
% s.AxesProperties.YLimits = [-0.01 0.6];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
%% Process :: Calibrate
%
clc
close all
S = timerange(datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)));
TZ = LopapTimeT(S,:);
S = timerange(datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)));
TN = LopapTimeT(S,:);
calib.mean = mean(TN{:,:},1) - mean(TZ{:,:},1);
calib.std = std(TN{:,:},1);
calib.mixingRatio = LopapParameters.concStandard./1000.*LopapParameters.flowReagents./1000./46.*6.022E23./...
(246E17.*LopapParameters.flowGas).*1E12;
[LopapTimeT, calib.fitresult] = fitLopapCalibration(LopapTimeT,calib);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
fprintf('Detection limit channel #1 is %.0f ppt \n',calib.fitresult.calibChannel1(3.*calib.std(3)))
fprintf('Detection limit channel #2 is %.0f ppt \n',calib.fitresult.calibChannel2(3.*calib.std(4)))
%% Plot calibrated data
%
clc
close all
% Plot
figure('Name', 'Lopap TimeTable Calibrated')
s=stackedplot(LopapTimeT,{{'Ch #1 Calib','Ch #2 Calib'}});
title(['HONO in ppt, Lopap' datestr(LopapParameters.MeasurementDate)])
% s.AxesProperties.YLimits = [-100 2500];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Save the data
%
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'LopapDataCh1','LopapDataCh2','LopapParameters','calib','LopapSpecCh1','LopapSpecCh2')
save(['Lopap' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'LopapParameters','LopapTimeT','TimeRange')
%% Local functions
function [fileLocation, myColorMap] = localInitializeMatlab
%INITIALIZEMATLAB Sets Lopap specific settings
% Adds processing functions to path and defines file locations for mac, windows,
% and Matlab online (which runs on unix)
% change working directory
path = matlab.desktop.editor.getActiveFilename;
[filepath,~,~] = fileparts(path);
cd(filepath);
% Data Server Address
if ismac
fileLocation.server = '/Volumes/';
fileLocation.dataPath = 'data_2411/Lopap/';
% add Lopap functions and general add-ons to path
addpath('/Users/bartels-rausch/SwitchDrive/PSI_Matlab/Add-Ons');
addpath('/Users/bartels-rausch/SwitchDrive/PSI_Matlab/Lopap/1.LOPAP_Restliche Daten');
elseif ispc
fileLocation.server = '\\fs00\'; %UNC name for file server that hosts Napp shared drive
fileLocation.dataPath = 'data_2411\Lopap\';
% add Lopap functions and general add-ons to path
addpath('C:\Users\vieill_l\switchdrive\Add-Ons');
addpath('C:\Users\vieill_l\switchdrive\Lopap\1.LOPAP_Restliche Daten');
elseif isunix
fileLocation.server = '/';
fileLocation.dataPath = 'Lopap/';
end
% Color maps
myColorMap.red = colormap(autumn(7));
myColorMap.blue = colormap(winter(7));
close all
end
BIN
View File
Binary file not shown.
File diff suppressed because it is too large. Load diff
BIN
View File
Binary file not shown.
BIN
View File
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
+615
View File
@@ -0,0 +1,615 @@
%% Lopap Template
% Uses data from 20.12.2020 13h47
%
%% Info: How to use this
%
% # Make a copy of this file for each Lopap measurement. Name it
% "Lopap"_date_time.m and store it in the folder with the Lopap data. %
% Remove all the general text and comments in your file and add notes and
% observations about the experiment.
%
% # You can run this script, publish it, or evaluate it section wise. For
% new data run it section wise to carefully update the settings and to
% observe the results.
% # Once the data is loaded, you can also run individual sections only.
% # Use the LabBook :: Chrono Notes section as lab book.
% # First, you need to give the date and time of the Lopap measurement to
% load and process the data.
% # All sections that can be modified are marked with the word "<Adapt>:".
% # You also need to give time regions such as for the background fit, as
% well as to calculate means. There are several options how to do this
%
% * While running the experiments and observing the signal in the Lopap
% software, you can enter the date and times directly in the LabBook :
% Select data ranges section.
% * You can manually select the regions in a plot. A window will pop up to
% select the region to modify. To modify several regions, re-run this
% section of the script. This overwrites the current settings. You need to
% copy the output the the script to save your modified settings for the
% next time.
%% ToDO
%
% * Analyse time delay of Lopap and and to exported variables. We need this
% when looking at flow tube data.
%% Info: Output & Variables
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
% * *LopapTimeT* This is the main results table:
%
% * time
% * raw signal channel 1
% * raw signal channel 2
% * background subtracted channel 1
% * background subtracted channel 2
% * calibrated signal channel 1 in ppt
% * calibrated signal channel 2 in ppt
% * HONO in ppt (channel 2 -channel 1)
%
%
% * *TimeRangesCalib.mat*: Time ranges after they have been manually modified. If
% this file exists, the time range definitions in this script are ignored.
% * *LopapMatlabDataDATE_TIMECalib.mat*: Raw data and experimental settings.
% Background subtraction and calibration are not saved in this file. When
% this file exists, data will not be read from the windows server.
% * *LopapDATE_TIMECalib.mat*: Processed data - absorption with time (column 1
% and 2: 'Ch #1' and 'Ch #2') after background (columns 3 and 4: 'Ch #1
% Bkgrd' and 'Ch #2 Bkgrd') subtraction and calibration (columns 5 and 6:
% 'Ch #1 Calib' and 'Ch #2 Calib') for one wavelength - that can be used
% for further comparison with other experiments and discussion.
%% Matlab: Initialise
%
% This starts Matlab and defines some common settings such as the path to
% the Lopap data. Any modifications can be done at the end of this file in
% where the local function _localInitializeMatlab_ is stored.
% clear variables in workspace, clear text in command window
clearvars
clc
% set the path to the data and define some custom colours for plotting
[fileLocation, myColorMap] = localInitializeMatlab;
% <Adapt>: Date when Lopap recording was started. Defines names of Lopap
% files and folders. Seconds are irrelevant, but format is important.
LopapParameters.MeasurementDate = datenum('20-Dez-2021 13:47:00');
%%
% Print info to command window
disp('The data is read from:');
disp(fullfile(fileLocation.server,fileLocation.dataPath,datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM')));
disp(' ')
disp('Date of experiment:')
disp(datestr(LopapParameters.MeasurementDate));
disp(' ')
disp('Date of analyis:')
disp(datestr(today));
%% Process: Import data
%
% Reads and imports data from Lopap *.dat, *.spec, and *.prj files into:
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
if isfile(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'])
load(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'])
else
[LopapParameters.instrumentReferenceWavelength, LopapDataCh1,LopapDataCh2,LopapSpecCh1,LopapSpecCh2] ...
= importLopap(LopapParameters.MeasurementDate,fileLocation);
end
%% Matlab: Run mode
%
% <Adapt>:
%
% * runs without user interaction and without ability to modify data
% ranges, when set to 0; use this to publish.
% * if set to 1: stops for user interaction,
% modify data ranges is enabled.Use this to run.
LopapParameters.Interactive = 0;
%% LabBook : Measurement Settings
%
% <Adapt>: list the Lopap settings here:
LopapParameters.concStandard = 0.01; % mg/l concentration of nitrite in standard used for calibration
LopapParameters.flowReagents = 5./[22+0/60 25+26/60]; % ml/min, enter time for 5 ml in minutes,
% channel 1 and channel 2, first number is minutes, second is seconds.
LopapParameters.flowGas = 396; % ml/min, at room temp and pressure, take average of 10.
LopapParameters.temperature = 298; % temp of gas flow in K
LopapParameters.pressure = 1013; % pressure of gas flow in mbar
LopapParameters.flowGas = LopapParameters.flowGas.*LopapParameters.pressure ./LopapParameters.temperature...
.*298.15./1013.25; % STP ml/min (as used in original LOPAP documents,1 atm (1013 hPa) and 298.15K)
% <Adapt>: Choose a reference wavelength. Default is 650 nm.
LopapParameters.wavelengthR = 650.3; %nm, reference intensity
% <Adapt>: Choose a signal wavelength to calculate the absorption with time
% from the spectra
LopapParameters.wavelengthS = 530; %nm, signal intensity
%% LabBook : Define data ranges
%
% The time settings can be modified by selecting data ranges in a plot later.
% So, you do not necessarily need to type the dates and times here.
%
% 11/2022 TBR: Removed the user choice. Data ranges are now ALWAYS read
% from the script.
% <Adapt>: % State the start and end time of the calibration. Usually
% this is included zero air with R1&2, zero air with nitrite standard,
% zero air with R1&2
TimeRange.calibrationRange = [datenum('20-Dec-2021 15:00:00') datenum('20-Dec-2021 16:40:00')];
% <Adapt>: State the time for the first spectrum that you like to
% compare (to confirm nitrite)
TimeRange.specCompareOne = datenum('20-Dec-2021 15:23:00');
% <Adapt>: State the time for one spectrum of the second spectrum that
% you like to compare (to confirm nitrite)
TimeRange.specCompareTwo = datenum('20-Dec-2021 16:15:00');
% <Adapt>: Choose a start and end time for one as many intervals that can
% be used to fit the background. Background is the signal when R1&R2 was
% attached to the Lopap with zero air entering the sampling coil.
TimeRange.backgroundFit = [datenum('20-Dec-2021 15:10:58') datenum('20-Dec-2021 15:31:18'); datenum('20-Dec-2021 16:30:51') datenum('20-Dec-2021 16:25:11')];
% <Adapt>: State the start and end time of data that can be used to
% calculate the mean of the zero air signal without nitrite
TimeRange.MeanZeroAir = [datenum('20-Dec-2021 15:25:00') datenum('20-Dec-2021 15:27:00')];
% <Adapt>: State the start and end time of data that can be used to
% calculate the mean of the zero air signal with nitrite
TimeRange.MeanNitriteStandard = [datenum('20-Dec-2021 15:47:00') datenum('20-Dec-2021 15:49:00')];
% confirm that all dates are reasonable
tmp = fieldnames(TimeRange);
for ii = 1:length(tmp)
mustBeLopapReasonableDate(TimeRange.(tmp{ii}))
end
clearvars tmp*
%% Process: Modify data ranges
%
% Select the range to modify in the pop-up window and the number of fit
% ranges for the background fit in the 2nd pop-up window. Then, a figure
% will show, select the data ranges in the figure:
%
% * a line for calibration range, mean zero air and mean
% nitrite standard. Select the range by keeping the mouse pushed while
% moving from start to end time.
% * a line for each background region, default is 2. Select the range by
% keeping the mouse pushed while moving from start to end time.
% * a point for spectrum with reagents and spectrum with nitrite. Click
% once to select a specific time.
%
% * LopapParameters.Interactive needs to be 1 for this to run.
%
% # Run this section several times to modify several time regions
% # Press cancel in the pop-up window if you do not want to modify any
% region.
TimeRange = modifyLopapTimeRangesCalib(TimeRange, LopapParameters, LopapDataCh1, LopapDataCh2);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
% save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'],'TimeRange');
%% Process: Find closest match of selected dates in data set
%
% The data is sampled at specific times with a certain interval. The
% selected or chosen dates might not necessarily exactly match the data.
% Here, closest match to the given dates are selected.
tmp.d1 = TimeRange.specCompareOne;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specCompareOne = tmp.d2(ind1,:);
tmp.d1 = TimeRange.specCompareTwo ;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specCompareTwo = tmp.d2(ind1,:);
% save current state of
clearvars -except Lopap* fileLocation myColorMap TimeRange*
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'LopapDataCh1','LopapDataCh2','LopapParameters','LopapSpecCh1','LopapSpecCh2')
% save(['TimeRanges' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'TimeRange');
%% LabBook: Data ranges used in this work
%
fprintf('Calibration range: %s to %s \n',datestr(TimeRange.calibrationRange(1)), datestr(TimeRange.calibrationRange(2)))
for ii=1:size(TimeRange.backgroundFit,1)
fprintf('Background range: %s to %s \n',datestr(TimeRange.backgroundFit(ii,1)), datestr(TimeRange.backgroundFit(ii,2)))
end
fprintf('Zero mean: %s to %s \n',datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)))
fprintf('Nitrite mean: %s to %s \n',datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)))
fprintf('Comparison Spectrum #1: %s \n',datestr(TimeRange.specCompareOne(1)))
fprintf('Comparison Spectrum #2: %s \n',datestr(TimeRange.specCompareTwo(1)))
%% Plot all data
%
% Here, we plot the data as shown in the "Measurement window" of the Lopap
% software. That is, the processed $I^{*}_{i}$ data for all wavelength as
% defined in the settings of the Lopap software.
%
%
% For every measurement wavelength I_i the program calculates the
% measurement (absorption) value I_i^* from the spectral intensity on the
% measurement wavelength I_i and a reference wavelength I_r with the
% following operation:
%
%
% $$I^{*}_{i} = log_{10} (\frac{I_r}{I_i})$$
%
%
%
% The plots show time on the x-axis and the ABS intensity on the y axis for
% each wavelength in separate graphs. One Figure for each Channel.
figure( 'Name', 'LopapChannel1','NumberTitle','off','Position',[10 100 550 700]);
stackedplot(LopapDataCh1);
title(LopapParameters.instrumentReferenceWavelength(1))
for ii = 1:width(LopapDataCh1)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch1s.AxesProperties(ii).YLimits = [-0.4 -0.1];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2','NumberTitle','off','Position',[600 100 550 700]);
stackedplot(LopapDataCh2);
title(LopapParameters.instrumentReferenceWavelength(2))
for ii = 1:width(LopapDataCh2)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties(ii).YLimits = [-0.2 0.3];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot all data :: selected time range
%
TR = timerange(datestr(TimeRange.calibrationRange(1)),datestr(TimeRange.calibrationRange(2)));
figure( 'Name', 'LopapChannel1','NumberTitle','off','Position',[10 100 550 700]);
stackedplot(LopapDataCh1(TR,:));
title(LopapParameters.instrumentReferenceWavelength(1))
for ii = 1:width(LopapDataCh1)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch1s.AxesProperties(ii).YLimits = [-0.4 -0.1];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2','NumberTitle','off','Position',[600 100 550 700]);
stackedplot(LopapDataCh2(TR,:));
title(LopapParameters.instrumentReferenceWavelength(2))
for ii = 1:width(LopapDataCh2)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties(ii).YLimits = [-0.2 0.3];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot selected data in selected time range
%
% Here we plot the $I^{*}_{i}$ data for one wavelength only.
TR = timerange(datestr(TimeRange.calibrationRange(1)),datestr(TimeRange.calibrationRange(2)));
figure( 'Name', 'LopapChannel1 550','NumberTitle','off','Position',[10 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure.
ch1s=stackedplot(LopapDataCh1(TR,:),{'550.7'});
title(LopapParameters.instrumentReferenceWavelength(1))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
%ch1s.AxesProperties.YLimits = [-0.4 -0.1];
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2 550','NumberTitle','off','Position',[600 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure
ch2s=stackedplot(LopapDataCh2(TR,:),{'550.5'});
title(LopapParameters.instrumentReferenceWavelength(2))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties.YLimits = [-0.2 0.3];
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap lopap* TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot selected spectra
%
% Here, we plot a spectrum without ('R1&R2') and with nitrite
% ('R1+nitrite&R2'') for each channel or which ever combination of spectra
% that we have selected.
clc
close all
figure('Name', 'LopapSelectedSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specCompareOne),:},...
'DisplayName',['#1 ' datestr(TimeRange.specCompareOne) ]) % a time where only R1&R2 was connected to the Lopap
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specCompareOne),:},...
'DisplayName',['#2 ' datestr(TimeRange.specCompareOne)],'LineStyle','--') % a time where only R1&R2 was connected to the Lopap
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specCompareTwo),:},...
'DisplayName',['#1 ' datestr(TimeRange.specCompareTwo)]) % a time where nitrite in R1 was connected to the Lopap
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specCompareTwo),:},...
'DisplayName',['#2 ' datestr(TimeRange.specCompareTwo)],'LineStyle','--') % a time where nitrite in R1 was connected to the Lopap
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('Intensity')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process: Absorption spectrum
%
% Here, we reprocess the data to derive an absorption spectrum. This should
% match the spectra shown for nitrite in the Lopap manual.
clc
close all
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
figure('Name', 'LopapAbsorbtionSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),log10(LopapSpecCh1{datestr(TimeRange.specCompareTwo),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specCompareTwo),:})-log10(LopapSpecCh1{datestr(TimeRange.specCompareOne),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specCompareOne),:}),'DisplayName','#1')
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),log10(LopapSpecCh2{datestr(TimeRange.specCompareTwo),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specCompareTwo),:})-log10(LopapSpecCh2{datestr(TimeRange.specCompareOne),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specCompareOne),:}),'DisplayName','#2','LineStyle','--')
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('ABS (log I_{r}/I_{i}')
title([datestr(TimeRange.specCompareOne) ' & ' datestr(TimeRange.specCompareTwo)])
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Absorption with time
%
% Here, we reprocess the data to get the intensity with time for 1 chosen
% signal wavelength and 1 reference wavelength. The data is identical to that
% shown in the "measurement window" of the Lopap software.
%
% The data is written to the _LopapTimeT_ time table. Rows hold the time,
% columns hold the processed signal for channel 1 and channel 2.
%
% It is important that the same signal wavelength and reference wavelength
% are taken for the calibration and for the experimental data.
clc
close all
indexS(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexS(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
LopapTimeT = timetable(LopapSpecCh1.Time,log10(LopapSpecCh1{:,indexR(1)}./LopapSpecCh1{:,indexS(1)}),...
log10(LopapSpecCh2{:,indexR(2)}./LopapSpecCh2{:,indexS(2)}));
LopapTimeT.Properties.VariableNames ={'Ch #1','Ch #2'};
figure('Name', 'Lopap TimeTable')
stackedplot(LopapTimeT,{{'Ch #1','Ch #2'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Background
%
% Currently a linear fit forced through 0.
clc
close all
LopapTimeT = fitLopapBkgrd(LopapTimeT,'BkgrdPeriod',TimeRange.backgroundFit,'TimeRangeToPlot',TimeRange.calibrationRange);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%%
%
clc
close all
% Plot
TR = timerange(datestr(TimeRange.calibrationRange(1)),datestr(TimeRange.calibrationRange(2)));
figure('Name', 'Lopap TimeTable Bkgrd')
stackedplot(LopapTimeT(TR,:),{{'Ch #1 Bkgrd','Ch #2 Bkgrd'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
% s.AxesProperties.YLimits = [-0.01 0.6];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
%% Process :: Calibrate
%
clc
close all
S = timerange(datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)));
TZ = LopapTimeT(S,:);
S = timerange(datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)));
TN = LopapTimeT(S,:);
calib.mean = mean(TN{:,:},1) - mean(TZ{:,:},1);
calib.std = std(TN{:,:},1);
calib.stdZero = std(TZ{:,:},1);
calib.mixingRatio = LopapParameters.concStandard./1000.*LopapParameters.flowReagents./1000./46.*6.022E23./...
(246E17.*LopapParameters.flowGas).*1E12; % 246E17 are molecules in 1 ml at 298.15 K and 1013 hPa thus flow needs to be in these standard conditions!
[LopapTimeT, calib.fitresult] = fitLopapCalibration(LopapTimeT,calib);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
disp("***************************************************************")
disp("***************************************************************")
fprintf('Detection limit channel #1 is %.0f ppt \n',calib.fitresult.calibChannel1(3.*calib.stdZero(3)))
fprintf('Detection limit channel #2 is %.0f ppt \n',calib.fitresult.calibChannel2(3.*calib.stdZero(4)))
disp("***************************************************************")
fprintf('Calibration range: %s to %s \n',datestr(TimeRange.calibrationRange(1)), datestr(TimeRange.calibrationRange(2)))
for ii=1:size(TimeRange.backgroundFit,1)
fprintf('Background range: %s to %s \n',datestr(TimeRange.backgroundFit(ii,1)), datestr(TimeRange.backgroundFit(ii,2)))
end
fprintf('Zero mean: %s to %s \n',datestr(TimeRange.MeanZeroAir(1)), datestr(TimeRange.MeanZeroAir(2)))
fprintf('Nitrite mean: %s to %s \n',datestr(TimeRange.MeanNitriteStandard(1)), datestr(TimeRange.MeanNitriteStandard(2)))
%% Plot calibrated data
%
close all
% Plot
TR = timerange(datestr(TimeRange.calibrationRange(1)),datestr(TimeRange.calibrationRange(2)));
figure('Name', 'Lopap TimeTable Calibrated')
s=stackedplot(LopapTimeT(TR,:),{{'Ch #1 Calib','Ch #2 Calib'}});
title(['HONO in ppt, Lopap' datestr(LopapParameters.MeasurementDate)])
% s.AxesProperties.YLimits = [-100 2500];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Save the data
%
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'LopapDataCh1','LopapDataCh2','LopapParameters','calib','LopapSpecCh1','LopapSpecCh2')
save(['Lopap' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') 'Calib.mat'],'LopapParameters','LopapTimeT','TimeRange')
%% Local functions
function [fileLocation, myColorMap] = localInitializeMatlab
%INITIALIZEMATLAB Sets Lopap specific settings
% Adds processing functions to path and defines file locations for mac, windows,
% and Matlab online (which runs on unix)
% change working directory
path = matlab.desktop.editor.getActiveFilename;
[filepath,~,~] = fileparts(path);
cd(filepath);
% Data Server Address
if ismac
fileLocation.server = '/Volumes/';
fileLocation.dataPath = 'data_2411/Lopap/';
% add Lopap functions and general add-ons to path
addpath('/Users/bartels-rausch/SwitchDrive/PSI_Matlab/Add-Ons');
addpath('/Users/bartels-rausch/SwitchDrive/PSI_Matlab/Lopap/1.LOPAP_Restliche Daten');
elseif ispc
fileLocation.server = '\\fs00\'; %UNC name for file server that hosts Napp shared drive
fileLocation.dataPath = 'data_2411\Lopap\';
% add Lopap functions and general add-ons to path
addpath('C:\Users\iezzi_l\switchdrive\Add-Ons');
addpath('C:\Users\vieill_l\switchdrive\Add-Ons');
addpath('C:\Users\bletry_c\switchdrive\Add-Ons');
addpath('C:\Users\iezzi_l\switchdrive\Lopap\1.LOPAP_Restliche Daten');
addpath('C:\Users\vieill_l\switchdrive\Lopap\1.LOPAP_Restliche Daten');
addpath('C:\Users\bletry_c\switchdrive\Lopap\1.LOPAP_Restliche Daten');
addpath('U:\SwitchDrive\Lopap\1.LOPAP_Restliche Daten');
addpath('U:\SwitchDrive\Add-Ons');
elseif isunix
fileLocation.server = '/';
fileLocation.dataPath = 'Lopap/';
end
% Color maps
myColorMap.red = colormap(autumn(7));
myColorMap.blue = colormap(winter(7));
close all
end
+785
View File
@@ -0,0 +1,785 @@
%% Lopap Template
% Uses data from 20211115_1431
%
%% Info: How to use this
%
% # Make a copy of this file for each Lopap measurement. Name it
% "Lopap"_date_time.m and store it in the folder with the Lopap data. %
% Remove all the general text and comments in your file and add notes and
% observations about the experiment.
%
% # You can run this script, publish it, or evaluate it section wise. For
% new data run it section wise to carefully update the settings and to
% observe the results.
% # Once the data is loaded, you can also run individual sections only.
% # You may use the LabBook :: Chrono Notes section as lab book.
% # All sections that can be modified are marked with the word "<Adapt>:".
% # You also need to give time regions such as for the background fit, as
% well as to calculate means. There are several options how to do this
%
% * While running the experiments and observing the signal in the Lopap
% software, you can enter the date and times directly in the LabBook :
% Select data ranges section.
% * You can manually select the regions in a plot. A window will pop up to
% select the region to modify. To modify several regions, re-run this
% section of the script. This overwrites the current settings.
% * The times of the regions are stored in the "TimeRanges.mat" file. If it
% exists, this is read; the LabBook : Select data ranges section is then
% ignored.
%% Info: Output & Variables
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
% * *LopapTimeT* This is the main results table:
%
% * time
% * raw signal channel 1
% * raw signal channel 2
% * background subtracted channel 1
% * background subtracted channel 2
% * calibrated signal channel 1 in ppt
% * calibrated signal channel 2 in ppt
% * HONO in ppt (channel 2 -channel 1)
%
%
% * *TimeRanges.mat*: Time ranges after they have been manually modified. If
% this file exists, the time range definitions in this script are ignored.
% * *LopapMatlabDataDATE_TIME.mat*: Raw data and experimental settings.
% Background subtraction and calibration are not saved in this file. When
% this file exists, data will not be read from the windows server.
% * *LopapDATE_TIME.mat*: Processed data - absorption with time (column 1
% and 2: 'Ch #1' and 'Ch #2') after background (columns 3 and 4: 'Ch #1
% Bkgrd' and 'Ch #2 Bkgrd') subtraction and calibration (columns 5 and 6:
% 'Ch #1 Calib' and 'Ch #2 Calib') for one wavelength - that can be used
% for further comparison with other experiments and discussion. Also
% exported is HONO in ppt (channel2-channel1 after calibration).
%
%% Matlab: Initialise
%
% This starts Matlab and defines some common settings such as the path to
% the Lopap data. Any modifications can be done at the end of this file in
% where the local function _localInitializeMatlab_ is stored.
% clear variables in workspace, clear text in command window
clearvars
clc
% set the path to the data and define some custom colours for plotting
[fileLocation, myColorMap] = localInitializeMatlab;
% <Adapt>: Date when Lopap recording was started. Defines names of Lopap
% files and folders. Seconds are irrelevant, but format is important.
LopapParameters.MeasurementDate = datenum('15-Nov-2021 14:31:00');
% <Adapt>: Date when Lopap calibration for this measurement was started.
% Defines names of Lopap files and folders. Seconds are irrelevant, but
% format is important.
LopapParameters.CalibrationDate = datenum('11-March-2021 09:51:00');
%%
% Print info to command window
disp('The data is read from:');
disp(fullfile(fileLocation.server,fileLocation.dataPath,datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM')));
disp(' ')
disp('Date of experiment:')
disp(datestr(LopapParameters.MeasurementDate));
try
disp(' ')
disp('Date of calibration:')
disp(datestr(LopapParameters.CalibrationDate));
catch
end
disp(' ')
disp('Date of analyis:')
disp(datestr(today));
%% Process: Import data
%
% Reads and imports data from Lopap *.dat, *.spec, and *.prj files and the
% calibration results and settings from a *.mat file (*.mat file needs to
% be in same folder as this file).
%
% * *lopapReferenceWavelength*: This _string variable_ holds text with
% stating the reference wavelength for channel 1 and 2 used to process the
% data by the Lopap software.
% * *LopapDataCh1 & LopapDataCh2*: These _time tables_ hold the data that
% have been processed by the Lopap Software during the data acquisition. It
% contains the absorption for a number
% of selected wavelength. The absorption is the log of the intensity at the
% reference wavelength relative to the intensity at each wavelength of
% interest. The timetables hold the time in each row and the wavelength in
% the columns. The content is thus identical to that shown in the
% "measurement window" of the Lopap software.
% * *LopapSpecCh1,LopapSpecCh2*: The _time tables_ hold the full spectrum
% for each time. Time is in the rows, wavelength in the columns. The
% wavelength are read from the *.prj file, the spectra data from the *.spec
% file. The information is identical to that shown in the "spectrum window"
% of the Lopap software.
% * *LopapParametersCalibration*: Holds the date, wavelength,
% concentrations, flows, temperature & pressure of the calibration
% measurement.
% * *LopapCalibCalibration* holds the mean, standard deviation, mixing ratio
% and fit results of the calibration. The fit results can be used to
% calibrate this data-set. It's a linear fit via 0 of signal intensity vs.
% HONO in ppt.
if isfile(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'])
load(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'])
disp('Lopap data loaded from mat file')
else
[LopapParameters.instrumentReferenceWavelength, LopapDataCh1,LopapDataCh2,LopapSpecCh1,LopapSpecCh2] ...
= importLopap(LopapParameters.MeasurementDate,fileLocation);
disp('Lopap data imported from intrument')
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'],'LopapDataCh1','LopapDataCh2','LopapParameters','LopapSpecCh1','LopapSpecCh2')
end
% if the calibration file exist, it will read the results and the
% reference wavelength and signal wavelength from that file
if isfile(['LopapMatlabData' datestr(LopapParameters.CalibrationDate,'yyyymmdd_HHMM') 'Calib.mat'])
LopapParametersCalibration = load(['LopapMatlabData' datestr(LopapParameters.CalibrationDate,'yyyymmdd_HHMM') 'Calib.mat'],'LopapParameters');
LopapCalibCalibration = load(['LopapMatlabData' datestr(LopapParameters.CalibrationDate,'yyyymmdd_HHMM') 'Calib.mat'],'calib');
LopapParameters.wavelengthR = LopapParametersCalibration.LopapParameters.wavelengthR;
LopapParameters.wavelengthS = LopapParametersCalibration.LopapParameters.wavelengthS;
disp('Calibration data loaded from mat file')
else
% otherwise, you can choose these wavelength here
disp('No calibration done yet!')
% <Adapt>: Choose a reference wavelength. Default is 650 nm.
LopapParameters.wavelengthR = 650; %nm, reference intensity
% <Adapt>: Choose a signal wavelength to calculate the absorption with time
% from the spectra
LopapParameters.wavelengthS = 530; %nm, signal intensity
end
%% Matlab: Run mode
%
% <Adapt>:
%
% * runs without user interaction and without ability to modify data
% ranges, when set to 0; use this to publish.
% * if set to 1: stops for user interaction,
% modify data ranges is enabled.Use this to run.
LopapParameters.Interactive = 1;
%% LabBook : Define data ranges
%%
% These settings can be modified by selecting data ranges in a plot later.
% So, you do not necessarily need to type the dates and times here.
% <Adapt> (Optional - comment if no needed) state the start and end
% time of the time that you want to look at (to exclude very noisy
% signal in the beginning and at the end of the data)
TimeRange.dataRange = [datenum('15-Nov-2021 16:30:00') datenum('15-Nov-2021 17:20:00')];
% <Adapt>: State the time for the first spectrum that you like to
% compare (to confirm nitrite)
TimeRange.specCompareOne = datenum('15-Nov-2021 16:49:30');
% <Adapt>: State the time for one spectrum of the second spectrum that
% you like to compare (to confirm nitrite)
TimeRange.specCompareTwo = datenum('15-Nov-2021 17:02:00');
% <Adapt>: Choose a start and end time for one as many intervals that can
% be used to fit the background. Background is the signal when R1&R2 was
% attached to the Lopap with zero air entering the sampling coil.
TimeRange.backgroundFit = [datenum('15-Nov-2021 15:10:00') datenum('15-Nov-2021 15:50:00')];% ;...
% datenum('16-Sep-2020 15:40:50') datenum('16-Sep-2020 15:51:06')];
% confirm that all dates are reasonable
tmp = fieldnames(TimeRange);
for ii = 1:length(tmp)
mustBeLopapReasonableDate(TimeRange.(tmp{ii}))
end
clearvars tmp*
%% Process: Modify data ranges
%
% Select the range to modify in the pop-up window and the number of fit
% ranges for the background fit in the 2nd pop-up window. Then, a figure
% will show, select the data ranges in the figure:
%
% * a line for calibration range, mean zero air and mean
% nitrite standard. Select the range by keeping the mouse pushed while
% moving from start to end time.
% * a line for each background region, default is 2. Select the range by
% keeping the mouse pushed while moving from start to end time.
% * a point for spectrum with reagents and spectrum with nitrite. Click
% once to select a specific time.
%
% * LopapParameters.Interactive needs to be 1 for this to run.
%
% # Run this section several times to modify several time regions
% # Press cancel in the pop-up window if you do not want to modify any
% region.
TimeRange = modifyLopapTimeRanges(TimeRange, LopapParameters, LopapDataCh1, LopapDataCh2);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
%% Process: Find closest match of selected dates in data set
%
% The data is sampled at specific times with a certain interval. The
% selected or chosen dates might not necessarily exactly match the data.
% Here, closest match to the given dates are selected.
tmp.d1 = TimeRange.specCompareOne;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specCompareOne = tmp.d2(ind1,:);
tmp.d1 = TimeRange.specCompareTwo ;
tmp.d2 = LopapSpecCh1.Time;
[~,ind1] = min(abs(bsxfun(@minus,datenum(tmp.d2),datenum(tmp.d1)')));
TimeRange.specCompareTwo = tmp.d2(ind1,:);
% save current state of
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'],'LopapDataCh1','LopapDataCh2','LopapParameters','LopapSpecCh1','LopapSpecCh2')
%% LabBook: Data ranges used in this work
%
% This plots the
for ii=1:size(TimeRange.backgroundFit,1)
fprintf('Background range: %s to %s \n',datestr(TimeRange.backgroundFit(ii,1)), datestr(TimeRange.backgroundFit(ii,2)))
end
fprintf('Comparision Spectrum 1: %s \n',datestr(TimeRange.specCompareOne(1)))
fprintf('Comparision Spectrum 2: %s \n',datestr(TimeRange.specCompareTwo(1)))
%% Plot all data
%
% Here, we plot the data as shown in the "Measurement window" of the Lopap
% software. That is, the processed $I^{*}_{i}$ data for all wavelength as
% defined in the settings of the Lopap software.
%
%
% For every measurement wavelength I_i the program calculates the
% measurement (absorption) value I_i^* from the spectral intensity on the
% measurement wavelength I_i and a reference wavelength I_r with the
% following operation:
%
%
% $$I^{*}_{i} = log_{10} (\frac{I_r}{I_i})$$
%
%
%
% The plots show time on the x-axis and the ABS intensity on the y axis for
% each wavelength in separate graphs. One Figure for each Channel.
figure( 'Name', 'LopapChannel1','NumberTitle','off','Position',[10 100 550 700]);
ch1s=stackedplot(LopapDataCh1);
title(LopapParameters.instrumentReferenceWavelength(1))
for ii = 1:width(LopapDataCh1)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch1s.AxesProperties(ii).YLimits = [-0.4 -0.1];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2','NumberTitle','off','Position',[600 100 550 700]);
ch2s=stackedplot(LopapDataCh2);
title(LopapParameters.instrumentReferenceWavelength(2))
for ii = 1:width(LopapDataCh2)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties(ii).YLimits = [-0.2 0.3];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot all data :: selected time range
%
if isfield(TimeRange,'dataRange')
TR = timerange(datestr(TimeRange.dataRange(1)),datestr(TimeRange.dataRange(2)));
figure( 'Name', 'LopapChannel1','NumberTitle','off','Position',[10 100 550 700]);
ch1s=stackedplot(LopapDataCh1(TR,:));
title(LopapParameters.instrumentReferenceWavelength(1))
for ii = 1:width(LopapDataCh1)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch1s.AxesProperties(ii).YLimits = [-0.4 -0.1];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
figure( 'Name', 'LopapChannel2','NumberTitle','off','Position',[600 100 550 700]);
ch2s=stackedplot(LopapDataCh2(TR,:));
title(LopapParameters.instrumentReferenceWavelength(2))
for ii = 1:width(LopapDataCh2)
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties(ii).YLimits = [-0.2 0.3];
end
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
end
%% Plot selected data with Notes
%
% Here we plot the $I^{*}_{i}$ data for one wavelength only.
if isfield(TimeRange,'dataRange')
TR = timerange(datestr(TimeRange.dataRange(1)),datestr(TimeRange.dataRange(2)));
end
ch1f = figure( 'Name', 'LopapChannel1 550','NumberTitle','off','Position',[10 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure.
if isfield(TimeRange,'dataRange')
ch1ax=plot(LopapDataCh1(TR,:).Time,LopapDataCh1(TR,:).('550.7'));
else
ch1ax=plot(LopapDataCh1.Time,LopapDataCh1.('550.7'));
end
title(LopapParameters.instrumentReferenceWavelength(1))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
%ch1s.AxesProperties.YLimits = [-0.4 -0.1];
grid ON
% Add title and axis labels
xlabel('Lab Time')
ch2f = figure( 'Name', 'LopapChannel2 550','NumberTitle','off','Position',[600 100 550 800]);
% <Adapt>: The numbers in the curled brackets give wavelength for which you'd
% like to display the data. The numbers need to be exactly match the
% labels in the previous figure
if isfield(TimeRange,'dataRange')
ch2ax=plot(LopapDataCh2(TR,:).Time,LopapDataCh2(TR,:).('550.5'));
else
ch2ax=plot(LopapDataCh2.Time,LopapDataCh2.('550.5'));
end
title(LopapParameters.instrumentReferenceWavelength(2))
% <Adapt>: The numbers in the brackets give the lower and upper limit
% of the y-axis. Adapt these to fit your data and remove the % sign.
% ch2s.AxesProperties.YLimits = [-0.2 0.3];
grid ON
% Add title and axis labels
xlabel('Lab Time')
clearvars -except Lopap* fileLocation myColorMap lopap* TimeRange* ch1f ch2f
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% LabBook :: Annotate Figure
%
%
% * 16:07 N2+water
% * 16:23 N2+R1+R2
% * 16:37 bubbles in the system
% * 16:46 N2/O2/RH(78%)+NO2+ flow tube in line + LIGHT
% * 17:05 N2/O2/RH(78%)+NO2+ flow tube in line è switch off light
% * 17:15 Close ABS channel and switch the tubes of R1+R2 to the water
% * 17:21 N2/O2/RH(78%) è open the ABS channel, bypass the reactor and close NO2 supply
% * 17:30 N2 è switch off O2 supply and the humidifier
% * 17:32 stop recording
% clc
% ax = ch1f.CurrentAxes;
% deltaT = 6;
% eTime = datetime("25-Feb-2021 16:07:00")+minutes(deltaT);
% xPos = find(LopapDataCh1.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh1(xPos,:).('550.7'),...
% '\leftarrow N2+water');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 16:23:00")+minutes(deltaT);
% xPos = find(LopapDataCh1.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh1(xPos,:).('550.7'),...
% '\leftarrow N2+R1&R2');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 16:37:00")+minutes(deltaT);
% xPos = find(LopapDataCh1.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh1(xPos,:).('550.7'),...
% '\leftarrow Bubbles!');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 16:46:00")+minutes(deltaT);
% xPos = find(LopapDataCh1.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh1(xPos,:).('550.7'),...
% {['\leftarrow N2/O2/RH(78%)+NO2+'],[' flow tube in line + LIGHT']});
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 17:05:00")+minutes(deltaT);
% xPos = find(LopapDataCh1.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh1(xPos,:).('550.7'),...
% '\leftarrow light off');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 17:15:00")+minutes(deltaT);
% xPos = find(LopapDataCh1.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh1(xPos,:).('550.7'),...
% '\leftarrow close ABS, water');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 17:21:00")+minutes(deltaT);
% xPos = find(LopapDataCh1.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh1(xPos,:).('550.7'),...
% '\leftarrow N2/O2/RH(78%) ABS open');
% set(ht,'Rotation',90);
%
% ax = ch2f.CurrentAxes;
% deltaT = 6;
% eTime = datetime("25-Feb-2021 16:07:00")+minutes(deltaT);
% xPos = find(LopapDataCh2.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh2(xPos,:).('550.5'),...
% '\leftarrow N2+water');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 16:23:00")+minutes(deltaT);
% xPos = find(LopapDataCh2.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh2(xPos,:).('550.5'),...
% '\leftarrow N2+R1&R2');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 16:37:00")+minutes(deltaT);
% xPos = find(LopapDataCh2.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh2(xPos,:).('550.5'),...
% '\leftarrow Bubbles!');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 16:46:00")+minutes(deltaT);
% xPos = find(LopapDataCh2.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh2(xPos,:).('550.5'),...
% {['\leftarrow N2/O2/RH(78%)+NO2+'],[' flow tube in line + LIGHT']});
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 17:05:00")+minutes(deltaT);
% xPos = find(LopapDataCh2.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh2(xPos,:).('550.5'),...
% '\leftarrow light off');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 17:15:00")+minutes(deltaT);
% xPos = find(LopapDataCh2.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh2(xPos,:).('550.5'),...
% '\leftarrow close ABS, water');
% set(ht,'Rotation',90);
%
% eTime = datetime("25-Feb-2021 17:21:00")+minutes(deltaT);
% xPos = find(LopapDataCh2.Time >= eTime,1);
% ht = text(ax,eTime,LopapDataCh2(xPos,:).('550.5'),...
% '\leftarrow N2/O2/RH(78%) ABS open');
% set(ht,'Rotation',90);
%% Plot selected spectra
%
% Compare to spectra to see if change is over whole wavelength or just
% where nitrite absorbs.
clc
close all
figure('Name', 'LopapSelectedSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specCompareOne),:},...
'DisplayName',['#1 ' datestr(TimeRange.specCompareOne) ]) % a time where only R1&R2 was connected to the Lopap
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specCompareOne),:},...
'DisplayName',['#2 ' datestr(TimeRange.specCompareOne)],'LineStyle','--') % a time where only R1&R2 was connected to the Lopap
plot(str2double(LopapSpecCh1.Properties.VariableNames),LopapSpecCh1{datestr(TimeRange.specCompareTwo),:},...
'DisplayName',['#1 ' datestr(TimeRange.specCompareTwo)]) % a time where nitrite in R1 was connected to the Lopap
plot(str2double(LopapSpecCh2.Properties.VariableNames),LopapSpecCh2{datestr(TimeRange.specCompareTwo),:},...
'DisplayName',['#2 ' datestr(TimeRange.specCompareTwo)],'LineStyle','--') % a time where nitrite in R1 was connected to the Lopap
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('Intensity')
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process: Absorption spectrum
%
% Here, we reprocess the data to derive an absorption spectrum. This should
% match the spectra shown for nitrite in the Lopap manual.
clc
close all
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
figure('Name', 'LopapAbsorbtionSpectra')
plot(str2double(LopapSpecCh1.Properties.VariableNames),log10(LopapSpecCh1{datestr(TimeRange.specCompareTwo),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specCompareTwo),:})-log10(LopapSpecCh1{datestr(TimeRange.specCompareOne),...
indexR(1)}./LopapSpecCh1{datestr(TimeRange.specCompareOne),:}),'DisplayName','#1')
hold on
plot(str2double(LopapSpecCh2.Properties.VariableNames),log10(LopapSpecCh2{datestr(TimeRange.specCompareTwo),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specCompareTwo),:})-log10(LopapSpecCh2{datestr(TimeRange.specCompareOne),...
indexR(2)}./LopapSpecCh2{datestr(TimeRange.specCompareOne),:}),'DisplayName','#2','LineStyle','--')
legend
box ON
grid ON
% Add title and axis labels
xlabel('wavelength [nm]')
ylabel('ABS (log I_{r}/I_{i}')
title([datestr(TimeRange.specCompareOne) ' & ' datestr(TimeRange.specCompareTwo)])
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Absorption with time
%
% Here, we reprocess the data to get the intensity with time for 1 chosen
% signal wavelength and 1 reference wavelength. The data is identical to that
% shown in the "measurement window" of the Lopap software.
%
% *The signal and the reference wavelength need to be identical to those
% used for the calibration* and are therefore taken from the calibration
% result file if present.
%
% The data is written to the _LopapTimeT_ time table. Rows hold the time,
% columns hold the processed signal for channel 1 and channel 2. Did I
% mention that it is important that the same signal wavelength and
% reference wavelength are taken for the calibration and for the
% experimental data.
clc
close all
indexS(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexS(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthS,1);
indexR(1) = find(str2double(LopapSpecCh1.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
indexR(2) = find(str2double(LopapSpecCh2.Properties.VariableNames) >= LopapParameters.wavelengthR,1);
LopapTimeT = timetable(LopapSpecCh1.Time,log10(LopapSpecCh1{:,indexR(1)}./LopapSpecCh1{:,indexS(1)}),...
log10(LopapSpecCh2{:,indexR(2)}./LopapSpecCh2{:,indexS(2)}));
LopapTimeT.Properties.VariableNames ={'Ch #1','Ch #2'};
TR = timerange(datestr(TimeRange.dataRange(1)),datestr(TimeRange.dataRange(2)));
figure('Name', 'Lopap TimeTable')
stackedplot(LopapTimeT(TR,:),{{'Ch #1','Ch #2'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Background
%
% A linear background is fitted and subtracted (may be changed to
% higher order). Background MUST be measured at the beginning
% and at the end of the experiment for at least 30 minutes, ideally every
% few hours. Background is R1&R2 with zero-air (or N2) entering the
% stripping coil.
clc
close all
LopapTimeT = fitLopapBkgrd(LopapTimeT,'BkgrdPeriod',TimeRange.backgroundFit,'TimeRangeToPlot',TimeRange.dataRange);
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Plot bkgrd corrected data
%
clc
close all
% Plot
TR = timerange(datestr(TimeRange.dataRange(1)),datestr(TimeRange.dataRange(2)));
figure('Name', 'Lopap TimeTable Bkgrd')
s=stackedplot(LopapTimeT(TR,:),{{'Ch #1 Bkgrd','Ch #2 Bkgrd'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
s.AxesProperties.YLimits = [-0.01 0.6];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange*
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
%% Process :: Calibrate
%
% The results from the calibration fits are used here to calibrate the
% signal (intensity -> HONO in ppt).
clc
close all
if exist('LopapCalibCalibration','var')
LopapTimeT.('Ch #1 Calib') = LopapCalibCalibration.calib.fitresult.calibChannel1(LopapTimeT.('Ch #1 Bkgrd'));
LopapTimeT.('Ch #2 Calib') = LopapCalibCalibration.calib.fitresult.calibChannel2(LopapTimeT.('Ch #2 Bkgrd'));
disp('Settings of the calibration:')
disp(' ')
fprintf('Flow Reagents: %.2f and %.2f ml/min \n',...
LopapParametersCalibration.LopapParameters.flowReagents(1),...
LopapParametersCalibration.LopapParameters.flowReagents(2));
fprintf('Flow Gas: %.0f ml/min \n',...
LopapParametersCalibration.LopapParameters.flowGas(1));
fprintf('Temperature: %.0f K \n',LopapParametersCalibration.LopapParameters.temperature);
fprintf('Pressure: %.0f mbar \n',LopapParametersCalibration.LopapParameters.pressure);
fprintf('Reference Wavelength: %.0f nm \n',LopapParametersCalibration.LopapParameters.wavelengthR);
fprintf('Sample Wavelength: %.0f nm \n',LopapParametersCalibration.LopapParameters.wavelengthS);
disp('Settings of this measurement')
disp('')
fprintf('Reference Wavelength: %.0f nm \n',LopapParameters.wavelengthR);
fprintf('Sample Wavelength: %.0f nm \n',LopapParameters.wavelengthS);
% Plot
TR = timerange(datestr(TimeRange.dataRange(1)),datestr(TimeRange.dataRange(2)));
figure('Name', 'Lopap TimeTable Calibrated')
s=stackedplot(LopapTimeT(TR,:),{{'Ch #1 Calib','Ch #2 Calib'}});
title(['Lopap' datestr(LopapParameters.MeasurementDate)])
s.AxesProperties.YLimits = [-100 2500];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
end
%% And Finally :: HONO
%
% * Channel 1 is from stripping coil 2. Its signal comes from NO2 and other
% nitrogen oxides that yield nitrite in the stripping solution and/or dye
% * Channel 2 is from stripping coil 1. Its signal comes from HONO plus
% the impurities that are also collected in coil 2.
%
% HONO is thus calibrated signal coil 2 (Channel 1) - that of coil 1
% (channel 2).
close all
if exist('LopapCalibCalibration','var')
LopapTimeT.('HONO') = LopapTimeT.('Ch #1 Calib')-LopapTimeT.('Ch #2 Calib');
% Plot
TR = timerange(datestr(TimeRange.dataRange(1)),datestr(TimeRange.dataRange(2)));
figure('Name', 'HONO')
s=stackedplot(LopapTimeT(TR,:),{'HONO'});
s.Title=['HONO in ppt with time ' datestr(LopapParameters.MeasurementDate)];
s.AxesProperties.YLimits = [-10 2500];
clearvars -except Lopap* fileLocation myColorMap LopapParameters TimeRange* fitresult calib
if LopapParameters.Interactive == 1
distFig
fprintf('press any key if you want to proceed\n')
pause
%%
%
close all
end
end
%% Save the data
%
save(['LopapMatlabData' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'],'LopapDataCh1','LopapDataCh2','LopapParameters','LopapSpecCh1','LopapSpecCh2')
save(['Lopap' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.mat'],'LopapParameters','LopapTimeT','TimeRange')
% export to csv
TT = LopapTimeT(:,{'HONO'});
writetimetable(TT,['Lopap' datestr(LopapParameters.MeasurementDate,'yyyymmdd_HHMM') '.csv'])
clearvars TT
%% Local functions
function [fileLocation, myColorMap] = localInitializeMatlab
%INITIALIZEMATLAB Sets Lopap specific settings
% Adds processing functions to path and defines file locations for mac, windows,
% and Matlab online (which runs on unix)
% change working directory
path = matlab.desktop.editor.getActiveFilename;
[filepath,~,~] = fileparts(path);
cd(filepath);
% Data Server Address
if ismac
fileLocation.server = '/Volumes/';
fileLocation.dataPath = 'data_2411/Lopap';
% add Lopap functions and general add-ons to path
addpath('/Users/bartels-rausch/SwitchDrive/PSI_Matlab/Add-Ons');
addpath('/Users/bartels-rausch/SwitchDrive/PSI_Matlab/Lopap/1.LOPAP_Restliche Daten');
elseif ispc
fileLocation.server = '\\fs00\'; %UNC name for file server that hosts Napp shared drive
fileLocation.dataPath = 'data_2411\Lopap';
% add Lopap functions and general add-ons to path
addpath('C:\Users\iezzi_l\switchdrive\Add-Ons');
addpath('C:\Users\bletry_c\switchdrive\Add-Ons');
addpath('C:\Users\vieill_l\switchdrive\Add-Ons');
addpath('C:\Users\iezzi_l\switchdrive\Lopap\1.LOPAP_Restliche Daten');
addpath('C:\Users\bletry_c\switchdrive\Lopap\1.LOPAP_Restliche Daten');
addpath('C:\Users\vieill_l\switchdrive\Lopap\1.LOPAP_Restliche Daten');
% addpath('\\fs02\5601_surface_chemistry\People\Marc\Switchdrive\PSI_Matlab\Add-Ons');
% addpath('\\fs02\5601_surface_chemistry\People\Marc\Switchdrive\PSI_Matlab\Lopap\1.LOPAP_Restliche Daten');
elseif isunix
fileLocation.server = '/';
fileLocation.dataPath = 'Lopap/';
end
% Color maps
myColorMap.red = colormap(autumn(7));
myColorMap.blue = colormap(winter(7));
close all
end