Changeset 9462 for trunk/MagicSoft


Ignore:
Timestamp:
06/20/09 10:14:33 (15 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
19 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/Changelog

    r9460 r9462  
    1818
    1919                                                 -*-*- END OF LINE -*-*-
     20
     21 2009/06/20 Thomas Bretz
     22
     23   * mbase/MLut.[h,cc]:
     24     - renamed Default to DefaultCol
     25     - added DefaultRow
     26     - added Print
     27     - added some const qualifiers
     28
     29   * mbase/MParList.cc:
     30     - make sure a parlist is not added to itself
     31     - improved some debug output
     32
     33   * melectronics/MAvalanchePhotoDiode.[h,cc]:
     34     - keep a time-stamp fTime
     35     - added functions to "relax/evolve" the apd over a given time
     36
     37   * mgeom/MGeomCamDwarf.cc:
     38     - fixed a broken comment
     39
     40   * mjobs/MJPedestal.cc:
     41     - improved output in case of error
     42
     43   * mjobs/MJSimulation.cc:
     44     - implemented a PreCut (to gain processing time)
     45     - implemented the correct type id in the file name
     46
     47   * mpedestal/MPedCalcPedRun.cc:
     48     - check for the run-number only if real data
     49
     50   * mraw/MRawRunHeader.[h,cc]:
     51     - added GetRuntypeChar
     52     - fixed a bug which returned the runtype instead of the
     53       telescopenumber
     54
     55   * mreport/MReport.cc:
     56     - added a new report version number 200905170 (preliminray!)
     57
     58   * mreport/MReportDrive.cc:
     59     - added a fix for the starguider sttaus in the report
     60
     61   * msimcamera/MSimAPD.cc:
     62     - instead of a full initialization of the APD before each event
     63       we now only simulate a time evolution which is needed to relax
     64       a cell to one permille error. This is much faster.
     65
     66   * msimcamera/MSimTrigger.[h,cc]:
     67     - allow switching off of the electronics trigger
     68     - replaced the real coincidence by better and faster algorithm
     69     - the new algorithm also supports multiplicity triggers
     70     - therefore added fMinMultiplicity
     71     - updated the handling and ouput in case of empty coincidence maps
     72
     73   * mtrigger/MTriggerPattern.h:
     74     - added some information in the bit description
     75
     76   * mtrigger/MTriggerPatternDecode.cc:
     77     - fixed the decoding for runs for which L1TPU was connected
     78       instead of L1
     79
     80
    2081
    2182 2009/06/18 Daniela Dorner
  • trunk/MagicSoft/Mars/mbase/MLut.cc

    r9318 r9462  
    7070//
    7171// Check if it is a default lut which would just map every entry to itself
    72 // An empty Lut is a default lut
    73 //
    74 Bool_t MLut::IsDefault() const
     72// An empty Lut is a default column
     73//
     74Bool_t MLut::IsDefaultCol() const
    7575{
    7676    if (IsEmpty())
     
    9292// Setup a default lut which just maps n-entris to themself
    9393//
    94 void MLut::SetDefault(UInt_t n)
     94void MLut::SetDefaultCol(UInt_t n)
    9595{
    9696    Delete();
     
    106106    fMinEntries = 1;
    107107    fMaxEntries = 1;
     108
     109    fMaxIndex = n;
     110}
     111
     112// --------------------------------------------------------------------------
     113//
     114// Check if it is a default lut which would just map all entries to one.
     115// An empty Lut is a default row
     116//
     117Bool_t MLut::IsDefaultRow() const
     118{
     119    if (IsEmpty())
     120        return kTRUE;
     121
     122    if (GetEntriesFast()!=1)
     123        return kFALSE;
     124
     125    const MArrayI &idx = GetRow(0);
     126
     127    // Loop over all rows
     128    for (UInt_t x=0; x<idx.GetSize(); x++)
     129        if (UInt_t(idx[x])!=x)
     130            return kFALSE;
     131
     132    return kTRUE;
     133}
     134
     135
     136// --------------------------------------------------------------------------
     137//
     138// Setup a default lut which maps all n-entris to one
     139//
     140void MLut::SetDefaultRow(UInt_t n)
     141{
     142    Delete();
     143
     144    MArrayI &idx = *new MArrayI(n);
     145
     146    for (UInt_t y=0; y<n; y++)
     147        idx[y] = y;
     148
     149    Add(&idx);
     150
     151    fMinEntries = n;
     152    fMaxEntries = n;
    108153
    109154    fMaxIndex = n;
     
    288333// Write a lut to a stream.
    289334//
    290 Int_t MLut::WriteStream(ostream &out)
     335Int_t MLut::WriteStream(ostream &out) const
    291336{
    292337    const Int_t n = GetEntriesFast();
     
    329374// Write a lut to a file
    330375//
    331 Int_t MLut::WriteFile(const char *fname)
     376Int_t MLut::WriteFile(const char *fname) const
    332377{
    333378    TString expname(fname);
     
    344389    return WriteStream(fout);
    345390}
     391
     392void MLut::Print(const Option_t *o) const
     393{
     394    gLog << all;
     395    WriteStream(gLog);
     396}
  • trunk/MagicSoft/Mars/mbase/MLut.h

    r9318 r9462  
    3636    Bool_t HasConstantLength() const { return fMinEntries==fMaxEntries; }
    3737    Bool_t IsEmpty() const { return fMaxEntries==0; }
    38     Bool_t IsDefault() const;
     38    Bool_t IsDefaultCol() const;
     39    Bool_t IsDefaultRow() const;
    3940
    4041    // MLut conversions
     
    4344
    4445    // Setter
    45     void SetDefault(UInt_t n);
     46    void SetDefaultCol(UInt_t n);
     47    void SetDefaultRow(UInt_t n);
    4648
    4749    // MLut I/O
    4850    Int_t ReadStream(istream &in);
    49     Int_t WriteStream(ostream &out);
     51    Int_t WriteStream(ostream &out) const;
    5052
    5153    Int_t ReadFile(const char *fname);
    52     Int_t WriteFile(const char *fname);
     54    Int_t WriteFile(const char *fname) const;
     55
     56    // TObject
     57    void Print(const Option_t *o="") const;
     58    void Print(const Option_t *o, const Option_t *o2) const { Print(o); }
    5359
    5460    ClassDef(MLut, 1) // A simple and fast easy-to-use look-up-table
  • trunk/MagicSoft/Mars/mbase/MParList.cc

    r9347 r9462  
    184184        return kFALSE;
    185185
     186    if (cont==this)
     187    {
     188        *fLog << err << dbginf << "Error: It is not allowed to add a parameter list to itself." << endl;
     189        return kFALSE;
     190    }
     191
    186192    //
    187193    // Get Name of new container
     
    223229        if (!fContainer->FindObject(where))
    224230        {
    225             *fLog << dbginf << "Error: Cannot find parameter container after which the new one should be added!" << endl;
     231            *fLog << err << dbginf << "Error: Cannot find parameter container after which the new one should be added!" << endl;
    226232            return kFALSE;
    227233        }
     
    230236    if (!cont->InheritsFrom(MParContainer::Class()))
    231237    {
    232         *fLog << dbginf << "Error: Cantainer MUST derive from MParContainer!" << endl;
     238        *fLog << err << dbginf << "Error: Cantainer MUST derive from MParContainer!" << endl;
    233239        return kFALSE;
    234240    }
  • trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.cc

    r9261 r9462  
    7777APD::APD(Int_t n, Float_t prob, Float_t dt, Float_t rt)
    7878    : fHist("APD", "", n, 0.5, n+0.5, n, 0.5, n+0.5),
    79     fCrosstalkProb(prob), fDeadTime(dt), fRecoveryTime(rt)
     79    fCrosstalkProb(prob), fDeadTime(dt), fRecoveryTime(rt), fTime(-1)
    8080{
    8181    fHist.SetDirectory(0);
     82}
     83
     84// --------------------------------------------------------------------------
     85//
     86// This is the time the cell needs after a hit until less than threshold
     87// (0.001 == one permille) of the signal is lost.
     88//
     89// If all cells of a G-APD have fired simultaneously and they are fired
     90// once more after the relaxation time (with the defaultthreshold of 0.001)
     91// the chip will show a signal which is one permille too small.
     92//
     93Float_t APD::GetRelaxationTime(Float_t threshold) const
     94{
     95    return fDeadTime-TMath::Log(threshold)*fRecoveryTime;
    8296}
    8397
     
    197211// The default time is 0.
    198212//
     213// If you want t w.r.t. fTime use HitRandomCellRelative istead.
     214//
    199215Float_t APD::HitRandomCell(Float_t t)
    200216{
     
    219235// If deadtime and recovery time are 0 then t-1 is set.
    220236//
     237// Sets fTime to t
     238//
    221239// The default time is 0.
    222240//
     
    231249
    232250    fHist.SetEntries(1);
     251
     252    fTime = t;
    233253}
    234254
     
    238258// by calling HitCellImp with time t-gRandom->Exp(n/rate) with n being
    239259// the total number of cells.
     260//
     261// Sets fTime to t
     262//
    240263// The default time is 0.
    241264//
     
    253276    for (int x=1; x<=nx; x++)
    254277        for (int y=1; y<=ny; y++)
    255         {
    256278            HitCellImp(x, y, t-MMath::RndmExp(f));
    257         }
     279
     280    fTime = t;
     281}
     282
     283// --------------------------------------------------------------------------
     284//
     285// Evolve the chip from fTime to fTime+dt if it with a given frequency
     286// freq. Returns the total signal "recorded".
     287//
     288// fTime is set to the fTime+dt
     289//
     290// If you want to evolve over a default relaxation time (relax the chip
     291// from a signal) use Relax instead.
     292//
     293Float_t APD::Evolve(Double_t freq, Double_t dt)
     294{
     295    const Double_t avglen = 1./freq;
     296
     297    const Double_t end  = fTime+dt;
     298
     299    Float_t hits = 0;
     300
     301    Double_t time = fTime;
     302    while (1)
     303    {
     304        time += MMath::RndmExp(avglen);
     305        if (time>end)
     306            break;
     307
     308        hits += HitRandomCell(time);
     309    }
     310
     311    fTime = end;
     312
     313    return hits;
    258314}
    259315
  • trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.h

    r9261 r9462  
    1515    Float_t fRecoveryTime;   // Recoverytime after Deadtime (1-exp(-t/fRecoveryTime)
    1616
     17    Float_t fTime;           // A user settable time of the system
     18
    1719    Float_t HitCellImp(Int_t x, Int_t y, Float_t t=0);
    1820
     
    2224    Float_t HitCell(Int_t x, Int_t y, Float_t t=0);
    2325    Float_t HitRandomCell(Float_t t=0);
     26    Float_t HitRandomCellRelative(Float_t t=0) { return HitRandomCell(fTime+t); }
    2427
    2528    void FillEmpty(Float_t t=0);
    2629    void FillRandom(Float_t rate, Float_t t=0);
     30
     31    void Init(Float_t rate) { if (fTime<0) FillRandom(rate); else Relax(rate); }
    2732
    2833    Int_t CountDeadCells(Float_t t=0) const;
     
    3540    Float_t GetDeadTime() const { return fDeadTime; }
    3641    Float_t GetRecoveryTime() const { return fRecoveryTime; }
     42    Float_t GetTime() const { return fTime; }
     43
     44    Float_t GetRelaxationTime(Float_t threshold=0.001) const;
     45
     46    Float_t GetLastHit() const { return fHist.GetMaximum(); }
     47
     48    void SetTime(Float_t tm) { fTime=tm; }
     49    void IncreaseTime(Float_t dt) { fTime += dt; }
     50
     51    Float_t Evolve(Double_t freq, Double_t dt);
     52    Float_t Relax(Double_t freq, Float_t threshold=0.001) { return Evolve(freq, GetRelaxationTime(threshold)); }
    3753
    3854    void Draw(Option_t *o="") { fHist.Draw(o); }
  • trunk/MagicSoft/Mars/mgeom/MGeomCamDwarf.cc

    r9385 r9462  
    3838// hexagonal camera, but the edges filled with additional pixels
    3939// inside a circle.
    40 //    MGeomCamDwarf cam(
     40//    MGeomCamDwarf cam(209.5, 13.2);
    4141//
    4242////////////////////////////////////////////////////////////////////////////
  • trunk/MagicSoft/Mars/mjobs/MJPedestal.cc

    r9317 r9462  
    11541154    if (taskenv.GetNumExecutions()<fMinEvents)
    11551155    {
    1156         *fLog << err << GetDescriptor() << ": Failed. Less than the required " << fMinEvents << " evts processed." << endl;
     1156        *fLog << err << GetDescriptor() << ": Failed. Less (" << taskenv.GetNumExecutions() << ") than the required " << fMinEvents << " pedestal events extracted." << endl;
    11571157        return kFALSE;
    11581158    }
     
    11601160    if (fIsPulsePosCheck && pulcam.GetNumEvents()<fMinCosmics)
    11611161    {
    1162         *fLog << err << GetDescriptor() << ": Failed. Less than the required " << fMinCosmics << " cosmics evts processed." << endl;
     1162        *fLog << err << GetDescriptor() << ": Failed. Less (" << pulcam.GetNumEvents() << ") than the required " << fMinCosmics << " cosmics evts processed." << endl;
    11631163        return kFALSE;
    11641164    }
     
    11661166    if (fPedestalCamOut.GetNumEvents()<fMinPedestals)
    11671167    {
    1168         *fLog << err << GetDescriptor() << ": Failed. Less than the required " << fMinPedestals << " pedetsl evts processed." << endl;
     1168        *fLog << err << GetDescriptor() << ": Failed. Less (" << fPedestalCamOut.GetNumEvents() << ") than the required " << fMinPedestals << " pedestal evts processed." << endl;
    11691169        return kFALSE;
    11701170    }
  • trunk/MagicSoft/Mars/mjobs/MJSimulation.cc

    r9461 r9462  
    551551    MSimCalibrationSignal simcal;
    552552    simcal.SetNameGeomCam("GeomCones");
     553
     554    // Sky Quality Meter: 21.82M = 2.02e-4cd/m^2
     555    // 1cd = 12.566 lm / 4pi sr
    553556
    554557    // FIXME: Simulate photons before cones and QE!
  • trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc

    r9157 r9462  
    11/* ======================================================================== *\
    2 ! $Name: not supported by cvs2svn $:$Id: MPedCalcPedRun.cc,v 1.54 2008-11-12 16:04:18 tbretz Exp $
     2! $Name: not supported by cvs2svn $:$Id: MPedCalcPedRun.cc,v 1.55 2009-06-20 09:14:33 tbretz Exp $
    33! --------------------------------------------------------------------------
    44!
     
    320320//
    321321Bool_t MPedCalcPedRun::IsPedBitSet()
    322 {
    323     if (fRunHeader->GetTelescopeNumber()==1 && fRunHeader->GetRunNumber()<gkFirstRunWithFinalBits)
     322{     
     323    if (!fRunHeader->IsMonteCarloRun() && fRunHeader->GetTelescopeNumber()==1 && fRunHeader->GetRunNumber()<gkFirstRunWithFinalBits)
    324324        return kFALSE;
    325325
  • trunk/MagicSoft/Mars/mraw/MRawRunHeader.cc

    r9433 r9462  
    940940// to print it as readable text.
    941941//
    942 const char *MRawRunHeader::GetRunTypeStr() const
     942const Char_t *MRawRunHeader::GetRunTypeStr() const
    943943{
    944944    switch (fRunType)
     
    960960    default:
    961961        return "<unknown>";
     962    }
     963}
     964
     965const Char_t MRawRunHeader::GetRunTypeChar() const
     966{
     967    switch (fRunType&0xff)
     968    {
     969    case kRTData:
     970    case kRTPointRun:
     971        return 'D';
     972    case kRTPedestal:
     973        return 'P';
     974    case kRTCalibration:
     975        return 'C';
     976    case kRTLinearity:
     977        return 'N';
     978    default:
     979        return ' ';
    962980    }
    963981}
  • trunk/MagicSoft/Mars/mraw/MRawRunHeader.h

    r9433 r9462  
    126126    UInt_t   GetFileID() const            { return fRunNumber>1000000?(fRunNumber%1000000)*1000+(fFileNumber%1000):fRunNumber; }
    127127    TString  GetStringID() const;
    128     UShort_t GetTelescopeNumber() const   { return fRunType; }
     128    UShort_t GetTelescopeNumber() const   { return fTelescopeNumber; }
    129129    UShort_t GetRunType() const           { return fRunType; }
    130130    const Char_t *GetRunTypeStr() const;
     131    const Char_t  GetRunTypeChar() const;
    131132    const Char_t *GetProjectName() const  { return fProjectName; }
    132133    const Char_t *GetSourceName() const   { return fSourceName; }
  • trunk/MagicSoft/Mars/mreport/MReport.cc

    r9423 r9462  
    169169//    200805190  | 54711.5 |         |  200809030
    170170//    200812040  | 54814.5 |         |  200812140
     171//    200902210  | 54968.5 |         |  200905170
    171172//
    172173Int_t MReport::Interprete(TString &str, const MTime &start, const MTime &stop, Int_t ver)
     
    208209        ver=200812140;
    209210
     211    if (ver==200902210 && GetMjd()>54968.5)
     212        ver=200905170;
     213
    210214    // Interprete body (contents) of report
    211215    const Int_t rc = InterpreteBody(str, ver);
  • trunk/MagicSoft/Mars/mreport/MReportDrive.cc

    r8885 r9462  
    135135    }
    136136
     137    if (ver>=200805170)
     138    {
     139        Int_t dummy; // Starguider switched on or not
     140        n=sscanf(str.Data(), "%d %n", &dummy, &len);
     141        if (n!=1)
     142        {
     143            *fLog << warn << "WARNING - Not enough arguments." << endl;
     144            return kCONTINUE;
     145        }
     146
     147        str.Remove(0, len);
     148        str = str.Strip(TString::kBoth);
     149    }
     150
    137151    return str.IsNull() ? kTRUE : kCONTINUE;
    138152}
  • trunk/MagicSoft/Mars/msimcamera/MSimAPD.cc

    r9427 r9462  
    206206        return kERROR;
    207207    }
    208 
     208/*
    209209    for (UInt_t idx=0; idx<npix; idx++)
    210210    {
     
    212212        static_cast<APD*>(fAPDs.UncheckedAt(idx))->FillRandom(freq, fStat->GetTimeFirst());
    213213    }
     214*/
     215
     216    // This tries to initialize dead and relaxing cells properly. If
     217    // the APD has not been initialized before the chip is randomsly
     218    // filled, otherwise a time window of the default relaxing time
     219    // is simulated, so that the previous influence is less than a permille.
     220    for (UInt_t idx=0; idx<npix; idx++)
     221    {
     222        const Double_t freq = fRates ? (*fRates)[idx].GetPedestal() : fFreq;
     223        static_cast<APD*>(fAPDs.UncheckedAt(idx))->Init(freq);
     224    }
    214225
    215226    // Get number of photons
     
    223234
    224235        // Get arrival time of photon and idx
    225         const Double_t t = ph.GetTime();
     236        const Double_t t = ph.GetTime()-fStat->GetTimeFirst();
    226237        const Int_t  idx = ph.GetTag();
    227238        if (idx<0)
     
    239250        // Simulate hitting the APD (the signal height in effective
    240251        // "number of photons" is returned)
    241         const Double_t hits = static_cast<APD*>(fAPDs.UncheckedAt(idx))->HitRandomCell(t);
     252        const Double_t hits = static_cast<APD*>(fAPDs.UncheckedAt(idx))->HitRandomCellRelative(t);
    242253
    243254        // FIXME: Make a proper simulation of the excess noise!!!
     
    248259    }
    249260
     261    // Now we have to shift the evolved time of all APDs to the end of our
     262    // simulated time.
     263    for (UInt_t idx=0; idx<npix; idx++)
     264        static_cast<APD*>(fAPDs.UncheckedAt(idx))->IncreaseTime(fStat->GetTimeLast());
     265
    250266    return kTRUE;
    251267}
  • trunk/MagicSoft/Mars/msimcamera/MSimTrigger.cc

    r9318 r9462  
    4242//
    4343// With a second look-up table fCoincidenceMap the analog channels are
    44 // checked for coincidences. The coincidense must at least be of the length
     44// checked for coincidences. The coincidence must at least be of the length
    4545// defined by fCoincidenceTime. The earliest coincide is then stored as
    4646// trigger position.
     47//
     48// If a minimum multiplicity m is given, m signals above threshold
     49// in the coincidence patterns are enough to emit a trigger signal.
    4750//
    4851//
     
    97100    fEvtHeader(0), fElectronicNoise(0), fGain(0),
    98101    fDiscriminatorThreshold(-1), fDigitalSignalLength(8), fCoincidenceTime(0.5),
    99     fShiftBaseline(kTRUE), fUngainSignal(kTRUE)
     102    fShiftBaseline(kTRUE), fUngainSignal(kTRUE), fSimulateElectronics(kTRUE),
     103    fMinMultiplicity(-1)
    100104{
    101105    fName  = name  ? name  : "MSimTrigger";
     
    107111// Take two TObjArrays with a collection of digital signals.
    108112// Every signal from one array is compared with any from the other array.
    109 // For all signals whihc overlaps and which have an overlap time >gate
     113// For all signals which overlap and which have an overlap time >gate
    110114// a new digital signal is created storing start time and length of overlap.
    111115// They are collected in a newly allocated TObjArray. A pointer to this array
    112116// is returned.
    113117//
    114 // Th euser gains owenership of the object, ie.e., the user is responsible of
     118// The user gains owenership of the object, i.e., the user is responsible of
    115119// deleting the memory.
    116120//
    117 TObjArray *MSimTrigger::CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2, Float_t gate) const
     121TObjArray *MSimTrigger::CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2/*, Float_t gate*/) const
    118122{
    119123    TObjArray *res = new TObjArray;
     
    131135        {
    132136            MDigitalSignal *ttl = new MDigitalSignal(*ttl1, *ttl2);
    133 
     137            /*
    134138            if (ttl->GetLength()<=gate)
    135139            {
     
    137141                continue;
    138142            }
    139 
     143            */
    140144            res->Add(ttl);
    141145        }
     
    143147
    144148    res->SetOwner();
     149
     150    return res;
     151}
     152
     153class Edge : public TObject
     154{
     155private:
     156    Double_t fEdge;
     157    Int_t    fRising;
     158
     159public:
     160    Edge(Double_t t, Int_t rising) : fEdge(t), fRising(rising) { }
     161    Bool_t IsSortable() const { return kTRUE; }
     162    Int_t Compare(const TObject *o) const { const Edge *e = static_cast<const Edge*>(o); if (e->fEdge<fEdge) return 1; if (e->fEdge>fEdge) return -1; return 0; }
     163
     164    Int_t IsRising() const { return fRising; }
     165    Double_t GetEdge() const { return fEdge; }
     166};
     167
     168// --------------------------------------------------------------------------
     169//
     170// Calculate a multiplicity trigger on the given array(s). The idx-array
     171// conatins all channels which should be checked for coincidences
     172// and the ttls array conatins the arrays with the digital signals.
     173//
     174// For the windows in which more or euqal than threshold channels have
     175// a high signal a new MDigitalSignal is created.  newly allocated
     176// array with a collection of these trigger signals is returned.
     177//
     178TObjArray *MSimTrigger::CalcMinMultiplicity(const MArrayI &idx, const TObjArray &ttls, Int_t threshold) const
     179{
     180    // Create a new array for the rising and falling edges of the signals
     181    TObjArray times;
     182    times.SetOwner();
     183
     184    // Fill the array with edges from all digital signals of all our channels
     185    for (UInt_t k=0; k<idx.GetSize(); k++)
     186    {
     187        TObjArray *arr = static_cast<TObjArray*>(ttls[idx[k]]);
     188
     189        TIter Next(arr);
     190        MDigitalSignal *ttl = 0;
     191        while ((ttl=static_cast<MDigitalSignal*>(Next())))
     192        {
     193            times.Add(new Edge(ttl->GetStart(),  1));
     194            times.Add(new Edge(ttl->GetEnd(),   -1));
     195        }
     196    }
     197
     198    // Sort them in time
     199    times.Sort();
     200
     201    // Start with no channel active
     202    Int_t lvl = 0;
     203
     204    TObjArray *res = new TObjArray;
     205    res->SetOwner();
     206
     207    // First remove all edges which do not change the status
     208    // "below threshold" or "above threshold"
     209    for (int i=0; i<times.GetEntriesFast(); i++)
     210    {
     211        // Get i-th edge
     212        const Edge &e = *static_cast<Edge*>(times.UncheckedAt(i));
     213
     214        // Claculate what the number of active channels after the edge is
     215        const Int_t lvl1 = lvl + e.IsRising();
     216
     217        // Remove edge if number of active channels before and after the
     218        // edge lower is lower than the threshold or higher than
     219        // the threshold
     220        if (lvl+1<threshold || lvl-1>=threshold)
     221            delete times.RemoveAt(i);
     222
     223        // keep the (now) "previous" level
     224        lvl = lvl1<0 ? 0 : lvl1;
     225    }
     226
     227    // Remove the empty slots from the array
     228    times.Compress();
     229
     230    //
     231    for (int i=0; i<times.GetEntriesFast()-1; i++)
     232    {
     233        // get the current edge
     234        const Edge &e0 = *static_cast<Edge*>(times.UncheckedAt(i));
     235
     236        // go ahead if this is a falling edge
     237        if (e0.IsRising()!=1)
     238            continue;
     239
     240        // get the following edge (must be a falling edge now)
     241        const Edge &e1 = *static_cast<Edge*>(times.UncheckedAt(i+1));
     242
     243        // calculate the length of the digital signal
     244        const Double_t len = e1.GetEdge()-e0.GetEdge();
     245
     246        // Create a digital trigger signal
     247        MDigitalSignal *ds = new MDigitalSignal(e0.GetEdge(), len);
     248        //ds->SetIndex(lvl);
     249        res->Add(ds);
     250    }
    145251
    146252    return res;
     
    231337    }
    232338
    233     if (fElectronicNoise && !fRouteAC.IsEmpty() && !fRouteAC.IsDefault())
     339    if (fElectronicNoise && !fRouteAC.IsEmpty() && !fRouteAC.IsDefaultCol())
    234340    {
    235341        // FIXME: Apply to analog channels when summing
     
    238344    }
    239345
    240     if (fGain && !fRouteAC.IsEmpty() && !fRouteAC.IsDefault())
     346    if (fGain && !fRouteAC.IsEmpty() && !fRouteAC.IsDefaultCol())
    241347    {
    242348        // FIXME: Apply to analog channels when summing
     
    255361        *fLog << "Using " << fNameRouteAC << " for re-routing/summing of analog channels before discriminator." << endl;
    256362
    257     if (fCoincidenceMap.IsEmpty())
     363    if (fCoincidenceMap.IsEmpty() && fMinMultiplicity==0)
    258364        *fLog << "No coincidences of digital channels will be checked. Signal-above-threshold trigger applied." << endl;
    259365    else
    260         *fLog << "Using " << fNameCoincidenceMap << " to check for coincidences of the digital channels." << endl;
     366    {
     367        *fLog << "Using ";
     368        if (fCoincidenceMap.IsEmpty())
     369            *fLog << "the whole camera";
     370        else
     371            *fLog << "patterns from " << fNameCoincidenceMap;
     372        *fLog << " to check for ";
     373        if (fMinMultiplicity==0)
     374            *fLog << "coincidences of the digital channels." << endl;
     375        else
     376            *fLog << fMinMultiplicity << " multiplicity." << endl;
     377        }
    261378
    262379    *fLog << "Using discriminator threshold of " << fDiscriminatorThreshold << endl;
     
    265382}
    266383
     384/*
     385class MDigitalChannel : public TObjArray
     386{
     387private:
     388    TObjArray fArray;
     389
     390public:
     391    MDigitalSignal *GetSignal(UInt_t i) { return static_cast<MDigitalSignal*>(fArray[i]); }
     392
     393};
     394*/
     395
     396#include "MCamEvent.h"
     397class MTriggerSignal : public MParContainer, public MCamEvent
     398{
     399private:
     400    TObjArray fSignals;
     401
     402public:
     403    MTriggerSignal() { fSignals.SetOwner(); }
     404
     405    void Add(MDigitalSignal *signal) { fSignals.Add(signal); }
     406
     407    MDigitalSignal *GetSignal(UInt_t i) { return static_cast<MDigitalSignal*>(fSignals[i]); }
     408
     409    void Sort() { fSignals.Sort(); }
     410
     411    Int_t GetNumSignals() const { return fSignals.GetEntriesFast(); }
     412
     413    Float_t GetFirstTrigger() const
     414    {
     415        MDigitalSignal *sig = static_cast<MDigitalSignal*>(fSignals[0]);
     416        return sig ? sig->GetStart() : -1;
     417    }
     418    Bool_t GetPixelContent(Double_t&, Int_t, const MGeomCam&, Int_t) const
     419    {
     420        switch (1)
     421        {
     422        case 1: // yes/no
     423        case 2: // first time
     424        case 3: // length
     425        case 4: // n
     426            break;
     427        }
     428
     429        return kTRUE;
     430    }
     431    void DrawPixelContent(Int_t) const
     432    {
     433    }
     434};
     435
     436
     437void MSimTrigger::SetTrigger(Double_t pos)
     438{
     439    // FIXME: Jitter! (Own class?)
     440    fTrigger->SetVal(pos);
     441    fTrigger->SetReadyToSave();
     442
     443    // Trigger pattern to be set
     444    // FIXME: Make flexible
     445    const UInt_t pat = (UInt_t)~(MTriggerPattern::kTriggerLvl1 | (MTriggerPattern::kTriggerLvl1<<8));
     446
     447    // Flag this event as triggered by the lvl1 trigger (FIXME?)
     448    fEvtHeader->SetTriggerPattern(pos<0 ? 0 : pat);
     449    fEvtHeader->SetReadyToSave();
     450}
     451
    267452// --------------------------------------------------------------------------
    268453//
     
    270455{
    271456    // Invalidate trigger
    272     fTrigger->SetVal(-1);
     457    //fTrigger->SetVal(-1);
     458    // Calculate the minimum and maximum time for a valid trigger
     459    const Double_t freq    = fRunHeader->GetFreqSampling()/1000.;
     460    const Float_t  nsamp   = fRunHeader->GetNumSamplesHiGain();
     461    const Float_t  pulspos = fPulsePos->GetVal()/freq;
     462
     463    // Valid range in units of bins
     464    const Float_t min = fCamera->GetValidRangeMin()+pulspos;
     465    const Float_t max = fCamera->GetValidRangeMax()-(nsamp-pulspos);
     466
     467    if (!fSimulateElectronics)
     468    {
     469        SetTrigger(min);
     470        return kTRUE;
     471    }
    273472
    274473    // ================== Simulate channel bundling ====================
     
    297496            {
    298497                const UInt_t idx = row[j];
    299                 (*patches)[i].AddSignal((*fCamera)[idx]);
     498
     499                // FIXME: Shrinking the mapping table earlier (e.g.
     500                //        ReInit) would avoid a lot of if's
     501                if (idx<fCamera->GetNumChannels())
     502                    (*patches)[i].AddSignal((*fCamera)[idx]);
    300503            }
    301504        }
     
    329532    //        omitted in the loop
    330533    if (fCoincidenceMap.IsEmpty())
    331         fCoincidenceMap.SetDefault(npatch);
    332 
    333     // Calculate the minimum and maximum time for a valid trigger
    334     const Double_t freq    = fRunHeader->GetFreqSampling()/1000.;
    335     const Float_t  nsamp   = fRunHeader->GetNumSamplesHiGain();
    336     const Float_t  pulspos = fPulsePos->GetVal()/freq;
    337 
    338     // Valid range in units of bins
    339     const Float_t min = fCamera->GetValidRangeMin()+pulspos;
    340     const Float_t max = fCamera->GetValidRangeMax()-(nsamp-pulspos);
     534    {
     535        if (fMinMultiplicity>0)
     536            fCoincidenceMap.SetDefaultRow(npatch);
     537        else
     538            fCoincidenceMap.SetDefaultCol(npatch);
     539    }
    341540
    342541    // Create an array for the individual triggers
    343     TObjArray triggers;
    344     triggers.SetOwner();
     542    MTriggerSignal triggers;
    345543
    346544    Int_t cnt  = 0;
     
    352550        const MArrayI &idx = fCoincidenceMap.GetRow(j);
    353551
    354         // Start with a copy of the first coincidence channel
    355         TObjArray *arr = static_cast<TObjArray*>(ttls[idx[0]]->Clone());
    356         arr->SetOwner();
    357 
    358         // compare to all other channels in this coincidence patch, one by one
    359         for (UInt_t k=1; k<idx.GetSize() && arr->GetEntriesFast()>0; k++)
    360         {
    361             TObjArray *res = CalcCoincidence(*arr, *static_cast<TObjArray*>(ttls[idx[k]]),
    362                                              fCoincidenceTime);
    363 
    364             // Delete the original array and keep the new one
    365             delete arr;
    366             arr = res;
    367         }
     552        TObjArray *arr = 0;
     553
     554        if (fMinMultiplicity>0)
     555        {
     556            arr = CalcMinMultiplicity(idx, ttls, fMinMultiplicity);
     557        }
     558        else
     559        {
     560            arr = CalcMinMultiplicity(idx, ttls, idx.GetSize());
     561            /*
     562            // Start with a copy of the first coincidence channel
     563            arr = static_cast<TObjArray*>(ttls[idx[0]]->Clone());
     564            arr->SetOwner();
     565
     566            // compare to all other channels in this coincidence patch, one by one
     567            for (UInt_t k=1; k<idx.GetSize() && arr->GetEntriesFast()>0; k++)
     568            {
     569                TObjArray *res = CalcCoincidence(*arr, *static_cast<TObjArray*>(ttls[idx[k]]));//, fCoincidenceTime);
     570
     571                // Delete the original array and keep the new one
     572                delete arr;
     573                arr = res;
     574            }*/
     575        }
     576
     577        // Count the number of totally emitted coincidence signals
     578        cnt += arr->GetEntriesFast();
    368579
    369580        // Remove all signals which are not in the valid digitization range
    370581        // (This is not the digitization window, but the region in which
    371582        //  the analog channels contain usefull data)
     583        // and which are shorter than the defined coincidence gate.
    372584        TIter Next(arr);
    373585        MDigitalSignal *ttl = 0;
    374586        while ((ttl=static_cast<MDigitalSignal*>(Next())))
    375587        {
     588            if (ttl->GetLength()<fCoincidenceTime)
     589            {
     590                delete arr->Remove(ttl);
     591                continue;
     592            }
     593
    376594            if (ttl->GetStart()<min)
    377595            {
    378596                delete arr->Remove(ttl);
    379597                rmlo++;
     598                continue;
    380599            }
    381600            if (ttl->GetStart()>max)
     
    383602                delete arr->Remove(ttl);
    384603                rmhi++;
     604                continue;
    385605            }
     606
     607            // Set trigger channel index
     608            ttl->SetIndex(j);
    386609        }
    387610
     
    389612        arr->Compress();
    390613
    391         cnt += arr->GetEntriesFast();
    392 
    393614        // If we have at least one trigger keep the earliest one.
    394         // FIXME: The triggers should be ordered in time automatically: To be checked!
     615        // FIXME: The triggers might be ordered in time automatically:
     616        //        To be checked!
    395617        // FIXME: Simulate trigger dead-time!
    396618        if (arr->GetEntriesFast()>0)
    397             triggers.Add(arr->RemoveAt(0));
     619            triggers.Add(static_cast<MDigitalSignal*>(arr->RemoveAt(0)));
    398620
    399621        // delete the allocated space
    400622        delete arr;
    401     }
    402 
    403     // No trigger issued. Go on.
    404     if (triggers.GetEntriesFast()==0)
    405     {
    406         if (rmlo>0 || rmhi>0)
    407             *fLog << inf2 << rmlo << "/" << rmhi << " trigger out of valid range. No trigger raised." << endl;
    408         return kTRUE;
    409623    }
    410624
     
    414628    //         is omitted
    415629    triggers.Sort();
    416 
    417     // FIXME: Jitter! (Own class?)
    418     fTrigger->SetVal(static_cast<MDigitalSignal*>(triggers[0])->GetStart());
    419     fTrigger->SetReadyToSave();
    420 
    421     // Flag this event as triggered by the lvl1 trigger (FIXME?)
    422     fEvtHeader->SetTriggerPattern((UInt_t)~(MTriggerPattern::kTriggerLvl1 | (MTriggerPattern::kTriggerLvl1<<8)));
    423     fEvtHeader->SetReadyToSave();
    424 
    425     // inf2?
    426     *fLog << inf;
    427     *fLog << cnt << " triggers left in " << triggers.GetEntriesFast() << " patches (" << rmlo << "/" << rmhi << " trigger out of valid range), T=" << fTrigger->GetVal();
     630    // FIXME: Store triggers! (+ Reversed pixels?)
     631
     632    SetTrigger(triggers.GetFirstTrigger());
     633
     634    // No trigger issued. Go on.
     635    if (triggers.GetNumSignals()==0)
     636    {
     637        if (rmlo>0 || rmhi>0)
     638            *fLog << inf2 << GetNumExecutions() << ": " << rmlo << "/" << rmhi << " trigger out of valid range. No trigger raised." << endl;
     639        return kTRUE;
     640    }
     641
     642    // Number of patches which have triggered out of the total number of
     643    // Coincidence signals emitted. (If the total number is higher than
     644    // the number of triggers either some triggers had to be removed or
     645    // or a patch has emitted more than one trigger signal)
     646    // FIXME: inf2?
     647    *fLog << inf << GetNumExecutions() << ": ";
     648    *fLog << setw(3) << triggers.GetNumSignals() << " triggers left out of ";
     649    *fLog << setw(3) << cnt << " (" << rmlo << "/" << rmhi << " trigger out of valid range), T=" << fTrigger->GetVal();
    428650    *fLog << endl;
     651
     652    //# Trigger characteristics: gate length (ns), min. overlapping time (ns),
     653    //# amplitude and FWHM of (gaussian) single phe response for trigger:
     654    //trigger_prop 3.0 0.25 1.0 2.0
    429655
    430656    return kTRUE;
     
    438664// DigitalSignalLength:     8
    439665// CoincidenceTime:         0.5
     666// SimulateElectronics:     Yes
    440667//
    441668Int_t MSimTrigger::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
     
    472699    }
    473700
     701    if (IsEnvDefined(env, prefix, "SimulateElectronics", print))
     702    {
     703        rc = kTRUE;
     704        fSimulateElectronics = GetEnvValue(env, prefix, "SimulateElectronics", fSimulateElectronics);
     705    }
     706
     707    if (IsEnvDefined(env, prefix, "MinMultiplicity", print))
     708    {
     709        rc = kTRUE;
     710        fMinMultiplicity = GetEnvValue(env, prefix, "MinMultiplicity", fMinMultiplicity);
     711    }
     712
    474713    return rc;
    475714}
  • trunk/MagicSoft/Mars/msimcamera/MSimTrigger.h

    r9318 r9462  
    4040    Bool_t  fShiftBaseline;             // Shift the baseline back to 0 for the threshold (needs ElectronicNoise [MPedestalCam])
    4141    Bool_t  fUngainSignal;              // "Remove" the gain from the signal (needs Gain [MPedestalCam])
     42    Bool_t  fSimulateElectronics;       // If the electronics is not simulated the trigger is set artificially to the first photon arrived
     43
     44    Int_t   fMinMultiplicity;           // N out of M
    4245
    4346    // MSimTrigger
    44     TObjArray *CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2, Float_t gate=0) const;
     47    TObjArray *CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2/*, Float_t gate=0*/) const;
     48    TObjArray *CalcMinMultiplicity(const MArrayI &idx, const TObjArray &ttls, Int_t threshold) const;
     49    TObjArray *CalcCoincidences(const MArrayI &idx, const TObjArray &ttls) const;
     50    void SetTrigger(Double_t pos);
    4551
    4652    // MTask
  • trunk/MagicSoft/Mars/mtrigger/MTriggerPattern.h

    r9326 r9462  
    1212public:
    1313    enum Pattern_t {
    14         kTriggerLvl1 = BIT(0), //  1:
    15         kCalibration = BIT(1), //  2: Pulse Trigger
    16         kTriggerLvl2 = BIT(2), //  4: LUT Pseudo Size selection
    17         kPedestal    = BIT(3), //  8:
    18         kPinDiode    = BIT(4), // 16:
    19         kSumTrigger  = BIT(5)  // 32: Flag for an event taken with sum trigger
    20         //kUnused      = BIT(6)
    21         //kUnused      = BIT(7)
     14//        kUndefined1  = BIT(0), //  1   1: Level 1 from L2 board
     15        kTriggerLvl1 = BIT(0), //  1   1: Level 1 from L2 board
     16        kCalibration = BIT(1), //  2   2: Pulse Trigger
     17        kTriggerLvl2 = BIT(2), //  4   4: LUT Pseudo Size selection
     18        kPedestal    = BIT(3), //  8   8: Artificial pedestal event
     19        kPinDiode    = BIT(4), // 10  16:
     20        kSumTrigger  = BIT(5), // 20  32: Flag for an event taken with sum trigger
     21//        kTriggerLvl1 = BIT(6), // 40  64: Trigger lvl1 directly from L1 without going through L2
     22        kUndefined1  = BIT(6), // 40  64: Trigger lvl1 directly from L1 without going through L2
     23        kUndefined2  = BIT(7)  // 80 128: Undefined?
    2224    };
    2325
  • trunk/MagicSoft/Mars/mtrigger/MTriggerPatternDecode.cc

    r7170 r9462  
    106106        return kTRUE;
    107107
    108     UInt_t pattern = ~fEvtHeader->GetTriggerID();
     108    const UInt_t pattern = ~fEvtHeader->GetTriggerID();
    109109
    110110    // The the trigger pattern is currently written with inverted bits,
     
    133133    fPattern->fUnprescaled = (pattern>>8) & 0xff;
    134134
     135    // This is a workaround for the new scheme in which L1TPU (the signal
     136    // comming directly from the L1 is connected, but the L1 (routed
     137    // over L2 is disconnected)
     138    if (!fRunHeader->IsMonteCarloRun() && fRunHeader->GetTelescopeNumber()==1 &&
     139        fRunHeader->GetRunNumber()>1006246)
     140    {
     141        fPattern->fPrescaled   |= (pattern>> 6)&1;
     142        fPattern->fUnprescaled |= (pattern>>14)&1;
     143    }
     144
    135145    return kTRUE;
    136146}
Note: See TracChangeset for help on using the changeset viewer.