Changeset 17893 for trunk/Mars/fact


Ignore:
Timestamp:
05/26/14 16:26:53 (10 years ago)
Author:
Daniela Dorner
Message:
new callisto not using the light pulser, using a new extractor and pixel saturation treatment, not using delays
File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/Mars/fact/analysis/callisto.C

    r17885 r17893  
    11#include "MLogManip.h"
    22
    3 int callisto(const char *seqfile="seq/2012/01/23/20120123_023.seq", const char *outpath = "output")
     3int callisto9(const char *seqfile="seq/2012/01/23/20120123_023.seq", const char *outpath = "output", bool use_delays=false)
    44{
    55    // ======================================================
    … …  
    1212
    1313    // map file to use (get that from La Palma!)
    14     const char *map = usemap ? "/home/fact/FACT++/FACTmap111030.txt" : NULL;
    15 /*
    16     Bool_t maximum = kTRUE;
    17 
    18     const char *lp_template    = maximum ?
    19         "template-lp-extractor-maximum.root" :
    20         "template-lp-extractor-leading-edge.root";
    21 
    22     const char *pulse_template = "template-pulse.root";
    23 */
     14//    const char *map = usemap ? "/home/fact/FACT++/FACTmap111030.txt" : NULL;
     15    const char *map = usemap ? "/scratch/fact/FACTmap111030.txt" : NULL;
     16
    2417    // ------------------------------------------------------
    25 
    26     //bool use_delays=false;
    27 
    28     int spike_removal=3;
    29 
    30     // The gain as extracted from our dark count spectra
    31     double gain = 258;
    32 
    33     // ------------------------------------------------------
    34 
    35     // Extraction range in slices. It will always(!) contain the full range
    36     // of integration
    37     const int first_slice =  25; //  10ns
    38     const int last_slice  = 225; // 125ns
    39 
    40     // Note that rise and fall time mean different things whether you use IntegralFixed or IntegraRel:
    41     //
    42     //  IntegralFixed:
    43     //    * fRiseTime: Number of slices left  from arrival time
    44     //    * fFallTime: Number of slices right from arrival time
    45     //  IntegralRel:
    46     //    * fRiseTime: Number of slices left  from maximum time
    47     //    * fFallTime: Number of slices right from maximum time
    48     //
    49 /*
    50     const int rise_time_cal = maximum ?  40 :  10; // was 13;   5ns
    51     const int fall_time_cal = maximum ? 120 : 160; // was 23;  80ns
    52 
    53     const int rise_time_dat = maximum ?  10 :   2; // was 13; was 10;   1ns
    54     const int fall_time_dat = maximum ?  40 :  48; // was 23; was 40;  24ns
    55 
    56     // Extraction type: Extract integral and half leading edge
    57 
    58     const int type = maximum ? (MExtralgoSpline::kAmplitudeRel) : (MExtralgoSpline::kIntegralFixed);
    59     //const int type = MExtralgoSpline::kIntegralFixed;
    60 
    61 
    62     const double heighttm   = 0.5; // IntegralAbs { 1.5pe * 9.6mV/pe } / IntegralRel { 0.5 }
    63 */
    64     Long_t max  =    0;  // All
    65     Long_t max0 =  max;  // Time marker
    66     //Long_t max1 =  max;  // Light pulser
    67     //Long_t max2 = 3000;  // Calibration ratio
    68     //Long_t max3 =  max;  // Pedestal Rndm
    69     //Long_t max4 =  max;  // Pedestal Ext
    70     Long_t max5 =  max;  // Data
    71 
    72     // ======================================================
    73 
    74     if (map && gSystem->AccessPathName(map, kFileExists))
    75     {
    76         gLog << err << "ERROR - Cannot access mapping file '" << map << "'" << endl;
    77         return 1;
    78     }
    79 
    80     // The sequence file which defines the files for the analysis
    81     MSequence seq(seqfile);
    82     if (!seq.IsValid())
    83     {
    84         gLog << err << "ERROR - Sequence '" << seqfile << "' invalid!" << endl;
    85         return 2;
    86     }
    87 
    88     // --------------------------------------------------------------------------------
    89 
    90     gLog.Separator("Callisto");
    91     gLog << all << "Calibrate data of sequence '" << seq.GetFileName() << "'" << endl;
    92     gLog << endl;
    93 
    94     // ------------------------------------------------------
    95 
    96     ostringstream drsname;
    97     drsname << gSystem->DirName(seqfile) << "/";
    98     drsname << seq.GetNight().GetNightAsInt() << "_";
    99     drsname << Form("%03d", seq.GetDrsSequence()) << ".drs.seq";
    100 
    101     MSequence drs(drsname.str().c_str());
    102     if (!drs.IsValid())
    103     {
    104         gLog << err << "ERROR - DRS sequence invalid!" << endl;
    105         return 3;
    106     }
    107 
    108     gLog << all << "DRS sequence file: " << drsname.str() << '\n' << endl;
    109 
    110     TString drsfile = seq.GetFileName(0, MSequence::kRawDrs);
    111     if (drsfile.IsNull())
    112     {
    113         gLog << err << "No DRS file available in sequence." << endl;
    114         return 4;
    115     }
    116 
    117     TString timfile = drs.GetFileName(0, MSequence::kFitsDat);
    118     TString drs1024 = drs.GetFileName(0, MSequence::kFitsDrs);
    119     TString calfile = seq.GetFileName(0, MSequence::kFitsCal);
    120     //TString pedfile = seq.GetFileName(0, MSequence::kFitsPed);
    121 
    122     gLog << all;
    123     gLog << "DRS calib     300: " << drsfile << '\n';
    124     gLog << "DRS calib    1024: " << drs1024 << "\n\n";
    125 
    126     MDrsCalibration drscalib300;
    127     if (!drscalib300.ReadFits(drsfile.Data()))
    128         return 5;
    129 
    130     MDrsCalibration drscalib1024;
    131     if (!drscalib1024.ReadFits(drs1024.Data()))
    132         return 6;
    133 
    134     gLog << all;
    135     gLog << "Time calibration : " << timfile << '\n';
    136     gLog << "Light Pulser file: " << calfile << '\n' << endl;
    137     //gLog << "Pedestal     file: " << pedfile << '\n';
    138 
    139     // ------------------------------------------------------
    140 
    141     MDirIter iter;
    142     if (seq.GetRuns(iter, MSequence::kFitsDat)<=0)
    143     {
    144         gLog << err << "ERROR - Sequence valid but without files." << endl;
    145         return 7;
    146     }
    147     iter.Print();
    148 
    149     // ======================================================
    150 
    151     MStatusArray arrt, arrp;
    152 
    153     TFile ft(lp_template);
    154     if (arrt.Read()<=0)
    155     {
    156         gLog << err << "ERROR - Reading LP template from " << lp_template << endl;
    157         return 100;
    158     }
    159 
    160     MHCamera *lpref = (MHCamera*)arrt.FindObjectInCanvas("ExtCalSig;avg", "MHCamera", "Cam");
    161     if (!lpref)
    162     {
    163         gLog << err << "ERROR - LP Template not found in " << lp_template << endl;
    164         return 101;
    165     }
    166     lpref->SetDirectory(0);
    167 
    168     MHCamera *gain = (MHCamera*)arrt.FindObjectInCanvas("gain", "MHCamera", "Gain");
    169     if (!gain)
    170     {
    171         gLog << err << "ERROR - Gain not found in " << lp_template << endl;
    172         return 101;
    173     }
    174     gain->SetDirectory(0);
    175 
    176     TFile fp(pulse_template);
    177     if (arrp.Read()<=0)
    178     {
    179         gLog << err << "ERROR - Reading Pulse template from " << pulse_template << endl;
    180         return 102;
    181     }
    182 
    183     TH1F *hpulse = (TH1F*)arrp.FindObjectInCanvas("hPixelEdgeMean0_0", "TH1F", "cgpPixelPulses0");
    184     if (!hpulse)
    185     {
    186         gLog << err << "ERROR - Pulse Template not found in " << pulse_template << endl;
    187         return 103;
    188     }
    189     hpulse->SetDirectory(0);
    190 
    191     // ======================================================
    19218
    19319    MStatusDisplay *d = new MStatusDisplay;
    … …  
    20127    badpixels[1208].SetUnsuitable(MBadPixelsPix::kUnsuitable);
    20228    badpixels[1399].SetUnsuitable(MBadPixelsPix::kUnsuitable);
    203 
    20429    //  Twin pixel
    20530    //     113
    … …  
    21035    //    1393
    21136
    212     MDrsCalibrationTime timecam;
     37    // ------------------------------------------------------
     38
     39    // ------------------------------------------------------
     40
     41    // Calib: 51 / 90 / 197 (20% TH)
     42    // Data:  52 / 64 / 104 (20% TH)
     43
     44    // Extraction range in slices. It will always(!) contain the full range
     45    // of integration
     46    const int first_slice =  25; //  10ns
     47    const int last_slice  = 225; // 125ns
     48
     49    const double heighttm   = 0.5; // IntegralAbs { 1.5pe * 9.6mV/pe } / IntegralRel { 0.5 }
     50
     51    Long_t max  =    0;  // All
     52    Long_t max0 =  max;  // Time marker
     53    Long_t max1 =  max;  // Light pulser
     54    //Long_t max2 = 3000;  // Calibration ratio
     55    Long_t max3 =  max;  // Pedestal Rndm
     56    Long_t max4 =  max;  // Pedestal Ext
     57    Long_t max5 =  max;  // Data
     58
     59    // ========================= Result ==================================
     60
     61    //double scale = 0.1;
     62    double scale = 0.1024;
     63
     64    // ======================================================
     65
     66    if (map && gSystem->AccessPathName(map, kFileExists))
     67    {
     68        gLog << err << "ERROR - Cannot access mapping file '" << map << "'" << endl;
     69        return 1;
     70    }
     71
     72    // The sequence file which defines the files for the analysis
     73    MSequence seq(seqfile);
     74    if (!seq.IsValid())
     75    {
     76        gLog << err << "ERROR - Sequence '" << seqfile << "' invalid!" << endl;
     77        return 2;
     78    }
     79
     80    // --------------------------------------------------------------------------------
     81
     82    gLog.Separator("Callisto");
     83    gLog << all << "Calibrate data of sequence '" << seq.GetFileName() << "'" << endl;
     84    gLog << endl;
     85
     86    // ------------------------------------------------------
     87
     88    ostringstream drsname;
     89    drsname << gSystem->DirName(seqfile) << "/";
     90    drsname << seq.GetNight().GetNightAsInt() << "_";
     91    drsname << Form("%03d", seq.GetDrsSequence()) << ".drs.seq";
     92
     93    MSequence drs(drsname.str().c_str());
     94    if (!drs.IsValid())
     95    {
     96        gLog << err << "ERROR - DRS sequence invalid!" << endl;
     97        return 3;
     98    }
     99
     100    gLog << all << "DRS sequence file: " << drsname.str() << '\n' << endl;
     101
     102    TString drsfile = seq.GetFileName(0, MSequence::kRawDrs);
     103    if (drsfile.IsNull())
     104    {
     105        gLog << err << "No DRS file available in sequence." << endl;
     106        return 4;
     107    }
     108
     109    TString timfile = drs.GetFileName(0, MSequence::kFitsDat);
     110    TString drs1024 = drs.GetFileName(0, MSequence::kFitsDrs);
     111    TString pedfile = seq.GetFileName(0, MSequence::kFitsPed);
     112    TString calfile = seq.GetFileName(0, MSequence::kFitsCal);
     113
     114    gLog << all;
     115    gLog << "DRS calib     300: " << drsfile << '\n';
     116    gLog << "DRS calib    1024: " << drs1024 << "\n\n";
     117
     118    MDrsCalibration drscalib300;
     119    if (!drscalib300.ReadFits(drsfile.Data()))
     120        return 5;
     121
     122    MDrsCalibration drscalib1024;
     123    if (!drscalib1024.ReadFits(drs1024.Data()))
     124        return 6;
     125
     126    gLog << all;
     127    gLog << "Time calibration : " << timfile << '\n';
     128    gLog << "Pedestal     file: " << pedfile << '\n';
     129    gLog << "Light Pulser file: " << calfile << '\n' << endl;
     130
     131    // ------------------------------------------------------
     132
     133    MDirIter iter;
     134    if (seq.GetRuns(iter, MSequence::kFitsDat)<=0)
     135    {
     136        gLog << err << "ERROR - Sequence valid but without files." << endl;
     137        return 7;
     138    }
     139    iter.Print();
     140
     141    // ======================================================
     142
     143/*
     144    MStatusArray arrt, arrp;
     145
     146    TFile ft(lp_template);
     147    if (arrt.Read()<=0)
     148    {
     149        gLog << err << "ERROR - Reading LP template from " << lp_template << endl;
     150        return 100;
     151    }
     152
     153    MHCamera *lpref = (MHCamera*)arrt.FindObjectInCanvas("ExtCalSig;avg", "MHCamera", "Cam");
     154    if (!lpref)
     155    {
     156        gLog << err << "ERROR - LP Template not found in " << lp_template << endl;
     157        return 101;
     158    }
     159    lpref->SetDirectory(0);
     160
     161    MHCamera *gain = (MHCamera*)arrt.FindObjectInCanvas("gain", "MHCamera", "Gain");
     162    if (!gain)
     163    {
     164        gLog << err << "ERROR - Gain not found in " << lp_template << endl;
     165        return 101;
     166    }
     167    gain->SetDirectory(0);
     168
     169    TFile fp(pulse_template);
     170    if (arrp.Read()<=0)
     171    {
     172        gLog << err << "ERROR - Reading Pulse template from " << pulse_template << endl;
     173        return 102;
     174    }
     175
     176    TH1F *hpulse = (TH1F*)arrp.FindObjectInCanvas("hPixelEdgeMean0_0", "TH1F", "cgpPixelPulses0");
     177    if (!hpulse)
     178    {
     179        gLog << err << "ERROR - Pulse Template not found in " << pulse_template << endl;
     180        return 103;
     181    }
     182    hpulse->SetDirectory(0);
     183*/
     184    // ======================================================
    213185
    214186    // Plot the trigger pattern rates vs. run-number
    … …  
    223195    // hrate.DefaultLabelY("ERROR");
    224196
    225     Bool_t isinteg =
    226         (type&MExtralgoSpline::kIntegral)    ||
    227         (type&MExtralgoSpline::kFixedWidth)  ||
    228         (type&MExtralgoSpline::kDynWidth)
    229         ? kTRUE : kFALSE;
     197    MDrsCalibrationTime timecam;
    230198
    231199    gStyle->SetOptFit(kTRUE);
    … …  
    241209    plist0.AddToList(&tlist0);
    242210    plist0.AddToList(&drscalib1024);
    243     plist0.AddToList(&badpixels);
    244211    plist0.AddToList(&timecam);
    245212
    … …  
    257224
    258225    MDrsCalibApply drsapply0;
    259     drsapply0.SetRemoveSpikes(spike_removal);
     226    //drsapply0.SetRemoveSpikes(4);
    260227
    261228    MFillH fill0("MHDrsCalibrationTime");
    … …  
    281248
    282249    // ======================================================
    283 /*
     250
    284251    gLog << endl;
    285252    gLog.Separator("Processing external light pulser run");
    … …  
    308275
    309276    MDrsCalibApply drsapply1;
    310     drsapply1.SetRemoveSpikes(spike_removal);
     277    //drsapply1.SetRemoveSpikes(4);
     278
     279    MFilterData filterdata1;
    311280
    312281    // ---
    313282
    314     MExtractTimeAndChargeSpline extractor1b("ExtractPulse");
     283    MExtractFACT extractor1b("ExtractPulse");
    315284    extractor1b.SetRange(first_slice, last_slice);
    316     extractor1b.SetRiseTimeHiGain(rise_time_cal);
    317     extractor1b.SetFallTimeHiGain(fall_time_cal);
    318     extractor1b.SetHeightTm(heighttm);
    319     extractor1b.SetChargeType(type);
    320     extractor1b.SetSaturationLimit(600000);
    321285    extractor1b.SetNoiseCalculation(kFALSE);
    322 
    323     MExtractTimeAndChargeSpline extractor1c("ExtractAmplitude");
    324     extractor1c.SetRange(first_slice, last_slice);
    325     extractor1c.SetChargeType(MExtralgoSpline::kAmplitude);
    326     extractor1c.SetSaturationLimit(600000);
    327     extractor1c.SetNoiseCalculation(kFALSE);
    328     extractor1c.SetNameSignalCam("Amplitude");
    329     extractor1c.SetNameTimeCam("AmplitudePos");
    330 
    331     // ---
    332 
    333     MHCamEvent evt1a(5, "CalRatio", "Ratio per slice between integrated signal and amplitude;; r [1/n]");
    334     evt1a.SetNameSub("Amplitude", kTRUE);
    335     MFillH fill1a(&evt1a, "MExtractedSignalCam", "FillRatio");
    336     fill1a.SetDrawOption("gaus");
    337 
    338     MParameterD ratio1a;
    339     ratio1a.SetVal(1./(fall_time_cal+rise_time_cal));
    340     fill1a.SetWeight(&ratio1a);
    341286
    342287    // ---
    … …  
    391336    tlist1.AddToList(&drsapply1);
    392337    tlist1.AddToList(&cont1);
     338    tlist1.AddToList(&filterdata1);
    393339    tlist1.AddToList(&extractor1b);
    394     if (isinteg)
    395     {
    396         tlist1.AddToList(&extractor1c);
    397         tlist1.AddToList(&fill1a);
    398     }
    399340    tlist1.AddToList(&calctm1a);
    400341    tlist1.AddToList(&treat1);
    … …  
    416357        timecam.SetDelays(*evt1h.GetHist());
    417358
    418     // ========================= Result ==================================
    419 
    420     Double_t avgS = evt1f.GetHist()->GetMean();
    421     Double_t medS = evt1f.GetHist()->GetMedian();
    422     Double_t rmsS = evt1f.GetHist()->GetRMS();
    423     Double_t maxS = evt1f.GetHist()->GetMaximum();
    424 
    425     MArrayF der1(hpulse->GetNbinsX());
    426     MArrayF der2(hpulse->GetNbinsX());
    427 
    428     MExtralgoSpline spline(hpulse->GetArray()+1, hpulse->GetNbinsX(),
    429                            der1.GetArray(), der2.GetArray());
    430     spline.SetRiseFallTime(rise_time_dat, fall_time_dat);
    431     spline.SetExtractionType(type);
    432     spline.SetHeightTm(heighttm);
    433 
    434     spline.Extract(hpulse->GetMaximumBin()-1);
    435 
    436     // The pulser signal is most probably around 400mV/9.5mV
    437     // IntegraFixed 2/48 corresponds to roughly 215mV*50slices
    438     Double_t scale = 1./spline.GetSignal();
    439 
    440     MArrayD calib(1440);
    441     for (int i=0; i<1440; i++)
    442     {
    443         Double_t g = gain->GetBinContent(i+1)>0.5 ? gain->GetBinContent(i+1) : 1;
    444         if (evt1f.GetHist()->GetBinContent(i+1)>0 && !badpixels[i].IsUnsuitable())
    445             calib[i] = lpref->GetBinContent(i+1) / evt1f.GetHist()->GetBinContent(i+1) / g;
    446     }
    447 
    448     gROOT->SetSelectedPad(0);
    449     d->AddTab("PulseTemp");
    450     gPad->SetGrid();
    451     hpulse->SetNameTitle("Pulse", "Single p.e. pulse template");
    452     hpulse->SetDirectory(0);
    453     hpulse->SetLineColor(kBlack);
    454     hpulse->DrawCopy();
    455 
    456     TAxis *ax = hpulse->GetXaxis();
    457 
    458     Double_t w = hpulse->GetBinWidth(1);
    459     Double_t T = w*(spline.GetTime()+0.5)       +ax->GetXmin();
    460     Double_t H = w*(hpulse->GetMaximumBin()+0.5)+ax->GetXmin();
    461 
    462     TLine line;
    463     line.SetLineColor(kRed);
    464     line.DrawLine(T-rise_time_dat*w, spline.GetHeight(),
    465                   T+fall_time_dat*w, spline.GetHeight());
    466     line.DrawLine(T, spline.GetHeight()/4, T, 3*spline.GetHeight()/4);
    467     line.DrawLine(T-rise_time_dat*w, 0,
    468                   T-rise_time_dat*w, spline.GetHeight());
    469     line.DrawLine(T+fall_time_dat*w, 0,
    470                   T+fall_time_dat*w, spline.GetHeight());
    471 
    472     TGraph gg;
    473     for (int ix=1; ix<=hpulse->GetNbinsX(); ix++)
    474         for (int i=0; i<10; i++)
    475         {
    476             Double_t x = hpulse->GetBinLowEdge(ix)+i*hpulse->GetBinWidth(ix)/10.;
    477             gg.SetPoint(gg.GetN(), x+w/2, spline.EvalAt(ix-1+i/10.));
    478         }
    479 
    480     gg.SetLineColor(kBlue);
    481     gg.SetMarkerColor(kBlue);
    482     gg.SetMarkerStyle(kFullDotMedium);
    483     gg.DrawClone("L");
    484 
    485     gROOT->SetSelectedPad(0);
    486     d->AddTab("CalConst");
    487     MGeomCamFACT fact;
    488     MHCamera hcalco(fact);
    489     hcalco.SetName("CalConst");
    490     hcalco.SetTitle(Form("Relative calibration constant [%.0f/pe]", 1./scale));
    491     hcalco.SetCamContent(calib);
    492     hcalco.SetAllUsed();
    493     //hcalco.Scale(scale);
    494     hcalco.DrawCopy();
    495 */
    496 
    497     TF1 f("f", "[0]*(1-1/(1+exp(x/[1])))*exp(-x/[2])", -15, 150);
    498     f.SetParameter(0, gain*0.05143);
    499     f.SetParameter(1, 1.075);
    500     f.SetParameter(2, 19.3);
    501     f.SetNpx(2*165); //2GHz
    502 
    503     // Convert to graph
    504     TGraph g(&f);
    505 
    506     // Convert to float
    507     MArrayF x(g.GetN());
    508     MArrayF y(g.GetN());
    509     for (int i=0; i<g.GetN(); i++)
    510     {
    511         x[i] = g.GetX()[i];
    512         y[i] = g.GetY()[i];
    513     }
    514 
    515     MFilterData filter;
    516     filter.Filter(1, g.GetN(), y.GetArray(), y.GetArray());
    517 
    518     // Define spline
    519     MArrayF der1(g.GetN());
    520     MArrayF der2(g.GetN());
    521 
    522     MExtralgoSpline spline(y.GetArray(), y.GetSize(),
    523                            der1.GetArray(), der2.GetArray());
    524 
    525     spline.SetExtractionType(MExtralgoSpline::kAmplitudeRel);
    526     spline.SetHeightTm(0.5);
    527 
    528     // Estimate where the maximum is and extract signal
    529     Long64_t maxi = TMath::LocMax(y.GetSize(), y.GetArray());
    530     spline.Extract(maxi);
    531 
    532     // Scale factor for signal extraction
    533     double scale = 1./spline.GetSignal();
    534 
    535359    // ======================================================
    536 /*
     360
    537361    gLog << endl;
    538362    gLog.Separator("Extracting random pedestal");
    … …  
    563387
    564388    MDrsCalibApply drsapply3;
    565     drsapply3.SetRemoveSpikes(spike_removal);
     389    //drsapply3.SetRemoveSpikes(4);
     390
     391    MFilterData filterdata3;
    566392
    567393    //---
    568394
    569     MExtractTimeAndChargeSpline extractor3;
     395    MExtractFACT extractor3;
    570396    extractor3.SetRange(first_slice, last_slice);
    571     extractor3.SetRiseTimeHiGain(rise_time_dat);
    572     extractor3.SetFallTimeHiGain(fall_time_dat);
    573     extractor3.SetHeightTm(heighttm);
    574     extractor3.SetChargeType(type);
    575     extractor3.SetSaturationLimit(600000);
    576397    extractor3.SetNoiseCalculation(kTRUE);
    577398
    578399    MCalibrateFact conv3;
    579400    conv3.SetScale(scale);
    580     conv3.SetCalibConst(calib);
     401    //conv3.SetCalibConst(calib);
    581402
    582403    MBadPixelsTreat treat3;
    … …  
    597418    tlist3.AddToList(&drsapply3);
    598419    tlist3.AddToList(&cont3);
     420    tlist3.AddToList(&filterdata3);
    599421    tlist3.AddToList(&extractor3);
    600422//    tlist3.AddToList(&fill3a);
    … …  
    637459
    638460    MDrsCalibApply drsapply4;
    639     drsapply4.SetRemoveSpikes(spike_removal);
    640 
    641     MExtractTimeAndChargeSpline extractor4;
     461    //drsapply4.SetRemoveSpikes(4);
     462
     463    MFilterData filterdata4;
     464
     465    MExtractFACT extractor4;
    642466    extractor4.SetRange(first_slice, last_slice);
    643     extractor4.SetRiseTimeHiGain(rise_time_dat);
    644     extractor4.SetFallTimeHiGain(fall_time_dat);
    645     extractor4.SetHeightTm(heighttm);
    646     extractor4.SetChargeType(type);
    647     extractor4.SetSaturationLimit(600000);
    648467    extractor4.SetNoiseCalculation(kFALSE);
    649468
    650469    MCalibrateFact conv4;
    651470    conv4.SetScale(scale);
    652     conv4.SetCalibConst(calib);
     471    //conv4.SetCalibConst(calib);
    653472
    654473    MBadPixelsTreat treat4;
    … …  
    668487    tlist4.AddToList(&drsapply4);
    669488    tlist4.AddToList(&cont4);
     489    tlist4.AddToList(&filterdata4);
    670490    tlist4.AddToList(&extractor4);
    671 //    tlist4.AddToList(&fill4a);
    672491    tlist4.AddToList(&conv4);
    673492    tlist4.AddToList(&treat4);
    … …  
    679498    if (!loop4.GetDisplay())
    680499        return 16;
    681 */
     500
    682501    // ===================================================================
    683502
    … …  
    708527
    709528    MDrsCalibApply drsapply5;
    710     drsapply5.SetRemoveSpikes(spike_removal);
     529    //drsapply5.SetRemoveSpikes(4);
     530
     531    MTreatSaturation treatsat5;
     532
     533    MFilterData filterdata5;
    711534
    712535    MFDataPhrase filterdat("(MRawEvtHeader.GetTriggerID&0xff00)==0",     "SelectDat");
    … …  
    719542    // ---
    720543
    721     MTreatSaturation treatsat5;
    722     MFilterData filterdata5;
    723 
    724544    MExtractFACT extractor5dat;
    725545    extractor5dat.SetRange(first_slice, last_slice);
    726 /*
    727     MExtractTimeAndChargeSpline extractor5cal;
     546    extractor5dat.SetNoiseCalculation(kFALSE);
     547
     548    MExtractFACT extractor5cal;
    728549    extractor5cal.SetRange(first_slice, last_slice);
    729     extractor5cal.SetRiseTimeHiGain(rise_time_cal);
    730     extractor5cal.SetFallTimeHiGain(fall_time_cal);
    731     extractor5cal.SetHeightTm(heighttm);
    732     extractor5cal.SetChargeType(type);
    733     extractor5cal.SetSaturationLimit(600000);
    734550    extractor5cal.SetNoiseCalculation(kFALSE);
    735551
    736     MExtractTimeAndChargeSpline extractor5tm("ExtractTM");
     552    MExtractFACT extractor5tm("ExtractTM");
    737553    extractor5tm.SetRange(last_slice, 294);
    738     extractor5tm.SetRiseTimeHiGain(1);
    739     extractor5tm.SetFallTimeHiGain(1);
    740     extractor5tm.SetHeightTm(0.5);
    741     extractor5tm.SetChargeType(MExtralgoSpline::kAmplitudeRel);
    742     extractor5tm.SetSaturationLimit(600000);
    743554    extractor5tm.SetNoiseCalculation(kFALSE);
    744555    extractor5tm.SetNameSignalCam("TimeMarkerAmplitude");
    745556    extractor5tm.SetNameTimeCam("TimeMarkerTime");
    746 */
     557
    747558    extractor5dat.SetFilter(&filterncl);
    748     //extractor5cal.SetFilter(&filtercal);
     559    extractor5cal.SetFilter(&filtercal);
    749560    //extractor4tm.SetFilter(&filtercal);
    750561
    … …  
    752563    MCalibrateFact conv5;
    753564    conv5.SetScale(scale);
    754     conv5.SetCalibConst(calib);
     565    //conv5.SetCalibConst(calib);
    755566
    756567    MCalibrateDrsTimes calctm5;
    757568    calctm5.SetNameUncalibrated("UncalibratedTimes");
    758 /*
     569
    759570    MCalibrateDrsTimes calctm5tm("CalibrateTimeMarker");
    760571    calctm5tm.SetNameArrivalTime("TimeMarkerTime");
    … …  
    806617    //fill5w.SetDrawOption("gaus");
    807618    //fill5x.SetDrawOption("gaus");
    808 */
     619
    809620
    810621    MBadPixelsTreat treat5;
    811622    treat5.SetProcessPedestalRun(kFALSE);
    812623    treat5.SetProcessPedestalEvt(kFALSE);
    813 /*
     624
    814625    MHSectorVsTime hist5cal("CalVsTm");
    815626    MHSectorVsTime hist5ped("PedVsTm");
    … …  
    829640    fill5cal.SetFilter(&filtercal);
    830641    fill5ped.SetFilter(&filterped);
    831 */
     642
    832643    MHCamEvent evt5b(0, "ExtSig",   "Extracted signal;;S [mV·sl]");
    833644    MHCamEvent evt5c(0, "CalSig",   "Calibrated and interpolated signal;;S [~phe]");
    … …  
    860671    // The second rule is for the case reading raw-files!
    861672    MWriteRootFile write5(2, fname, "RECREATE", "Calibrated Data");
    862     write5.AddContainer("MRawRunHeader",       "RunHeaders");
    863     write5.AddContainer("MGeomCam",            "RunHeaders");
    864     write5.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE);
    865     write5.AddContainer("MCorsikaRunHeader",   "RunHeaders", kFALSE);
    866     write5.AddContainer("MMcRunHeader",        "RunHeaders", kFALSE);
    867 
    868     // Common events
    869     write5.AddContainer("MCorsikaEvtHeader",   "Events", kFALSE);
    870     write5.AddContainer("MMcEvt",              "Events", kFALSE);
    871     write5.AddContainer("IncidentAngle",       "Events", kFALSE);
    872     write5.AddContainer("MPointingPos",        "Events", kFALSE);
    873     write5.AddContainer("MSignalCam",          "Events");
    874     write5.AddContainer("MTime",               "Events", kFALSE);
    875     write5.AddContainer("MRawEvtHeader",       "Events");
     673    write5.AddContainer("MRawRunHeader", "RunHeaders");
     674    write5.AddContainer("MGeomCam",      "RunHeaders");
     675    write5.AddContainer("MSignalCam",    "Events");
     676    write5.AddContainer("MTime",         "Events");
     677    write5.AddContainer("MRawEvtHeader", "Events");
     678    //write.AddContainer("MTriggerPattern", "Events");
    876679
    877680    // ------------------ Setup histograms and fill tasks ----------------
    … …  
    879682    MContinue test;
    880683    test.SetFilter(&filterncl);
    881 /*
     684
    882685    MTaskList tlist5tm;
    883686    tlist5tm.AddToList(&extractor5tm);
    … …  
    894697    tlist5tm.AddToList(&fill5x);
    895698    tlist5tm.SetFilter(&filtercal);
    896 */
     699
    897700    MTaskList tlist5dat;
    898701    tlist5dat.AddToList(&fill5b);
    … …  
    914717    tlist5.AddToList(&filterdata5);
    915718    tlist5.AddToList(&extractor5dat);
    916     //tlist5.AddToList(&extractor5cal);
     719    tlist5.AddToList(&extractor5cal);
    917720    tlist5.AddToList(&calctm5);
    918     //tlist5.AddToList(&tlist5tm);
     721    tlist5.AddToList(&tlist5tm);
    919722    tlist5.AddToList(&conv5);
    920723    tlist5.AddToList(&treat5);
    921     //tlist5.AddToList(&fill5ped);
    922     //tlist5.AddToList(&fill5cal);
     724    tlist5.AddToList(&fill5ped);
     725    tlist5.AddToList(&fill5cal);
    923726    tlist5.AddToList(&tlist5dat);
    924727    tlist5.AddToList(&write5);
    925728
    926     if (!loop5.Eventloop(max4))
     729    if (!loop5.Eventloop(max5))
    927730        return 18;
    928731
    … …  
    945748    return 0;
    946749}
    947 
    948 int callisto(const ULong64_t seqnum, const char *outpath = "output")
    949 {
    950     UInt_t night = seqnum/1000;
    951     UInt_t num   = seqnum%1000;
    952 
    953     TString file = Form("/scratch/fact/sequences/%04d/%02d/%02d/%06d_%03d.seq",
    954                         night/10000, (night/100)%100, night%100, num);
    955 
    956     return callisto(file.Data(), outpath);
    957 }
Note: See TracChangeset for help on using the changeset viewer.