<< NumericalMath`ListIntegrate`
<< Statistics`DataSmoothing`
<< NumericalMath`PolynomialFit`
hide = DisplayFunction -> Identity
SetDirectory["/afs/cern.ch/user/a/afabich/public/MERIT"]
<< "prolog.m"

    
SetDirectory["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/\
MERIT_libera"]
Files = FileNames["2007*result.txt"];

\(\(\( (*\ 
      Read\ run\ list\ of\ Marcus\ and\ modify\ *) \)\(
\)\(\(ReadList["/afs/cern.ch/user/a/afabich/public/MERIT/RunList.data", String]\)[\([1]\)]

    \(test = \(RunList = 
Transpose[
Map[
              StringSplit, 
Drop[
ReadList["/afs/cern.ch/user/a/afabich/public/MERIT/RunList.data", String]
, 
                1]
]
]\);\)

    Harmonic[x_] := 
      If[x == "1.31e-007" || 
          x == "1.31288e-007", "16", 
        If[x == "2.62576e-007" || 
            x == "2.62e-007", "8", 
          If[x == "5.25e-007", "5.25", 
            If[x == "5.27e-007", "6", 
              If[x == "7.95e-007", "8", 
                x]]]]]

    \(RunList = 
        ReplacePart[RunList, Map[Harmonic, RunList[\([2]\)]], 
          2];\)

    \(RunList = 
        Select[
ToExpression[
            Transpose[RunList]
], #[\([1]\)] ≠ 
              0 &];\)

    \(RunListOff = 
        Select[
ToExpression[
            RunList]
, #[\([2]\)] ≠ 
                16 && #[\([2]\)] ≠ 8 &];\)\ 

    \(RunList8 = 
        Select[
ToExpression[
            RunList]
, #[\([2]\)] == 
              8 &];\)\ 

    \(RunList = 
        Select[
ToExpression[
            RunList]
, #[\([2]\)] == 
              16 &];\)\ 

    RunList[\([1]\)]
 (*\ 
      Check\ availability\ of\ input\ files\ *) 

    Length[RunList]\)\)\)

\(Do[
If[
        1 ≠ 
FileNames[
"/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/MERIT_libera/" <> 
ToString[RunList[\([k, 4]\)]] <> "-" <> 
              ToString[RunList[\([k, 5]\)]] <> "-" <> 
              ToString[RunList[\([k, 6]\)]] <> "_" <> 
              ToString[RunList[\([k, 7]\)]] <> "-" <> 
              ToString[RunList[\([k, 8]\)]] <> "-" <> 
              ToString[
                RunList[\([k, 
                    9]\)]] <> "_MEAS_result.txt"]
, 
        Print["missing"]
]
, {k, 
        Length[RunList]}]\)
Print["if not <missing> printed, every file is available"]

SetDirectory["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_\
MERIT/MERIT_libera"]

      
data162 = ReadList[Files[[162]], Number];    
FullData162 = 
      Transpose[{Table[m, {m, Length[data162]}], 
          data162}];
Pulse = Take[data162, {873, 3000}];
func = 
        Interpolation[Pulse, InterpolationOrder -> 5];
func2 = 
        Interpolation[
          Transpose[{Table[m, {m, 60, 70}], Take[Pulse, {60, 70}]}], 
          InterpolationOrder -> 5];
PulseMaxI = 
        First[First[
            Position[Pulse, Max[Take[Pulse, {60, 70}]]]]];
PulseMax = 
        Table[{x, func2[x]}, {x, 60, 70, 0.05}];

    
PulseH = PulseMax[[First[
              First[Position[PulseMax, 
                  Max[PulseMax]]]]]][[1]]
PulseMax = func2[PulseH]

    
DisplayTogether[
ListPlot[Pulse, 
        PlotRange -> {{1, 100}, All}], 
Plot[
        func[x], {x, 1, 1000}, DisplayFunction -> $DisplayFunction, 
        PlotStyle -> 
          RGBColor[1, 0, 0]]
]

    
StartOsc = Select[
Position[
            Pulse, _?(# < 0 &)]
, #[[1]] > 
              10 &][[1, 1]]

    
Pulse = 
        ReplacePart[Pulse, 0, 
          Table[{m}, {m, 1, StartOsc}]];
func = 
        Interpolation[Pulse, 
          InterpolationOrder -> 5];

    
DisplayTogether[
ListPlot[Pulse, 
        PlotRange -> {{1, 100}, All}], 
Plot[
        func[x], {x, 1, 1000}, DisplayFunction -> $DisplayFunction, 
        PlotStyle -> RGBColor[1, 0, 0]], 
Plot[
        func2[x], {x, 60, 70}, DisplayFunction -> $DisplayFunction, 
        PlotStyle -> 
          RGBColor[0, 1, 0]]
]

\(\(\( (*\ 
      Verify\ coherence\ of\ Mercus\ input\ file\ and\ outfile\ report2\ *) \)\(
\)\(\(report2 = 
        ReadList["/afs/cern.ch/user/a/afabich/public/MERIT/transformer.dat", Table[Number, {40}]];\)

    \(check = 5;\)

    Do[{
If[
          Length[Position[\(Transpose[report2]\)[\([1]\)], 
                report2[\([k, 1]\)]]] ≠ 
            1, {Print[k, ":", 
              report2[\([k, 
                  1]\)], " several entries for single run"], \(check++\)}], 
pos = 
          First[First[
              Position[\(Transpose[RunList]\)[\([1]\)], 
                report2[\([k, 1]\)]]]], 
If[\(Transpose[
                RunList]\)[\([10, pos]\)] ≠ \(Transpose[
                report2]\)[\([8, k]\)], {Print[
              k, ": ", \(Transpose[RunList]\)[\([10, 
                  pos]\)], "<->", \(Transpose[report2]\)[\([8, 
                  k]\)], " different bunch number"], \(check++\)}]
}, {k, Length[report2]}]

    If[check ≠ 0, 
      Print["differences in databases!!!!"]]\)\)\)

report2[[310]]

SetDirectory["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/\
MERIT_libera"]
debug = 5;
report = {};
report2 = 
        ReadList["/afs/cern.ch/user/a/afabich/public/MERIT/transformer.dat\
", Table[Number, {40}]];
\(Do[{
Print[k, ": Working on run: ", 
          RunList[\([k, 1]\)]], 
data = 
          ReadList["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/MERIT_libera/" <> 
ToString[
                RunList[\([k, 4]\)]] <> "-" <> 
              ToString[RunList[\([k, 5]\)]] <> "-" <> 
              ToString[RunList[\([k, 6]\)]] <> "_" <> 
              ToString[RunList[\([k, 7]\)]] <> "-" <> 
              ToString[RunList[\([k, 8]\)]] <> "-" <> 
              ToString[RunList[\([k, 9]\)]] <> "_MEAS_result.txt", 
            Number], 
If[debugg ≥ 1, 
          Print["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/MERIT_libera/" <> 
ToString[
                RunList[\([k, 4]\)]] <> "-" <> 
              ToString[RunList[\([k, 5]\)]] <> "-" <> 
              ToString[RunList[\([k, 6]\)]] <> "_" <> 
              ToString[RunList[\([k, 7]\)]] <> "-" <> 
              ToString[RunList[\([k, 8]\)]] <> "-" <> 
              ToString[
                RunList[\([k, 
                    9]\)]] <> "_MEAS_result.txt"]], 
FullData = 
          Transpose[{Table[m, {m, Length[data]}], 
              data}], 
 (*\ \(ListPlot[
            data]\)\(,\)\ *) 
Integral = {data[\([1]\)]}, 
Do[
          AppendTo[Integral, 
            data[\([m]\)] + Integral[\([\(-1\)]\)]]
, {m, 
            2, Length[data]}], 
Hpeak = 
          Max[data], 

LowBound = {}, 
Do[
If[
            data[\([m - 1]\)] < Hpeak/2 ≤ data[\([m]\)], 
            AppendTo[LowBound, {m, data[\([m]\)]}]]
, {m, 
            2, Length[
              
              data]}], 
HiBound = {}, 
Do[
If[
            data[\([m - 1]\)] > Hpeak/2 ≥ data[\([m]\)], 
            AppendTo[HiBound, {m, data[\([m]\)]}]]
, {m, 2, 
            Length[data]}], 

PeakPos = 
          IntegerPart[
            N[\((\(Transpose[LowBound]\)[\([1]\)] + \(Transpose[
                        HiBound]\)[\([1]\)])\)/
                2]], 
BunchStart = {}, 
        BunchDel = 
          5, 
Do[
Do[
If[data[\([n - 1]\)] > data[\([n]\)] || 
                n == LowBound[\([m, 1]\)] - BunchDel, {AppendTo[
                  BunchStart, {n, data[\([n]\)]}], 
                Break[]}]
, {n, LowBound[\([m, 1]\)], 
              LowBound[\([m, 1]\)] - 
                BunchDel, \(-1\)}]
, {m, 
            Length[LowBound]}], 
Bunch = {}, 
Do[
AppendTo[
            Bunch, 
Table[{n, data[\([n]\)]}, {n, 
                BunchStart[\([m, 1]\)], 
                BunchStart[\([m, 1]\)] + 
                  22}]
]
, {m, 
            Length[BunchStart]}], 
Bunch = 
          Flatten[Bunch, 1], 

If[
          debug ≥ 2, 
          Export["/afs/cern.ch/user/a/afabich/public/MERIT/pic1.eps", 
            Show[GraphicsArray[{
DisplayTogether[
ListPlot[data, 
                      PlotRange -> {{LowBound[\([1, 1]\)] - 100, 
                            HiBound[\([\(-1\), 1]\)] + 300}, All}, 
                      PlotStyle -> RGBColor[0, 0, 0], 
                      PlotJoined -> False, 
                      hide], 
ListPlot[LowBound, 
                      PlotStyle -> RGBColor[1, 0, 0], 
                      PlotJoined -> False], 
ListPlot[
                      HiBound, PlotStyle -> RGBColor[1, 0, 0], 
                      PlotJoined -> False]
, 
                    DisplayFunction -> 
                      Identity], 
ListPlot[Integral, hide, 
                    FrameTicks -> {True, False, False, 
                        False}]
}, 
                ImageSize -> 
                  800]]]], 

ShowPeak = 
          1, 
If[debug ≥ 2, 
          Show[GraphicsArray[{
DisplayTogether[
ListPlot[data, 
                    PlotRange -> {{LowBound[\([ShowPeak, 1]\)] - 30, 
                          HiBound[\([ShowPeak, 1]\)] + 50}, All}, 
                    PlotStyle -> RGBColor[0, 0, 0], 
                    PlotJoined -> False, 
                    hide], 
ListPlot[Bunch, 
                    PlotStyle -> RGBColor[0, 1, 0], 
                    PlotJoined -> False], 
ListPlot[
                    LowBound, PlotStyle -> RGBColor[1, 0, 0], 
                    PlotJoined -> False], 
ListPlot[
                    HiBound, PlotStyle -> RGBColor[1, 0, 0], 
                    PlotJoined -> False], 
ListPlot[
                    Transpose[{PeakPos, data[\([PeakPos]\)]}], 
                    PlotStyle -> RGBColor[1, 0, 0], 
                    PlotJoined -> False], 
ListPlot[
                    BunchStart, PlotStyle -> RGBColor[1, 0, 0], 
                    PlotJoined -> False]
, 
                  DisplayFunction -> 
                    Identity], 
DisplayTogether[
ListPlot[data, 
                    PlotRange -> {{LowBound[\([\(-1\), 1]\)] - 30, 
                          HiBound[\([\(-1\), 1]\)] + 50}, All}, 
                    PlotStyle -> RGBColor[0, 0, 0], 
                    PlotJoined -> False, 
                    hide], 
ListPlot[LowBound, 
                    PlotStyle -> RGBColor[1, 0, 0], 
                    PlotJoined -> False], 
ListPlot[
                    HiBound, PlotStyle -> RGBColor[1, 0, 0], 
                    PlotJoined -> False], 
ListPlot[
                    BunchStart, PlotStyle -> RGBColor[1, 0, 0], 
                    PlotJoined -> False]
, 
                  DisplayFunction -> Identity]
}, 
              ImageSize -> 
                900]]], 

dataM = {}, 
dataH = {}, 

DataMod = data, 

Do[{
dataP = 
              Take[DataMod, {IntegerPart[\(LowBound[\([m, 1]\)] + HiBound[\([m, 1]\)]\)\/2] + 46, 
                  IntegerPart[\(LowBound[\([m, 1]\)] + HiBound[\([m, 1]\)]\)\/2] + 56}], 
funcP = 
              PolynomialFit[dataP, 5], 
dataMi = 
              Table[{x, funcP[x]}, {x, 1, Length[dataP], 
                  0.05}], 
dataMi = \(dataMi[\([First[
                      First[Position[dataMi, 
                          Max[dataMi]]]]]\)]\)[\([1]\)], 
AppendTo[dataM, 
              dataMi + \(LowBound[\([m, 1]\)] + HiBound[\([m, 1]\)]\)\/2 + 
                46 - 1], 
AppendTo[dataH, 
              funcP[dataMi]], 

FirstFit = {}, 
Do[{
DataMod = 
                  ReplacePart[DataMod, 
                    DataMod[\([p]\)] - 
                      func[p - dataM[\([m]\)] + PulseH]\ /PulseMax*
                        dataH[\([m]\)], p], 
AppendTo[
                  FirstFit, {p, 
                    func[p - dataM[\([m]\)] + PulseH]\ /PulseMax*
                      dataH[\([m]\)]}]
}, {p, 
                Ceiling[1 + dataM[\([m]\)] - PulseH], 
                Ceiling[1 + dataM[\([m]\)] - PulseH] + 
                  800}], 
If[
              Mod[m, 100] == 0, 
Show[
                GraphicsArray[{
MultipleListPlot[data, 
                      DataMod, 
                      FirstFit, 
PlotRange -> {{850, 
                            1250}, {\(-1000\), 1000}}, 
                      PlotStyle -> {RGBColor[0, 0, 0], 
                          RGBColor[1, 0, 0], RGBColor[0, 1, 0]}, 
                      hide], 
DisplayTogether[
ListPlot[dataP, 
                        PlotRange -> {All, 
                            All}], 
Plot[
                        funcP[x], {x, 1, Length[dataP]}, 
                        DisplayFunction -> $DisplayFunction, 
                        PlotStyle -> 
                          RGBColor[1, 0, 0]]
, 
                      DisplayFunction -> Identity]
}, 
                  ImageSize -> 
                    1000]]
]
}, {m, 1, 
            Length[LowBound]}], 

If[
          1 < RunList[\([k, 1]\)] < 
            20000, {\  (*\ \(If\)\([\)\(8999 < RunList[\([k, 1]\)] < 
                  12000, \(\({\)\(\ \)\)\)*) 
DataModFull = 
              Transpose[{Table[m, {m, Length[DataMod]}], 
                  DataMod}], 
Back3 = {
N[
                  Mean[Take[
                      Select[DataModFull, #[\([1]\)] < 
                            LowBound[\([1, 1]\)] - 10 &], 
                      50]]], 
N[
                  Mean[Take[
                      Select[DataModFull, #[\([1]\)] < 
                            LowBound[\([1, 1]\)] - 
                              10 &], \(-50\)]]]
}, 
Do[{
AppendTo[
                  Back3, 
N[
                    Mean[Select[DataModFull, 
                        HiBound[\([m - 1, 1]\)] + 10 < #[\([1]\)] < 
                            LowBound[\([m, 1]\)] - 
                              10 &]]]
]
}, {m, 2, Length[LowBound]}], 
Back3 = 
              Join[Back3, 
                N[Map[Mean, 
                    Partition[
                      Select[DataModFull, #[\([1]\)] > 
                            HiBound[\([\(-1\), 1]\)] + 10 &], 
                      50]]]], 
func4 = 
              Interpolation[Back3, 
                InterpolationOrder -> 1], 
DataMod = 
              Take[\(Transpose[DataModFull]\)[\([2]\)], {50, \(-1000\)}] - 
                func4[Take[\(Transpose[
                        DataModFull]\)[\([1]\)], {50, \(-1000\)}]]
}], 
NewFullIntegral = {DataMod[\([1]\)]}, 
Do[
AppendTo[
            NewFullIntegral, 
NewFullIntegral[\([\(-1\)]\)] + DataMod[\([m]\)]
]
, {m, 1, 
            Length[DataMod]}], 

level = {}, 
Do[
If[
            NewFullIntegral[\([m]\)] - 
                NewFullIntegral[\([m - 1]\)] ≤ 
              Max[NewFullIntegral]/
                  Length[LowBound]\ 0.003
, 
            AppendTo[
              level, {m, NewFullIntegral[\([m]\)]}]]
, {m, 
            2, HiBound[\([\(-1\), 1]\)] + 
              50}], 

levels = {}, 
levelpos = {}, 
temp = {level[\([1]\)]}, 
Do[
If[
            level[\([m, 2]\)] - temp[\([\(-1\), 2]\)] ≤ 
              Max[NewFullIntegral]/Length[LowBound]\ 0.2
, 
            AppendTo[temp, 
              level[\([m]\)]]
, {AppendTo[levels, temp], 
              AppendTo[levelpos, level[\([m - 1, 1]\)]], 
              temp = {level[\([m]\)]}}
]
, {m, 2, Length[level]}], 
AppendTo[levels, 
          temp], 
levels = 
          ReplacePart[levels, Take[levels[\([1]\)], \(-20\)], 
            1], 
AppendTo[levelpos, 
          level[\([\(-1\), 
              1]\)]], 

GrepLo[x_] := 
          If[Length[x] > 100, Take[x, 100], x], 
GrepHi[
            x_] := If[Length[x] > 100, Take[x, \(-100\)], 
            x], 
LoLevel = 
          Map[GrepLo, levels], 
HiLevel = 
          Map[GrepHi, levels], 
CalcAverage[x_] := 
          N[Mean[\(Transpose[x]\)[\([2]\)]]], 
CalcError[
            x_] := N[
            StandardDeviation[\(Transpose[
                  x]\)[\([2]\)]]], 
Intensities = 
          Drop[Map[CalcAverage, LoLevel], 1] - 
            Drop[Map[CalcAverage, 
                HiLevel], \(-1\)], 
Error = 
          Sqrt[Drop[Map[CalcError, LoLevel], 1]\^2 + 
              Drop[Map[CalcError, HiLevel], \(-1\)]\^2], 

If[
          debug ≥ 
            2, {
Export["/afs/cern.ch/user/a/afabich/public/MERIT/pic6.eps", 
Show[
                GraphicsArray[{
DisplayTogether[
ListPlot[data, 
                        PlotRange -> {{820, 2200}, {\(-500\), 500}}, 
                        PlotLabel -> "Error plot"], 
ListPlot[DataMod, 
                        PlotStyle -> 
                          RGBColor[1, 0, 0]]
, 
                      DisplayFunction -> 
                        Identity], 
DisplayTogether[
ListPlot[data, 
                        PlotRange -> {{820, 1050}, {\(-500\), 500}}, 
                        PlotLabel -> "Error plot"], 
ListPlot[DataMod, 
                        PlotStyle -> 
                          RGBColor[1, 0, 0]]
, 
                      DisplayFunction -> Identity]
, 
                    DisplayTogether[
ListPlot[data, 
                        PlotRange -> {{LowBound[\([\(-1\), 1]\)] - 100, 
                              HiBound[\([\(-1\), 1]\)] + 20}, {\(-500\), 
                              1500}}, 
                        PlotLabel -> "Error plot"], 
ListPlot[DataMod, 
                        PlotStyle -> 
                          RGBColor[1, 0, 0]]
, 
                      DisplayFunction -> Identity]
}, 
                  ImageSize -> 
                    900]]
]}], 

If[debug ≥ 1, 
          Export["/afs/cern.ch/user/a/afabich/public/MERIT/pic7.eps", 
            Show[GraphicsArray[{
MultipleListPlot[
                    NewFullIntegral, 
                    level, 
Flatten[HiLevel, 1], 
                    Flatten[LoLevel, 1], 
Join[
                      Transpose[{\(Transpose[HiBound]\)[\([1]\)] - 50, 
                          Drop[Map[CalcAverage, LoLevel], 
                            1]}], 
Transpose[{\(Transpose[
                                LowBound]\)[\([1]\)] - 50, 
                          Drop[Map[CalcAverage, 
                              HiLevel], \(-1\)]}]]
, 
                    PlotRange -> {{850 - 50, 
                          HiBound[\([\(-1\), 1]\)] + 200}, All}, 
                    PlotJoined -> {True, False, False, False, False}, 
                    SymbolShape -> {None, PlotSymbol[Triangle], 
                        PlotSymbol[Triangle], PlotSymbol[Triangle], 
                        PlotSymbol[Box]}, 
                    SymbolStyle -> {RGBColor[0, 0, 0], 
                        RGBColor[1, 0, 0], RGBColor[0, 0, 1], 
                        RGBColor[0, 0, 1], RGBColor[0, 1, 0]}, 
                    hide], 
MultipleListPlot[
                    NewFullIntegral, 
                    level, 
Flatten[HiLevel, 1], 
                    Flatten[LoLevel, 1], 
Join[
                      Transpose[{\(Transpose[HiBound]\)[\([1]\)] - 50, 
                          Drop[Map[CalcAverage, LoLevel], 
                            1]}], 
Transpose[{\(Transpose[
                                LowBound]\)[\([1]\)] - 50, 
                          Drop[Map[CalcAverage, 
                              HiLevel], \(-1\)]}]]
, 
                    PlotRange -> {{LowBound[\([\(-1\), 1]\)] - 100 - 50, 
                          HiBound[\([\(-1\), 1]\)] + 100 - 50}, All}, 
                    PlotJoined -> {True, False, False, False, False}, 
                    SymbolShape -> {None, PlotSymbol[Triangle], 
                        PlotSymbol[Triangle], PlotSymbol[Triangle], 
                        PlotSymbol[Box]}, 
                    SymbolStyle -> {RGBColor[0, 0, 0], 
                        RGBColor[1, 0, 0], RGBColor[0, 0, 1], 
                        RGBColor[0, 0, 1], RGBColor[0, 1, 0]}, 
                    hide]
}, 
                ImageSize -> 
                  1100]]
]], 

If[
          debug ≥ 2, 
ListPlot[Intensities, 
            PlotLabel -> "Bunch intensities"]
], 

result = 
          Join[{RunList[\([k, 1]\)]}, 
            Take[RunList[\([k]\)], {4, 9}], {Length[
                LowBound]}, 
Intensities, 
            Table[0, {16 - Length[LowBound]}], 
N[
              IntegerPart[Error/Intensities\ 100\ 100]/100], 
            Table[0, {16 - 
                  Length[LowBound]}]
], 
If[
          result[\([8]\)] > 16 || 
            RunList[\([k, 10]\)] ≠ result[\([8]\)] || 
            Length[result] ≠ 40, {
result = 
              Join[{result[\([1]\)]}, Table[0, {39}]], 
            Print["FAULTY RUN ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"]
}], 
PosRep = 
          Position[\(Transpose[report2]\)[\([1]\)], 
            RunList[\([k, 1]\)]], 
If[PosRep == {}, 
          AppendTo[report2, result], 
          report2 = 
            ReplacePart[report2, result, 
              First[First[
                  PosRep]]]], 
 (*\ \(Print["Errors: ", 
            N[IntegerPart[Error/Intensities\ 100\ 100]/
                100]]\)\(,\)*) 
Print["Maximum error: ", Max[N[IntegerPart[Error/Intensities\ 100\ 100]/
                100]]]
, 
        Export["/afs/cern.ch/user/a/afabich/public/MERIT/transformer.dat", 
report2
], 
Export["/afs/cern.ch/user/a/afabich/public/MERIT/runmod/" <> ToString[RunList[\([k, 1]\)]] <> ".dat", 
          N[DataMod]]
}, {k, 160, 
        160}]
\)
Print["++++++"]
Print["The report contains ", Dimensions[
          report2][[1]], " entries for beam runs. Does this \
correspond with the requested evaluation?"]

Take[report2, {275, 276}]
N[%[[2]]/%[[1]]]

Join[
      Transpose[{Transpose[HiBound][[1]], 
          Drop[Map[CalcAverage, LoLevel], 
            1]}], 
Transpose[{Transpose[
              LowBound][[1]], 
          Drop[Map[CalcAverage, HiLevel], -1]}]]

data16 = 
        ReadList["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/\
MERIT_libera/" <> 
ToString[
              RunList[[141, 4]]] <> "-" <> 
            ToString[RunList[[141, 5]]] <> "-" <> 
            ToString[RunList[[141, 6]]] <> "_" <> 
            ToString[RunList[[141, 7]]] <> "-" <> 
            ToString[RunList[[141, 8]]] <> "-" <> 
            ToString[RunList[[141, 9]]] <> "_MEAS_result.txt", 
          Number];
Export["/afs/cern.ch/user/a/afabich/public/MERIT/Run_" <> 
        ToString[
          RunList[[141, 
              1]]] <> "_Clearly16.gif", 
ListPlot[
        data16, PlotRange -> {{800, 1200}, All}, 
        PlotLabel -> "PS in h16"]
]
data8 = 
        ReadList["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/\
MERIT_libera/" <> 
ToString[
              RunList8[[14, 4]]] <> "-" <> 
            ToString[RunList8[[14, 5]]] <> "-" <> 
            ToString[RunList8[[14, 6]]] <> "_" <> 
            ToString[RunList8[[14, 7]]] <> "-" <> 
            ToString[RunList8[[14, 8]]] <> "-" <> 
            ToString[RunList8[[14, 9]]] <> "_MEAS_result.txt", 
          Number];
Export["/afs/cern.ch/user/a/afabich/public/MERIT/Run_" <> 
        ToString[
          RunList8[[14, 
              1]]] <> "_Clearly8.gif", 
ListPlot[
        data8, PlotRange -> {{800, 1200}, All}, 
        PlotLabel -> "PS in h8"]
]

\(Do[{
dataOff = 
          ReadList["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/MERIT_libera/" <> 
ToString[
                RunListOff[\([k, 4]\)]] <> "-" <> 
              ToString[RunListOff[\([k, 5]\)]] <> "-" <> 
              ToString[RunListOff[\([k, 6]\)]] <> "_" <> 
              ToString[RunListOff[\([k, 7]\)]] <> "-" <> 
              ToString[RunListOff[\([k, 8]\)]] <> "-" <> 
              ToString[RunListOff[\([k, 9]\)]] <> "_MEAS_result.txt", 
            Number], 
Export["/afs/cern.ch/user/a/afabich/public/MERIT/Run_" <> ToString[RunListOff[\([k, 
                  1]\)]] <> ".gif", 
ListPlot[dataOff, 
            PlotRange -> {{800, 1200}, All}, 
            PlotLabel -> "PS in harminc unknown"]
]
}, {k, Length[RunListOff]}]\)

SetDirectory["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/\
MERIT_libera"]
debug = 0;
report = {};
report2 = 
        ReadList["/afs/cern.ch/user/a/afabich/public/MERIT/transformer.dat\
", Table[Number, {40}]];
saturation = {};
\(Do[{
data = 
          ReadList["/nfs/cs-ccr-nfs4/data/slopsdata/2007/MD_DATA/CPS_MERIT/MERIT_libera/" <> 
ToString[
                RunList[\([k, 4]\)]] <> "-" <> 
              ToString[RunList[\([k, 5]\)]] <> "-" <> 
              ToString[RunList[\([k, 6]\)]] <> "_" <> 
              ToString[RunList[\([k, 7]\)]] <> "-" <> 
              ToString[RunList[\([k, 8]\)]] <> "-" <> 
              ToString[RunList[\([k, 9]\)]] <> "_MEAS_result.txt", 
            Number], 
If[
          Max[data] ≥ 
            7900, 
{
If[
              debug ≥ 
                1, {
Print[k, ": Working on run: ", 
                  RunList[\([k, 1]\)], " : ", 
                  Max[data]], 
maxpos = 
                  First[First[
                      Position[data, 
                        Max[data]]]], 
ListPlot[
                  Take[data, {maxpos - 100, 
                      maxpos + 
                        300}]]
}], 
AppendTo[saturation, 
              Join[report2[\([k]\)], {1}]]
}, 
          AppendTo[saturation, 
            Join[report2[\([k]\)], {0}]]
]
}, {k, 1, 310}]\)
\(Export["/afs/cern.ch/user/a/afabich/public/MERIT/transformerWithSaturation.dat", 
saturation
]\)

report2 = 
        ReadList["/afs/cern.ch/user/a/afabich/public/MERIT/transformer.dat\
", Table[Number, {40}]];
rep = 
        Select[report2, #[[8]] ≠ 
              0 &];
index = 
        Select[
Transpose[{
Table[
                m, {m, Length[Transpose[
                        rep][[1]]]}], 
Transpose[
                  rep][[1]]
}]
 , 
          Mod[#[[2]] - 1, 1000] ≤ 
              1 &];
index = 
        Transpose[{Transpose[
                index][[1]], 
IntegerPart[Transpose[
                    index][[2]]/1000]}];
ListPlot[
Transpose[
          rep][[8]]
, PlotJoined -> False, 
      FrameTicks -> {index, True, False, False}, ImageSize -> 900, 
      FrameLabel -> {"# day index", "# of bunches"}]