SPS-BCT lifetime calculation for MD4/2004

<< Graphics`Colors`
<< Graphics`MultipleListPlot`
<< Graphics`Legend`
<< Statistics`LinearRegression`

AllColors;

ROSALI Package version 0.1

Data Formats

Whatever the syntax, the structure  of all BDI datasets is always the same with a header and an arbitrary number of data tables; each table includes a table header followed by the data proper.
      There are three types of tables (like in Excel): Line (the abscissa is implicitly the natural numbers), XYScatter (the abscissa is the first column of the table) and Surface for 2D information, typically images.

The standard BDI data format on external files  is the Rosali XML format (in 2004). However, most SPS instruments produce adhoc formats and require a specific read procedure. The XML and BCT formats are handeled in Rosali v0.1.
In Mathematica, a Rosali object is a replacement list made of a header and a table list. Each entry in the table list is a replacement rule for a table header and the table data. The data structure and the available procedures for access to the object follow as much as possible those of the MFS package by J. Jowett.

Available Functionalities

Read BDI files

- rOsRead[infile]: read in an XML Rosali file and produces a Rosali object; XML comments are presently not supported
- rOsReadBCT[infile]: read in a SPS BCT file and produces a Rosali object; returns {-1} if the file does not exist

Summary of information in an object

- rOsHeader[dataset]: returns the header of the Rosali object in form of a list
- rOsTableDirectory[dataset]: returns a list containing the numerical index and the dataset description for each table

Access to the header information

- rOsHeaderKeyNames[dataset]:returns the list of key names in the header
- rOsHeaderKeyValue[dataset,keyname]:returns the value associated to a header key name
- rOsHeaderKeyUnit[dataset,keyname]:returns the unit associated to a header key name

Access to the table information

- rOsTableIndex[dataset,tableID]: for clarity, the tables be either pointed at by their numerical index or their dataset description (string); this procedures returns the numerical index whether tableID is the numerical index or the dataset description.
- rOsTableType[dataset,tableID]: returns the Rosali table type (Line, XYScatter, Surface)
- rOsTableHeader[dataset,tableID]: returns the table header in form of a list
- rOsTableHeaderKeyNames[dataset,tableID]: returns a list of the key names in the table header
- rOsTableHeaderKeyValue[dataset,tableID,keyname]: returns the value associated to a table key name
- rOsTableHeaderKeyUnit[dataset,tableID,keyname]: returns the unit associated to a table key name
- rOsTableData[dataset,tableID]: returns the data contained in the table including the column names and units if relevant
- rOsTableDataColumnNames[dataset,tableID]: returns the list of column names
- rOsTableDataColumnValue[dataset,tableID,columnname]: returns the data of a given column

Plotting

rOsPlotPalette: selection of colors and styles for default plotting
rOsTableDataColumnTitles[dataset,tableID]: prepares axes lables for scatter and line plots (column name + unit)
rOsTableSurfaceAxesLabels[dataset,tableID]: prepares axes labels for Surface plots
rOsPlot[dataset,tableID]: general plotting procedure calling the three following procedures for each table type
rOsXYScatterPlot[dataset,tableID]: 
rOsLinePlot[dataset,tableID]
rOsSurfacePlot[dataset,tableID]

Data processing

rOsDataList[ dataset,tableID, col1,col2]: returns a time series using col1 as the independent variable and col2 for the function
rOsExponentialFit[datalist,startT, stopT]: fits an exponential decay to a time series between startT and stopT

Warning

Exception situations are not tested but handled directly by Mathematica. This should be gradually improved to provide clear diagnostic messages.

Procedures

rOsRead[infile_]:=\({
ref = OpenRead[infile]; 

        inputBuffer = 
          ReadList[ref, Record, 
            RecordSeparators -> {{"<header>
", "<table>
"}, {"
</header>", "</table>
"}}]; 
        headerbuffer = 
          ReadList[StringToStream[inputBuffer[\([1]\)]], Word, 
            WordSeparators -> "	", RecordLists -> True]; 
        headercontent = 
          Table[{\(headerbuffer[\([i]\)]\)[\([1]\)], 
              If[StringMatchQ[\(headerbuffer[\([i]\)]\)[\([2]\)], ""*""], StringDrop[
                  StringDrop[\(headerbuffer[\([i]\)]\)[\([2]\)], 1], \(-1\)], 
                ToExpression[\(headerbuffer[\([i]\)]\)[\([2]\)]]], 
              StringDrop[
                StringDrop[\(headerbuffer[\([i]\)]\)[\([3]\)], 
                  1], \(-1\)]}, {i, 1, 
              Length[headerbuffer]}]; 

tablecontent = {}; 
For[itab = 1, 
          itab ≤ 
            Length[inputBuffer] - 1, \(itab++\), 
buffer1 = 
            Read[StringToStream[inputBuffer[\([itab + 1]\)]], Record, 
              RecordSeparators -> {{"<tableheader>
"}, {"
</tableheader>"}}]; 
          tableheaderbuffer = 
            ReadList[StringToStream[buffer1], Word, 
              WordSeparators -> "	", RecordLists -> True]; 
          tableheadercontent = 
            Table[{\(tableheaderbuffer[\([i]\)]\)[\([1]\)], 
                If[StringMatchQ[\(tableheaderbuffer[\([i]\)]\)[\([2]\)], ""*""], StringDrop[
                    StringDrop[\(tableheaderbuffer[\([i]\)]\)[\([2]\)], 
                      1], \(-1\)], 
                  ToExpression[\(tableheaderbuffer[\([i]\)]\)[\([2]\)]]], 
                StringDrop[
                  StringDrop[\(tableheaderbuffer[\([i]\)]\)[\([3]\)], 
                    1], \(-1\)]}, {i, 1, Length[tableheaderbuffer]}]; 
          buffer2 = 
            Read[StringToStream[inputBuffer[\([itab + 1]\)]], Record, 
              RecordSeparators -> {{"<tabledata>
"}, {"
</tabledata>"}}]; 
          tabledatabuffer = 
            ReadList[StringToStream[buffer2], Word, 
              WordSeparators -> "	", RecordLists -> True]; 
          tabledatacontent = 
            Table[If[
                StringMatchQ[tabledatabuffer[\([j, i]\)], ""*""] || 
                  StringMatchQ[tabledatabuffer[\([j, i]\)], "[*]"], 
                StringDrop[
                  StringDrop[tabledatabuffer[\([j, i]\)], 1], \(-1\)], 
                ToExpression[tabledatabuffer[\([j, i]\)]]], {j, 1, 
                Length[tabledatabuffer]}, {i, 1, 
                Length[tabledatabuffer[\([j]\)]]}]; 
          tablecontent = 
            Append[tablecontent, {tableheader -> tableheadercontent, 
                tabledata -> tabledatacontent}]
]; 
        Close[infile]; header -> headercontent, 
        table -> tablecontent}\)
rOsReadBCT[infile_]:=\((ref = OpenRead[infile]; 
        If[ref == $Failed, 
          Return[{\(-1\)}]]; 
headercontent = {}; 
        Read[ref, String]; buf = Read[ref, {Word, Word, Word}]; 
        headercontent = 
          Append[headercontent, {"InstrName", 
              buf[\([1]\)] <> buf[\([2]\)], ""}]; 
        headercontent = 
          Append[headercontent, {"DevName", buf[\([3]\)], ""}]; 
        ReadList[ref, String, 2]; buf = Read[ref, {Word, Word, Word}]; 
        headercontent = 
          Append[headercontent, {"CycName", buf[\([3]\)], ""}]; 
        headercontent = 
          Append[headercontent, {"AcqDesc", "Beam current evolution in cycle", ""}]; 
        headercontent = 
          Append[headercontent, {"CycDate", "Missing", ""}]; 
        headercontent = 
          Append[headercontent, {"CycTime", "Missing", ""}]; 
        headercontent = 
          Append[headercontent, {"CycTS", "Missing", "s"}]; 
        headercontent = 
          Append[headercontent, {"AcqDate", "Missing", ""}]; 
        headercontent = 
          Append[headercontent, {"AcqTime", "Missing", ""}]; 
        headercontent = 
          Append[headercontent, {"AcqTS", "Missing", "s"}]; 
        headercontent = 
          Append[headercontent, {"MsInCyc", "Missing", "ms"}]; 
buf = Read[ref, {Word, Word, Word, Word, Word}]; 
        headercontent = 
          Append[headercontent, {"SuperCycNb", buf[\([5]\)], ""}]; 
        buf2 = Read[ref, {Word, Word, Word, Word, Word, Word, Word, Word}]; 
        buf3 = Read[ref, {Word, Word, Word}]; 
        buf4 = Read[ref, {Word, Word, Word, Word}]; 
        buf = Read[ref, {Word, Word, Word, Word, Word}]; 
        headercontent = 
          Append[headercontent, {SCLength, buf[\([4]\)], buf[\([5]\)]}]; 
        tablecontent = {}; tableheadercontent = {}; tabledatacontent = {}; 
        tableheadercontent = 
          Append[tableheadercontent, {"DataSetType", \ "XYScatter", ""}]; tableheadercontent = 
          Append[tableheadercontent, {"DataSetDesc", "Beam current evolution in cycle"\ , ""}]; 
        tableheadercontent = 
          Append[tableheadercontent, {"NbOfColumns", 2\ , ""}]; 
        tableheadercontent = 
          Append[tableheadercontent, {"NbOfRaws", 
              ToExpression[buf3[\([3]\)]] + 2\ , ""}]; 
        tableheadercontent = 
          Append[tableheadercontent, {"DataFormat", "Double"\ , ""}]; 
        tableheadercontent = 
          Append[tableheadercontent, {"NbChannels", 1\ , ""}]; 
        tableheadercontent = 
          Append[tableheadercontent, {"NbAcqPoints", 
              buf3[\([3]\)]\ , ""}]; 
buf2 = 
          Read[ref, {Word, Word, Word, Word}]; 
        tablecontent = 
          Append[tablecontent, 
            tableheader -> tableheadercontent]; 
buf = 
          Read[ref, {Word, Word, Word, \ Word}]; 
buf3 = 
          Transpose[ReadList[ref, {Number, Number, Number}]]; 
        tabledatacontent = Transpose[{buf3[\([2]\)], buf3[\([3]\)]}]; 
        tabledatacontent = 
          Prepend[tabledatacontent, {StringDrop[buf[\([2]\)], 5], 
              buf2[\([4]\)] <> " " <> buf[\([4]\)]}]; 
        tabledatacontent = 
          Prepend[tabledatacontent, {StringTake[buf[\([2]\)], 4], 
              StringTake[buf[\([3]\)], 9]}]; 
        tablecontent = 
          Append[tablecontent, 
            tabledata -> 
              tabledatacontent]; 
{header -> 
            headercontent, table -> {tablecontent}})\)
rOsHeader[dataset_]:=header /. dataset
rOsHeaderKeyNames[dataset_]:=Transpose[
          header /. dataset][[1]]
simplifyUnitList[list_] := 
      If[Length[list] == 1, list[[1]], list]
rOsHeaderKeyValue[dataset_, keyname_]:=simplifyUnitList[
        Part[Transpose[header /. dataset][[2]], 
          Flatten[Table[
              Position[Transpose[header /. dataset][[1]], 
                ToString[Flatten[{keyname}][[ijk]]]], {ijk, 1, 
                Length[Flatten[{keyname}]]}]]]]
rOsHeaderKeyUnit[dataset_, keyname_]:=simplifyUnitList[
        Part[Transpose[header /. dataset][[3]], 
          Flatten[Table[
              Position[Transpose[header /. dataset][[1]], 
                ToString[Flatten[{keyname}][[ijk]]]], {ijk, 1, 
                Length[Flatten[{keyname}]]}]]]]
rOsTableDirectory[dataset_]:=Table[{i, rOsTableHeaderKeyValue[dataset, i, DataSetDesc]}, {i, 
          1, Length[table /. dataset]}]
rOsTableIndex[dataset_, tablename_]:=If[NumberQ[tablename], 
        tablename, Position[Transpose[
                  rOsTableDirectory[dataset]][[2]], 
              tablename][[1]][[1]]]
rOsTableType[dataset_, tableID_]:=rOsTableHeaderKeyValue[dataset, tableID, 
        DataSetType]
rOsTableHeader[dataset_, tableID_]:=tableheader /. (table /. dataset)[[rOsTableIndex[
              dataset, tableID]]]
rOsTableHeaderKeyNames[dataset_, tableID_]:=Transpose[
          tableheader /. (table /. dataset)[[rOsTableIndex[dataset, 
                  tableID]]]][[1]]
rOsTableHeaderKeyValue[dataset_, tableID_, 
        keyname_]:=simplifyUnitList[
        Part[Transpose[
              tableheader /. (table /. dataset)[[rOsTableIndex[dataset, 
                      tableID]]]][[2]], 
          Flatten[Table[
              Position[Transpose[
                    tableheader /. (table /. dataset)[[rOsTableIndex[
                            dataset, tableID]]]][[1]], 
                ToString[Flatten[{keyname}][[ijk]]]], {ijk, 1, 
                Length[Flatten[{keyname}]]}]]]]
rOsTableHeaderKeyUnit[dataset_, tableID_, 
        keyname_]:=simplifyUnitList[
        Part[Transpose[
              tableheader /. (table /. dataset)[[rOsTableIndex[dataset, 
                      tableID]]]][[3]], 
          Flatten[Table[
              Position[Transpose[
                    tableheader /. (table /. dataset)[[rOsTableIndex[
                            dataset, tableID]]]][[1]], 
                ToString[Flatten[{keyname}][[ijk]]]], {ijk, 1, 
                Length[Flatten[{keyname}]]}]]]]
rOsTableData[dataset_, tableID_]:=tabledata /. (table /. dataset)[[rOsTableIndex[dataset, 
              tableID]]]
rOsTableDataColumnNames[dataset_, tableID_]:=
If[
        rOsTableType[dataset, 
              rOsTableIndex[dataset, tableID]] == "XYScatter" || 
          rOsTableType[dataset, 
              rOsTableIndex[dataset, 
                tableID]] == "Line", (tabledata /. (table /. 
                    dataset)[[rOsTableIndex[dataset, 
                    tableID]]])[[1]], "{1.." <> 
          ToString[
            Length[(tabledata /. (table /. dataset)[[rOsTableIndex[
                          dataset, 
                          tableID]]])[[1]]]] <> "}"]
rOsTableDataColumnValue[dataset_, tableID_, 
        columnname_]:=\(If[rOsTableType[dataset, tableID] != "Surface", 
        Table[\(Transpose[
              Drop[rOsTableData[dataset, tableID], 
                2]]\)[\([\(\(Position[
                    rOsTableDataColumnNames[dataset, 
                      tableID], \(Flatten[{columnname}]\)[\([i]\)]]\)[\([1]\)]\)[\([1]\)]]\)], {i, 1, 
            Length[Flatten[{columnname}]]}], 
        Part[Transpose[rOsTableData[dataset, tableID]], 
          columnname]]\)
rOsPlotPalette={Red, Blue, Green, Yellow, Magenta, Orange, 
          LightBlue, {Violet, Dashing[{ .01}]}};
rOsTableDataColumnTitles[dataset_, tableID_]:=
If[
        rOsTableType[dataset, 
              rOsTableIndex[dataset, tableID]] == "XYScatter" || 
          rOsTableType[dataset, 
              rOsTableIndex[dataset, tableID]] == "Line", 
        buffer = 
          Transpose[{(tabledata /. (table /. dataset)[[rOsTableIndex[
                          dataset, 
                          tableID]]])[[1]], (tabledata /. (table /. 
                          dataset)[[rOsTableIndex[dataset, 
                          tableID]]])[[2]]}]; 
        Table[ToString[buffer[[i, 1]]] <> " [" <> 
            ToString[buffer[[i, 2]]] <> "]", {i, 1, 
            Length[buffer]}], "{1.." <> 
          ToString[
            Length[(tabledata /. (table /. dataset)[[rOsTableIndex[
                          dataset, 
                          tableID]]])[[1]]]] <> "}"]
rOsPlot[dataset_, tableID_]:=Which[
rOsTableType[dataset, 
            tableID] == "XYScatter", 
        rOsXYScatterPlot[dataset, 
          tableID], 
rOsTableType[dataset, 
            tableID] == "Line", 
        rOsLinePlot[dataset, 
          tableID], 
rOsTableType[dataset, 
            tableID] == "Surface", 
        rOsSurfacePlot[dataset, tableID]
]
rOsXYScatterPlot[dataset_, tableID_]:=\(MultipleListPlot[
        Table[
Transpose[{\(rOsTableDataColumnValue[
                  dataset, 
                  tableID, \(rOsTableDataColumnNames[dataset, 
                      tableID]\)[\([1]\)]]\)[\([1]\)], \(rOsTableDataColumnValue[dataset, tableID, 
                  Drop[rOsTableDataColumnNames[dataset, tableID], 
                    i - 1]]\)[\([1]\)]}], 
{i, 2, 
            Length[rOsTableDataColumnNames[dataset, tableID]]}], 
        PlotJoined -> True, PlotStyle -> rOsPlotPalette, 
        Frame -> True, FrameStyle -> Red, 
        GridLines -> Automatic, 
        PlotLegend -> Drop[rOsTableDataColumnTitles[dataset, tableID], 1], 
        LegendTextSpace -> 4, LegendShadow -> Automatic, 
        PlotLabel -> rOsTableHeaderKeyValue[dataset, tableID, DataSetDesc], 
        FrameLabel -> {\(rOsTableDataColumnTitles[dataset, 
                tableID]\)[\([1]\)], ""}]\)
rOsLinePlot[dataset_, tableID_]:=MultipleListPlot[
        Drop[Transpose[Drop[rOsTableData[dataset, tableID], 2]], 1], 
        PlotJoined -> True, PlotStyle -> rOsPlotPalette, 
        Frame -> True, FrameStyle -> Red, 
        GridLines -> Automatic, 
        PlotLegend -> Drop[rOsTableDataColumnTitles[dataset, tableID], 1], 
        LegendTextSpace -> 4, LegendShadow -> Automatic, 
        PlotLabel -> rOsTableHeaderKeyValue[dataset, tableID, DataSetDesc], 
        FrameLabel -> {"Index", ""}, 
        FrameTicks -> {Transpose[{Transpose[
                      Drop[rOsTableData[dataset, tableID], 2]][[2]] + 
                  1, Transpose[
                    Drop[rOsTableData[dataset, tableID], 2]][[1]]}], 
            Automatic}]
rOsTableSurfaceAxesLabels[dataset_, tableID_]:={ToString[rOsTableHeaderKeyValue[testdataset, 
                tableID, {RowTag, ColTag}][[1]]] <> " [" <> 
          ToString[rOsTableHeaderKeyUnit[testdataset, 
                tableID, {RowTag, ColTag}][[1]]] <> "]", 
        ToString[rOsTableHeaderKeyValue[testdataset, 
                tableID, {RowTag, ColTag}][[2]]] <> " [" <> 
          ToString[rOsTableHeaderKeyUnit[testdataset, 
                tableID, {RowTag, ColTag}][[2]]] <> "]"}
rOsSurfacePlot[dataset_, tableID_]:=ShowLegend[
        ListDensityPlot[rOsTableData[dataset, tableID], 
          FrameStyle -> Red, 
          PlotLabel -> rOsTableHeaderKeyValue[dataset, tableID, DataSetDesc], 
          FrameLabel -> rOsTableSurfaceAxesLabels[dataset, tableID], 
          DisplayFunction -> Identity], {GrayLevel[1 - #] &, 10, 
          ToString[Max[rOsTableData[dataset, tableID]]], 
          ToString[Min[rOsTableData[dataset, tableID]]], 
          LegendPosition -> {1.1, - .4}}]
rOsDataList[ dataset_, tableID_,  
        col1_, col2_]:=Transpose[{Flatten[
            rOsTableDataColumnValue[dataset, tableID, col1]], 
          Flatten[rOsTableDataColumnValue[dataset, tableID, 
              col2]]}]
rOsExponentialFit[ data_, startT_,  
        stopT_]:=\({x = \(Transpose[data]\)[\([1]\)]; 
        y = Log[\(Transpose[data]\)[\([2]\)]]; 

        debut = \(\(Position[x, startT]\)[\([1]\)]\)[\([1]\)]; 

        fin = \(\(Position[x, stopT]\)[\([1]\)]\)[\([1]\)]; 

        regress = 
          Regress[Take[Transpose[{x, y}], {debut, fin}], {1, t}, t, 
            RegressionReport -> {BestFit, ParameterTable}]; 
        a0 = \((BestFit /. regress)\) /. t -> 0; 
        k = \(-Coefficient[\((BestFit /. regress)\), t]\); τ = 1/k; 
        da0 = \(\(\((ParameterTable /. 
                    regress)\)[\([1]\)]\)[\([1]\)]\)[\([2]\)]; 

        dk = \(\(\((ParameterTable /. 
                    regress)\)[\([1]\)]\)[\([2]\)]\)[\([2]\)]; 

        equation -> Exp[a0]\ Exp[\(-t\)/τ], 
        coefficients -> {Exp[a0], τ}, 
        rmscoefficients -> {\((Exp[a0 + da0] - Exp[a0 - da0])\)/
              2, \((1/\((k - dk)\) - 1/\((k + dk)\))\)/2}}\)

Needed fixes or upgrades for next version

- Check the validity of the start and end times in rOsExponentialFit;

Specific procedures to the analysis of the BBLR MD 2004-3

\(showDecay[] := {p1 = 
          ListPlot[datalist, PlotStyle -> Red, PlotJoined -> True, 
            Frame -> True, GridLines -> Automatic, 
            FrameLabel -> rOsTableDataColumnTitles[bctdata, 1], 
            PlotLabel -> "current decay at cycle " <> 
                rOsHeaderKeyValue[bctdata, SuperCycNb] <> "; lifetime=" <> 
                ToString[\((coefficients /. fit)\)[\([2]\)]/
                    1000] <> "±" <> 
                ToString[\((rmscoefficients /. fit)\)[\([2]\)]/
                    1000] <> " s", DisplayFunction -> Identity]; 
        p2 = Plot[equation /. fit, {t, startfit, stopfit}, 
            PlotStyle -> Blue, DisplayFunction -> Identity]; 
        decaycalcperiod = {{startfit, 
              0}, {startfit, \(equation /. fit\) /. 
                t -> startfit}, {startfit, 0}, {stopfit, 
              0}, {stopfit, \(equation /. fit\) /. 
                t -> stopfit}, {stopfit, 0}}; 
        p3 = ListPlot[decaycalcperiod, PlotStyle -> Blue, 
            PlotJoined -> True, 
            DisplayFunction -> 
              Identity]; 
wirerampfordisplay = 
          Transpose[{\(Transpose[wireramp]\)[\([1]\)], 
              0.7  Max[\(Transpose[datalist]\)[\([2]\)]]\ \(Transpose[
                    wireramp]\)[\([2]\)]}]; 
p4 = 
          ListPlot[wirerampfordisplay, PlotStyle -> Green, 
            PlotJoined -> True, DisplayFunction -> Identity]; 
        pp1 = Show[p1, p2, p3, p4]; 
        pp2 = MultipleListPlot[
            Table[{{\(lifetimelog[\([i]\)]\)[\([1]\)], \(lifetimelog[\([i]\)]\)[\([2]\)]}, ErrorBar[\(lifetimelog[\([i]\)]\)[\([3]\)]]}, {i, 1, 
                Length[lifetimelog]}], PlotStyle -> Red, 
            PlotJoined -> True, GridLines -> Automatic, 
            Frame -> True, 
            PlotLabel -> "Lifetime versus cycle number", \ 
            FrameLabel -> {"cycle number", "lifetime in s"}, 
            PlotRange -> {Automatic, {0, 7200}}, 
            DisplayFunction -> Identity]; 
pp3 = 
          ListPlot[
            Transpose[{\(Transpose[
                    logDif]\)[\([1]\)], \(Transpose[\(Transpose[
                        logDif]\)[\([2]\)]]\)[\([1]\)]}], 
            PlotJoined -> True, PlotStyle -> Red, 
            Frame -> True, GridLines -> Automatic, 
            PlotLabel -> "Injected beam current", 
            FrameLabel -> {"cycle number", \((rOsTableDataColumnTitles[bctdata, 1])\)[\([2]\)]
}, 
            DisplayFunction -> Identity]; 
        pp4 = MultipleListPlot[
            Sequence @@ 
              Table[Transpose[{\(Transpose[
                        logDif]\)[\([1]\)], \(Transpose[\(Transpose[
                            logDif]\)[\([2]\)]]\)[\([i]\)]}], {i, 2, 
                  Length[\(logDif[\([1]\)]\)[\([2]\)]]}], 
            PlotJoined -> True, PlotStyle -> rOsPlotPalette, 
            Frame -> True, GridLines -> Automatic, 
            PlotLabel -> "Fractional intensity losses:
R:" <> 
                ToString[diftiming[\([2]\)]] <> " B:" <> 
                ToString[diftiming[\([3]\)]] <> " G:" <> 
                ToString[diftiming[\([4]\)]], 
            FrameLabel -> {"cycle number", \((rOsTableDataColumnTitles[bctdata, 1])\)[\([2]\)]}, 
            DisplayFunction -> Identity]; 
        Show[GraphicsArray[{{pp1, pp2}, {pp3, pp4}}], 
          DisplayFunction -> $DisplayFunction]}\)

rOsListDifferences[datalist_, intervals_] := 
      Table[xx = Position[Transpose[
                    datalist][[1]], intervals[[i]][[1]]][[1]][[1]]; 
        yy = Position[Transpose[
                    datalist][[1]], intervals[[i]][[2]]][[1]][[1]]; (Transpose[
                  datalist][[2]])[[xx]] - (Transpose[
                  datalist][[2]])[[yy]], {i, 1, 
          Length[intervals]}]

Analysis program

Access to MD files require a bridge to the local store:
Run on a PC: \\hpdepot\opdata3, with opera and sps2lhc
data is in data\SPS_MD\BBLR2004-2608 and BBLR2004-2608-2

SetDirectory["\\\\hpdepot\\opdata3\\SPS_MD\\BBLR04-02Sep04"]

wireramp = {{0, 0}, {0, 1}, {4300, 1}, {4600, 0}};

noiselevel = 20;

taulogNobb = lifetimelog; deltalogNobb = logDif;

taulogExc = lifetimelog; deltalogExc = logDif;

taulogComp = lifetimelog; deltalogComp = logDif;

lifetimelog = {}; logDif = {}; 
    diftiming = {{0, 4650}, {0, 100}, {100, 3450}, {3450, 3670}};

cyclenumber = 178452;

bctfilename = "sps_bct_4_MTG#" <> 
          ToString[cyclenumber] <> "_xrx34.dat";
bctdata = rOsReadBCT[bctfilename]; 
If[bctdata == {-1}, 
      Print["TRY AGAIN OR CHANGE FILE NAME"], (datalist = 
          rOsDataList[ bctdata, 1,  "time", "intensity"]; 
        maxlevel = Max[Transpose[datalist][[2]]];)]
\(\(If[
        bctdata != {\(-1\)}\  && \ 
          maxlevel > noiselevel, 
startfit = 100; \ 
        stopfit = 3450; 
        fit = rOsExponentialFit[datalist, startfit, stopfit]; 

        lifetimelog = 
          Append[lifetimelog, {ToExpression[
                rOsHeaderKeyValue[bctdata, 
                  SuperCycNb]], \((coefficients /. fit)\)[\([2]\)]/1000, 
              1.  \((rmscoefficients /. fit)\)[\([2]\)]/1000}]; 
        logDif = 
          Append[logDif, {ToExpression[
                rOsHeaderKeyValue[bctdata, SuperCycNb]], 
              rOsListDifferences[datalist, diftiming]}]; 

        showDecay[]; \(cyclenumber++\);];\)\)
If[bctdata != {-1}  &&  maxlevel <= noiselevel, showDecay[]; 
      Print["NO BEAM IN CYCLE ", 
        cyclenumber]; cyclenumber++]

logDif

Show[pp1, 
      Graphics[{Blue, 
          Text["Lifetime calculation\nwindow", {1800, 4800}]}], 
      DisplayFunction -> $DisplayFunction]

plot1 = 
      Show[pp2, Graphics[{Red, Text["compensation", {178392, 6200}]}], 
        Graphics[{Red, Text["excitation", {178415, 3800}]}], 
        Graphics[{Red, Text["no beam-beam", {178450, 5200}]}], 
        DisplayFunction -> $DisplayFunction]

Export["G:\\Users\\j\\jpk\\MyDocs\\Math\\BBLR2004-3\\fig1.eps",pp2,"EPS"]

\(pp3 = 
      ListPlot[Transpose[{\(Transpose[
                logDif]\)[\([1]\)], \(Transpose[\(Transpose[
                    logDif]\)[\([2]\)]]\)[\([1]\)]}], 
        PlotJoined -> True, PlotStyle -> Red, Frame -> True, 
        GridLines -> Automatic, 
        PlotLabel -> "Injected beam current", 
        FrameLabel -> {"cycle number", \((rOsTableDataColumnTitles[
                  bctdata, 1])\)[\([2]\)]
}]\)

Drop[diftiming, 1]

pp4 = 
      MultipleListPlot[
        Sequence @@ 
          Table[Transpose[{Transpose[
                    logDif][[1]], Transpose[Transpose[
                        logDif][[2]]][[i]]}], {i, 2, 
              Length[logDif[[1]][[2]]]}], 
        PlotJoined -> True, PlotStyle -> rOsPlotPalette, 
        Frame -> True, GridLines -> Automatic, 
        PlotLabel -> "Fractional intensity losses:\nR:" <> 
            ToString[diftiming[[2]]] <> " B:" <> 
            ToString[diftiming[[3]]] <> " G:" <> 
            ToString[diftiming[[4]]], 
        FrameLabel -> {"cycle number", (rOsTableDataColumnTitles[
                  bctdata, 1])[[2]]}, 
        PlotLegend -> Drop[diftiming, 1], LegendTextSpace -> 4, 
        LegendShadow -> Automatic]

pp4 = 
        MultipleListPlot[
          Sequence @@ 
            Table[Transpose[{Transpose[
                      logDif][[1]], Transpose[Transpose[
                          logDif][[2]]][[i]]}], {i, 1, 
                Length[logDif[[1]][[2]]]}], 
          PlotJoined -> True, PlotStyle -> rOsPlotPalette, 
          Frame -> True, GridLines -> Automatic, 
          PlotLabel -> "Fractional intensity losses", 
          FrameLabel -> {"cycle number", "intensity loss in "}, 
          PlotLegend -> diftiming, LegendTextSpace -> 4, 
          LegendShadow -> Automatic, 
          PlotRange -> {0, 1000}];

buffer = 
        Transpose[{Transpose[logDif][[1]], 
            Transpose[{Transpose[Transpose[logDif][[2]]][[1]]/
                  10. , Transpose[Transpose[
                        logDif][[2]]][[2]], Transpose[Transpose[
                        logDif][[2]]][[3]], Transpose[Transpose[logDif][[2]]][[4]]}]}];

pp3 = 
        MultipleListPlot[
          Sequence @@ 
            Table[Transpose[{Transpose[
                      buffer][[1]], Transpose[Transpose[
                          buffer][[2]]][[i]]}], {i, 1, 
                Length[buffer[[1]][[2]]]}], 
          PlotJoined -> True, PlotStyle -> rOsPlotPalette, 
          Frame -> True, GridLines -> Automatic, 
          PlotLabel -> "Fractional intensity losses", 
          FrameLabel -> {"cycle number", "intensity loss in "}, 
          PlotLegend -> diftiming, LegendTextSpace -> 4, 
          LegendShadow -> Automatic];

plot2 = 
      Show[pp3, Graphics[{Red, Text["compensation", {178392, 400}]}], 
        Graphics[{Red, Text["excitation", {178415, 450}]}], 
        Graphics[{Red, Text["no beam-beam", {178450, 460}]}], 
        DisplayFunction -> $DisplayFunction]

plot2 = 
      Show[pp3, Graphics[{Red, Text["compensation", {-0.6, 0}]}], 
        Graphics[{Red, Text["excitation", {-0.09, 0.1}]}], 
        Graphics[{Red, Text["no beam-beam", {0.75, 0.12}]}], 
        DisplayFunction -> $DisplayFunction]

taulogNobb = lifetimelog; deltalogNobb = logDif;

taulogExc = lifetimelog; deltalogExc = logDif;

taulogComp = lifetimelog; deltalogComp = logDif;

{Mean[Transpose[taulogComp][[2]]], 
      StandardDeviation[Transpose[taulogComp][[2]]]}

i = 1;

i = 4; 
    N[{Mean[Transpose[Transpose[deltalogNobb][[2]]][[i]]], 
        StandardDeviation[Transpose[Transpose[
                  deltalogNobb][[2]]][[i]]]}]