Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 9461) +++ trunk/MagicSoft/Mars/Changelog (revision 9462) @@ -18,4 +18,65 @@ -*-*- END OF LINE -*-*- + + 2009/06/20 Thomas Bretz + + * mbase/MLut.[h,cc]: + - renamed Default to DefaultCol + - added DefaultRow + - added Print + - added some const qualifiers + + * mbase/MParList.cc: + - make sure a parlist is not added to itself + - improved some debug output + + * melectronics/MAvalanchePhotoDiode.[h,cc]: + - keep a time-stamp fTime + - added functions to "relax/evolve" the apd over a given time + + * mgeom/MGeomCamDwarf.cc: + - fixed a broken comment + + * mjobs/MJPedestal.cc: + - improved output in case of error + + * mjobs/MJSimulation.cc: + - implemented a PreCut (to gain processing time) + - implemented the correct type id in the file name + + * mpedestal/MPedCalcPedRun.cc: + - check for the run-number only if real data + + * mraw/MRawRunHeader.[h,cc]: + - added GetRuntypeChar + - fixed a bug which returned the runtype instead of the + telescopenumber + + * mreport/MReport.cc: + - added a new report version number 200905170 (preliminray!) + + * mreport/MReportDrive.cc: + - added a fix for the starguider sttaus in the report + + * msimcamera/MSimAPD.cc: + - instead of a full initialization of the APD before each event + we now only simulate a time evolution which is needed to relax + a cell to one permille error. This is much faster. + + * msimcamera/MSimTrigger.[h,cc]: + - allow switching off of the electronics trigger + - replaced the real coincidence by better and faster algorithm + - the new algorithm also supports multiplicity triggers + - therefore added fMinMultiplicity + - updated the handling and ouput in case of empty coincidence maps + + * mtrigger/MTriggerPattern.h: + - added some information in the bit description + + * mtrigger/MTriggerPatternDecode.cc: + - fixed the decoding for runs for which L1TPU was connected + instead of L1 + + 2009/06/18 Daniela Dorner Index: trunk/MagicSoft/Mars/mbase/MLut.cc =================================================================== --- trunk/MagicSoft/Mars/mbase/MLut.cc (revision 9461) +++ trunk/MagicSoft/Mars/mbase/MLut.cc (revision 9462) @@ -70,7 +70,7 @@ // // Check if it is a default lut which would just map every entry to itself -// An empty Lut is a default lut -// -Bool_t MLut::IsDefault() const +// An empty Lut is a default column +// +Bool_t MLut::IsDefaultCol() const { if (IsEmpty()) @@ -92,5 +92,5 @@ // Setup a default lut which just maps n-entris to themself // -void MLut::SetDefault(UInt_t n) +void MLut::SetDefaultCol(UInt_t n) { Delete(); @@ -106,4 +106,49 @@ fMinEntries = 1; fMaxEntries = 1; + + fMaxIndex = n; +} + +// -------------------------------------------------------------------------- +// +// Check if it is a default lut which would just map all entries to one. +// An empty Lut is a default row +// +Bool_t MLut::IsDefaultRow() const +{ + if (IsEmpty()) + return kTRUE; + + if (GetEntriesFast()!=1) + return kFALSE; + + const MArrayI &idx = GetRow(0); + + // Loop over all rows + for (UInt_t x=0; xFindObject(where)) { - *fLog << dbginf << "Error: Cannot find parameter container after which the new one should be added!" << endl; + *fLog << err << dbginf << "Error: Cannot find parameter container after which the new one should be added!" << endl; return kFALSE; } @@ -230,5 +236,5 @@ if (!cont->InheritsFrom(MParContainer::Class())) { - *fLog << dbginf << "Error: Cantainer MUST derive from MParContainer!" << endl; + *fLog << err << dbginf << "Error: Cantainer MUST derive from MParContainer!" << endl; return kFALSE; } Index: trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.cc =================================================================== --- trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.cc (revision 9461) +++ trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.cc (revision 9462) @@ -77,7 +77,21 @@ APD::APD(Int_t n, Float_t prob, Float_t dt, Float_t rt) : fHist("APD", "", n, 0.5, n+0.5, n, 0.5, n+0.5), - fCrosstalkProb(prob), fDeadTime(dt), fRecoveryTime(rt) + fCrosstalkProb(prob), fDeadTime(dt), fRecoveryTime(rt), fTime(-1) { fHist.SetDirectory(0); +} + +// -------------------------------------------------------------------------- +// +// This is the time the cell needs after a hit until less than threshold +// (0.001 == one permille) of the signal is lost. +// +// If all cells of a G-APD have fired simultaneously and they are fired +// once more after the relaxation time (with the defaultthreshold of 0.001) +// the chip will show a signal which is one permille too small. +// +Float_t APD::GetRelaxationTime(Float_t threshold) const +{ + return fDeadTime-TMath::Log(threshold)*fRecoveryTime; } @@ -197,4 +211,6 @@ // The default time is 0. // +// If you want t w.r.t. fTime use HitRandomCellRelative istead. +// Float_t APD::HitRandomCell(Float_t t) { @@ -219,4 +235,6 @@ // If deadtime and recovery time are 0 then t-1 is set. // +// Sets fTime to t +// // The default time is 0. // @@ -231,4 +249,6 @@ fHist.SetEntries(1); + + fTime = t; } @@ -238,4 +258,7 @@ // by calling HitCellImp with time t-gRandom->Exp(n/rate) with n being // the total number of cells. +// +// Sets fTime to t +// // The default time is 0. // @@ -253,7 +276,40 @@ for (int x=1; x<=nx; x++) for (int y=1; y<=ny; y++) - { HitCellImp(x, y, t-MMath::RndmExp(f)); - } + + fTime = t; +} + +// -------------------------------------------------------------------------- +// +// Evolve the chip from fTime to fTime+dt if it with a given frequency +// freq. Returns the total signal "recorded". +// +// fTime is set to the fTime+dt +// +// If you want to evolve over a default relaxation time (relax the chip +// from a signal) use Relax instead. +// +Float_t APD::Evolve(Double_t freq, Double_t dt) +{ + const Double_t avglen = 1./freq; + + const Double_t end = fTime+dt; + + Float_t hits = 0; + + Double_t time = fTime; + while (1) + { + time += MMath::RndmExp(avglen); + if (time>end) + break; + + hits += HitRandomCell(time); + } + + fTime = end; + + return hits; } Index: trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.h =================================================================== --- trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.h (revision 9461) +++ trunk/MagicSoft/Mars/melectronics/MAvalanchePhotoDiode.h (revision 9462) @@ -15,4 +15,6 @@ Float_t fRecoveryTime; // Recoverytime after Deadtime (1-exp(-t/fRecoveryTime) + Float_t fTime; // A user settable time of the system + Float_t HitCellImp(Int_t x, Int_t y, Float_t t=0); @@ -22,7 +24,10 @@ Float_t HitCell(Int_t x, Int_t y, Float_t t=0); Float_t HitRandomCell(Float_t t=0); + Float_t HitRandomCellRelative(Float_t t=0) { return HitRandomCell(fTime+t); } void FillEmpty(Float_t t=0); void FillRandom(Float_t rate, Float_t t=0); + + void Init(Float_t rate) { if (fTime<0) FillRandom(rate); else Relax(rate); } Int_t CountDeadCells(Float_t t=0) const; @@ -35,4 +40,15 @@ Float_t GetDeadTime() const { return fDeadTime; } Float_t GetRecoveryTime() const { return fRecoveryTime; } + Float_t GetTime() const { return fTime; } + + Float_t GetRelaxationTime(Float_t threshold=0.001) const; + + Float_t GetLastHit() const { return fHist.GetMaximum(); } + + void SetTime(Float_t tm) { fTime=tm; } + void IncreaseTime(Float_t dt) { fTime += dt; } + + Float_t Evolve(Double_t freq, Double_t dt); + Float_t Relax(Double_t freq, Float_t threshold=0.001) { return Evolve(freq, GetRelaxationTime(threshold)); } void Draw(Option_t *o="") { fHist.Draw(o); } Index: trunk/MagicSoft/Mars/mgeom/MGeomCamDwarf.cc =================================================================== --- trunk/MagicSoft/Mars/mgeom/MGeomCamDwarf.cc (revision 9461) +++ trunk/MagicSoft/Mars/mgeom/MGeomCamDwarf.cc (revision 9462) @@ -38,5 +38,5 @@ // hexagonal camera, but the edges filled with additional pixels // inside a circle. -// MGeomCamDwarf cam( +// MGeomCamDwarf cam(209.5, 13.2); // //////////////////////////////////////////////////////////////////////////// Index: trunk/MagicSoft/Mars/mjobs/MJPedestal.cc =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJPedestal.cc (revision 9461) +++ trunk/MagicSoft/Mars/mjobs/MJPedestal.cc (revision 9462) @@ -1154,5 +1154,5 @@ if (taskenv.GetNumExecutions()GetTelescopeNumber()==1 && fRunHeader->GetRunNumber()IsMonteCarloRun() && fRunHeader->GetTelescopeNumber()==1 && fRunHeader->GetRunNumber()"; + } +} + +const Char_t MRawRunHeader::GetRunTypeChar() const +{ + switch (fRunType&0xff) + { + case kRTData: + case kRTPointRun: + return 'D'; + case kRTPedestal: + return 'P'; + case kRTCalibration: + return 'C'; + case kRTLinearity: + return 'N'; + default: + return ' '; } } Index: trunk/MagicSoft/Mars/mraw/MRawRunHeader.h =================================================================== --- trunk/MagicSoft/Mars/mraw/MRawRunHeader.h (revision 9461) +++ trunk/MagicSoft/Mars/mraw/MRawRunHeader.h (revision 9462) @@ -126,7 +126,8 @@ UInt_t GetFileID() const { return fRunNumber>1000000?(fRunNumber%1000000)*1000+(fFileNumber%1000):fRunNumber; } TString GetStringID() const; - UShort_t GetTelescopeNumber() const { return fRunType; } + UShort_t GetTelescopeNumber() const { return fTelescopeNumber; } UShort_t GetRunType() const { return fRunType; } const Char_t *GetRunTypeStr() const; + const Char_t GetRunTypeChar() const; const Char_t *GetProjectName() const { return fProjectName; } const Char_t *GetSourceName() const { return fSourceName; } Index: trunk/MagicSoft/Mars/mreport/MReport.cc =================================================================== --- trunk/MagicSoft/Mars/mreport/MReport.cc (revision 9461) +++ trunk/MagicSoft/Mars/mreport/MReport.cc (revision 9462) @@ -169,4 +169,5 @@ // 200805190 | 54711.5 | | 200809030 // 200812040 | 54814.5 | | 200812140 +// 200902210 | 54968.5 | | 200905170 // Int_t MReport::Interprete(TString &str, const MTime &start, const MTime &stop, Int_t ver) @@ -208,4 +209,7 @@ ver=200812140; + if (ver==200902210 && GetMjd()>54968.5) + ver=200905170; + // Interprete body (contents) of report const Int_t rc = InterpreteBody(str, ver); Index: trunk/MagicSoft/Mars/mreport/MReportDrive.cc =================================================================== --- trunk/MagicSoft/Mars/mreport/MReportDrive.cc (revision 9461) +++ trunk/MagicSoft/Mars/mreport/MReportDrive.cc (revision 9462) @@ -135,4 +135,18 @@ } + if (ver>=200805170) + { + Int_t dummy; // Starguider switched on or not + n=sscanf(str.Data(), "%d %n", &dummy, &len); + if (n!=1) + { + *fLog << warn << "WARNING - Not enough arguments." << endl; + return kCONTINUE; + } + + str.Remove(0, len); + str = str.Strip(TString::kBoth); + } + return str.IsNull() ? kTRUE : kCONTINUE; } Index: trunk/MagicSoft/Mars/msimcamera/MSimAPD.cc =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimAPD.cc (revision 9461) +++ trunk/MagicSoft/Mars/msimcamera/MSimAPD.cc (revision 9462) @@ -206,5 +206,5 @@ return kERROR; } - +/* for (UInt_t idx=0; idx(fAPDs.UncheckedAt(idx))->FillRandom(freq, fStat->GetTimeFirst()); } +*/ + + // This tries to initialize dead and relaxing cells properly. If + // the APD has not been initialized before the chip is randomsly + // filled, otherwise a time window of the default relaxing time + // is simulated, so that the previous influence is less than a permille. + for (UInt_t idx=0; idx(fAPDs.UncheckedAt(idx))->Init(freq); + } // Get number of photons @@ -223,5 +234,5 @@ // Get arrival time of photon and idx - const Double_t t = ph.GetTime(); + const Double_t t = ph.GetTime()-fStat->GetTimeFirst(); const Int_t idx = ph.GetTag(); if (idx<0) @@ -239,5 +250,5 @@ // Simulate hitting the APD (the signal height in effective // "number of photons" is returned) - const Double_t hits = static_cast(fAPDs.UncheckedAt(idx))->HitRandomCell(t); + const Double_t hits = static_cast(fAPDs.UncheckedAt(idx))->HitRandomCellRelative(t); // FIXME: Make a proper simulation of the excess noise!!! @@ -248,4 +259,9 @@ } + // Now we have to shift the evolved time of all APDs to the end of our + // simulated time. + for (UInt_t idx=0; idx(fAPDs.UncheckedAt(idx))->IncreaseTime(fStat->GetTimeLast()); + return kTRUE; } Index: trunk/MagicSoft/Mars/msimcamera/MSimTrigger.cc =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimTrigger.cc (revision 9461) +++ trunk/MagicSoft/Mars/msimcamera/MSimTrigger.cc (revision 9462) @@ -42,7 +42,10 @@ // // With a second look-up table fCoincidenceMap the analog channels are -// checked for coincidences. The coincidense must at least be of the length +// checked for coincidences. The coincidence must at least be of the length // defined by fCoincidenceTime. The earliest coincide is then stored as // trigger position. +// +// If a minimum multiplicity m is given, m signals above threshold +// in the coincidence patterns are enough to emit a trigger signal. // // @@ -97,5 +100,6 @@ fEvtHeader(0), fElectronicNoise(0), fGain(0), fDiscriminatorThreshold(-1), fDigitalSignalLength(8), fCoincidenceTime(0.5), - fShiftBaseline(kTRUE), fUngainSignal(kTRUE) + fShiftBaseline(kTRUE), fUngainSignal(kTRUE), fSimulateElectronics(kTRUE), + fMinMultiplicity(-1) { fName = name ? name : "MSimTrigger"; @@ -107,13 +111,13 @@ // Take two TObjArrays with a collection of digital signals. // Every signal from one array is compared with any from the other array. -// For all signals whihc overlaps and which have an overlap time >gate +// For all signals which overlap and which have an overlap time >gate // a new digital signal is created storing start time and length of overlap. // They are collected in a newly allocated TObjArray. A pointer to this array // is returned. // -// Th euser gains owenership of the object, ie.e., the user is responsible of +// The user gains owenership of the object, i.e., the user is responsible of // deleting the memory. // -TObjArray *MSimTrigger::CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2, Float_t gate) const +TObjArray *MSimTrigger::CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2/*, Float_t gate*/) const { TObjArray *res = new TObjArray; @@ -131,5 +135,5 @@ { MDigitalSignal *ttl = new MDigitalSignal(*ttl1, *ttl2); - + /* if (ttl->GetLength()<=gate) { @@ -137,5 +141,5 @@ continue; } - + */ res->Add(ttl); } @@ -143,4 +147,106 @@ res->SetOwner(); + + return res; +} + +class Edge : public TObject +{ +private: + Double_t fEdge; + Int_t fRising; + +public: + Edge(Double_t t, Int_t rising) : fEdge(t), fRising(rising) { } + Bool_t IsSortable() const { return kTRUE; } + Int_t Compare(const TObject *o) const { const Edge *e = static_cast(o); if (e->fEdgefEdge>fEdge) return -1; return 0; } + + Int_t IsRising() const { return fRising; } + Double_t GetEdge() const { return fEdge; } +}; + +// -------------------------------------------------------------------------- +// +// Calculate a multiplicity trigger on the given array(s). The idx-array +// conatins all channels which should be checked for coincidences +// and the ttls array conatins the arrays with the digital signals. +// +// For the windows in which more or euqal than threshold channels have +// a high signal a new MDigitalSignal is created. newly allocated +// array with a collection of these trigger signals is returned. +// +TObjArray *MSimTrigger::CalcMinMultiplicity(const MArrayI &idx, const TObjArray &ttls, Int_t threshold) const +{ + // Create a new array for the rising and falling edges of the signals + TObjArray times; + times.SetOwner(); + + // Fill the array with edges from all digital signals of all our channels + for (UInt_t k=0; k(ttls[idx[k]]); + + TIter Next(arr); + MDigitalSignal *ttl = 0; + while ((ttl=static_cast(Next()))) + { + times.Add(new Edge(ttl->GetStart(), 1)); + times.Add(new Edge(ttl->GetEnd(), -1)); + } + } + + // Sort them in time + times.Sort(); + + // Start with no channel active + Int_t lvl = 0; + + TObjArray *res = new TObjArray; + res->SetOwner(); + + // First remove all edges which do not change the status + // "below threshold" or "above threshold" + for (int i=0; i(times.UncheckedAt(i)); + + // Claculate what the number of active channels after the edge is + const Int_t lvl1 = lvl + e.IsRising(); + + // Remove edge if number of active channels before and after the + // edge lower is lower than the threshold or higher than + // the threshold + if (lvl+1=threshold) + delete times.RemoveAt(i); + + // keep the (now) "previous" level + lvl = lvl1<0 ? 0 : lvl1; + } + + // Remove the empty slots from the array + times.Compress(); + + // + for (int i=0; i(times.UncheckedAt(i)); + + // go ahead if this is a falling edge + if (e0.IsRising()!=1) + continue; + + // get the following edge (must be a falling edge now) + const Edge &e1 = *static_cast(times.UncheckedAt(i+1)); + + // calculate the length of the digital signal + const Double_t len = e1.GetEdge()-e0.GetEdge(); + + // Create a digital trigger signal + MDigitalSignal *ds = new MDigitalSignal(e0.GetEdge(), len); + //ds->SetIndex(lvl); + res->Add(ds); + } return res; @@ -231,5 +337,5 @@ } - if (fElectronicNoise && !fRouteAC.IsEmpty() && !fRouteAC.IsDefault()) + if (fElectronicNoise && !fRouteAC.IsEmpty() && !fRouteAC.IsDefaultCol()) { // FIXME: Apply to analog channels when summing @@ -238,5 +344,5 @@ } - if (fGain && !fRouteAC.IsEmpty() && !fRouteAC.IsDefault()) + if (fGain && !fRouteAC.IsEmpty() && !fRouteAC.IsDefaultCol()) { // FIXME: Apply to analog channels when summing @@ -255,8 +361,19 @@ *fLog << "Using " << fNameRouteAC << " for re-routing/summing of analog channels before discriminator." << endl; - if (fCoincidenceMap.IsEmpty()) + if (fCoincidenceMap.IsEmpty() && fMinMultiplicity==0) *fLog << "No coincidences of digital channels will be checked. Signal-above-threshold trigger applied." << endl; else - *fLog << "Using " << fNameCoincidenceMap << " to check for coincidences of the digital channels." << endl; + { + *fLog << "Using "; + if (fCoincidenceMap.IsEmpty()) + *fLog << "the whole camera"; + else + *fLog << "patterns from " << fNameCoincidenceMap; + *fLog << " to check for "; + if (fMinMultiplicity==0) + *fLog << "coincidences of the digital channels." << endl; + else + *fLog << fMinMultiplicity << " multiplicity." << endl; + } *fLog << "Using discriminator threshold of " << fDiscriminatorThreshold << endl; @@ -265,4 +382,72 @@ } +/* +class MDigitalChannel : public TObjArray +{ +private: + TObjArray fArray; + +public: + MDigitalSignal *GetSignal(UInt_t i) { return static_cast(fArray[i]); } + +}; +*/ + +#include "MCamEvent.h" +class MTriggerSignal : public MParContainer, public MCamEvent +{ +private: + TObjArray fSignals; + +public: + MTriggerSignal() { fSignals.SetOwner(); } + + void Add(MDigitalSignal *signal) { fSignals.Add(signal); } + + MDigitalSignal *GetSignal(UInt_t i) { return static_cast(fSignals[i]); } + + void Sort() { fSignals.Sort(); } + + Int_t GetNumSignals() const { return fSignals.GetEntriesFast(); } + + Float_t GetFirstTrigger() const + { + MDigitalSignal *sig = static_cast(fSignals[0]); + return sig ? sig->GetStart() : -1; + } + Bool_t GetPixelContent(Double_t&, Int_t, const MGeomCam&, Int_t) const + { + switch (1) + { + case 1: // yes/no + case 2: // first time + case 3: // length + case 4: // n + break; + } + + return kTRUE; + } + void DrawPixelContent(Int_t) const + { + } +}; + + +void MSimTrigger::SetTrigger(Double_t pos) +{ + // FIXME: Jitter! (Own class?) + fTrigger->SetVal(pos); + fTrigger->SetReadyToSave(); + + // Trigger pattern to be set + // FIXME: Make flexible + const UInt_t pat = (UInt_t)~(MTriggerPattern::kTriggerLvl1 | (MTriggerPattern::kTriggerLvl1<<8)); + + // Flag this event as triggered by the lvl1 trigger (FIXME?) + fEvtHeader->SetTriggerPattern(pos<0 ? 0 : pat); + fEvtHeader->SetReadyToSave(); +} + // -------------------------------------------------------------------------- // @@ -270,5 +455,19 @@ { // Invalidate trigger - fTrigger->SetVal(-1); + //fTrigger->SetVal(-1); + // Calculate the minimum and maximum time for a valid trigger + const Double_t freq = fRunHeader->GetFreqSampling()/1000.; + const Float_t nsamp = fRunHeader->GetNumSamplesHiGain(); + const Float_t pulspos = fPulsePos->GetVal()/freq; + + // Valid range in units of bins + const Float_t min = fCamera->GetValidRangeMin()+pulspos; + const Float_t max = fCamera->GetValidRangeMax()-(nsamp-pulspos); + + if (!fSimulateElectronics) + { + SetTrigger(min); + return kTRUE; + } // ================== Simulate channel bundling ==================== @@ -297,5 +496,9 @@ { const UInt_t idx = row[j]; - (*patches)[i].AddSignal((*fCamera)[idx]); + + // FIXME: Shrinking the mapping table earlier (e.g. + // ReInit) would avoid a lot of if's + if (idxGetNumChannels()) + (*patches)[i].AddSignal((*fCamera)[idx]); } } @@ -329,18 +532,13 @@ // omitted in the loop if (fCoincidenceMap.IsEmpty()) - fCoincidenceMap.SetDefault(npatch); - - // Calculate the minimum and maximum time for a valid trigger - const Double_t freq = fRunHeader->GetFreqSampling()/1000.; - const Float_t nsamp = fRunHeader->GetNumSamplesHiGain(); - const Float_t pulspos = fPulsePos->GetVal()/freq; - - // Valid range in units of bins - const Float_t min = fCamera->GetValidRangeMin()+pulspos; - const Float_t max = fCamera->GetValidRangeMax()-(nsamp-pulspos); + { + if (fMinMultiplicity>0) + fCoincidenceMap.SetDefaultRow(npatch); + else + fCoincidenceMap.SetDefaultCol(npatch); + } // Create an array for the individual triggers - TObjArray triggers; - triggers.SetOwner(); + MTriggerSignal triggers; Int_t cnt = 0; @@ -352,30 +550,51 @@ const MArrayI &idx = fCoincidenceMap.GetRow(j); - // Start with a copy of the first coincidence channel - TObjArray *arr = static_cast(ttls[idx[0]]->Clone()); - arr->SetOwner(); - - // compare to all other channels in this coincidence patch, one by one - for (UInt_t k=1; kGetEntriesFast()>0; k++) - { - TObjArray *res = CalcCoincidence(*arr, *static_cast(ttls[idx[k]]), - fCoincidenceTime); - - // Delete the original array and keep the new one - delete arr; - arr = res; - } + TObjArray *arr = 0; + + if (fMinMultiplicity>0) + { + arr = CalcMinMultiplicity(idx, ttls, fMinMultiplicity); + } + else + { + arr = CalcMinMultiplicity(idx, ttls, idx.GetSize()); + /* + // Start with a copy of the first coincidence channel + arr = static_cast(ttls[idx[0]]->Clone()); + arr->SetOwner(); + + // compare to all other channels in this coincidence patch, one by one + for (UInt_t k=1; kGetEntriesFast()>0; k++) + { + TObjArray *res = CalcCoincidence(*arr, *static_cast(ttls[idx[k]]));//, fCoincidenceTime); + + // Delete the original array and keep the new one + delete arr; + arr = res; + }*/ + } + + // Count the number of totally emitted coincidence signals + cnt += arr->GetEntriesFast(); // Remove all signals which are not in the valid digitization range // (This is not the digitization window, but the region in which // the analog channels contain usefull data) + // and which are shorter than the defined coincidence gate. TIter Next(arr); MDigitalSignal *ttl = 0; while ((ttl=static_cast(Next()))) { + if (ttl->GetLength()Remove(ttl); + continue; + } + if (ttl->GetStart()Remove(ttl); rmlo++; + continue; } if (ttl->GetStart()>max) @@ -383,5 +602,9 @@ delete arr->Remove(ttl); rmhi++; + continue; } + + // Set trigger channel index + ttl->SetIndex(j); } @@ -389,22 +612,13 @@ arr->Compress(); - cnt += arr->GetEntriesFast(); - // If we have at least one trigger keep the earliest one. - // FIXME: The triggers should be ordered in time automatically: To be checked! + // FIXME: The triggers might be ordered in time automatically: + // To be checked! // FIXME: Simulate trigger dead-time! if (arr->GetEntriesFast()>0) - triggers.Add(arr->RemoveAt(0)); + triggers.Add(static_cast(arr->RemoveAt(0))); // delete the allocated space delete arr; - } - - // No trigger issued. Go on. - if (triggers.GetEntriesFast()==0) - { - if (rmlo>0 || rmhi>0) - *fLog << inf2 << rmlo << "/" << rmhi << " trigger out of valid range. No trigger raised." << endl; - return kTRUE; } @@ -414,17 +628,29 @@ // is omitted triggers.Sort(); - - // FIXME: Jitter! (Own class?) - fTrigger->SetVal(static_cast(triggers[0])->GetStart()); - fTrigger->SetReadyToSave(); - - // Flag this event as triggered by the lvl1 trigger (FIXME?) - fEvtHeader->SetTriggerPattern((UInt_t)~(MTriggerPattern::kTriggerLvl1 | (MTriggerPattern::kTriggerLvl1<<8))); - fEvtHeader->SetReadyToSave(); - - // inf2? - *fLog << inf; - *fLog << cnt << " triggers left in " << triggers.GetEntriesFast() << " patches (" << rmlo << "/" << rmhi << " trigger out of valid range), T=" << fTrigger->GetVal(); + // FIXME: Store triggers! (+ Reversed pixels?) + + SetTrigger(triggers.GetFirstTrigger()); + + // No trigger issued. Go on. + if (triggers.GetNumSignals()==0) + { + if (rmlo>0 || rmhi>0) + *fLog << inf2 << GetNumExecutions() << ": " << rmlo << "/" << rmhi << " trigger out of valid range. No trigger raised." << endl; + return kTRUE; + } + + // Number of patches which have triggered out of the total number of + // Coincidence signals emitted. (If the total number is higher than + // the number of triggers either some triggers had to be removed or + // or a patch has emitted more than one trigger signal) + // FIXME: inf2? + *fLog << inf << GetNumExecutions() << ": "; + *fLog << setw(3) << triggers.GetNumSignals() << " triggers left out of "; + *fLog << setw(3) << cnt << " (" << rmlo << "/" << rmhi << " trigger out of valid range), T=" << fTrigger->GetVal(); *fLog << endl; + + //# Trigger characteristics: gate length (ns), min. overlapping time (ns), + //# amplitude and FWHM of (gaussian) single phe response for trigger: + //trigger_prop 3.0 0.25 1.0 2.0 return kTRUE; @@ -438,4 +664,5 @@ // DigitalSignalLength: 8 // CoincidenceTime: 0.5 +// SimulateElectronics: Yes // Int_t MSimTrigger::ReadEnv(const TEnv &env, TString prefix, Bool_t print) @@ -472,4 +699,16 @@ } + if (IsEnvDefined(env, prefix, "SimulateElectronics", print)) + { + rc = kTRUE; + fSimulateElectronics = GetEnvValue(env, prefix, "SimulateElectronics", fSimulateElectronics); + } + + if (IsEnvDefined(env, prefix, "MinMultiplicity", print)) + { + rc = kTRUE; + fMinMultiplicity = GetEnvValue(env, prefix, "MinMultiplicity", fMinMultiplicity); + } + return rc; } Index: trunk/MagicSoft/Mars/msimcamera/MSimTrigger.h =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimTrigger.h (revision 9461) +++ trunk/MagicSoft/Mars/msimcamera/MSimTrigger.h (revision 9462) @@ -40,7 +40,13 @@ Bool_t fShiftBaseline; // Shift the baseline back to 0 for the threshold (needs ElectronicNoise [MPedestalCam]) Bool_t fUngainSignal; // "Remove" the gain from the signal (needs Gain [MPedestalCam]) + Bool_t fSimulateElectronics; // If the electronics is not simulated the trigger is set artificially to the first photon arrived + + Int_t fMinMultiplicity; // N out of M // MSimTrigger - TObjArray *CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2, Float_t gate=0) const; + TObjArray *CalcCoincidence(const TObjArray &arr1, const TObjArray &arr2/*, Float_t gate=0*/) const; + TObjArray *CalcMinMultiplicity(const MArrayI &idx, const TObjArray &ttls, Int_t threshold) const; + TObjArray *CalcCoincidences(const MArrayI &idx, const TObjArray &ttls) const; + void SetTrigger(Double_t pos); // MTask Index: trunk/MagicSoft/Mars/mtrigger/MTriggerPattern.h =================================================================== --- trunk/MagicSoft/Mars/mtrigger/MTriggerPattern.h (revision 9461) +++ trunk/MagicSoft/Mars/mtrigger/MTriggerPattern.h (revision 9462) @@ -12,12 +12,14 @@ public: enum Pattern_t { - kTriggerLvl1 = BIT(0), // 1: - kCalibration = BIT(1), // 2: Pulse Trigger - kTriggerLvl2 = BIT(2), // 4: LUT Pseudo Size selection - kPedestal = BIT(3), // 8: - kPinDiode = BIT(4), // 16: - kSumTrigger = BIT(5) // 32: Flag for an event taken with sum trigger - //kUnused = BIT(6) - //kUnused = BIT(7) +// kUndefined1 = BIT(0), // 1 1: Level 1 from L2 board + kTriggerLvl1 = BIT(0), // 1 1: Level 1 from L2 board + kCalibration = BIT(1), // 2 2: Pulse Trigger + kTriggerLvl2 = BIT(2), // 4 4: LUT Pseudo Size selection + kPedestal = BIT(3), // 8 8: Artificial pedestal event + kPinDiode = BIT(4), // 10 16: + kSumTrigger = BIT(5), // 20 32: Flag for an event taken with sum trigger +// kTriggerLvl1 = BIT(6), // 40 64: Trigger lvl1 directly from L1 without going through L2 + kUndefined1 = BIT(6), // 40 64: Trigger lvl1 directly from L1 without going through L2 + kUndefined2 = BIT(7) // 80 128: Undefined? }; Index: trunk/MagicSoft/Mars/mtrigger/MTriggerPatternDecode.cc =================================================================== --- trunk/MagicSoft/Mars/mtrigger/MTriggerPatternDecode.cc (revision 9461) +++ trunk/MagicSoft/Mars/mtrigger/MTriggerPatternDecode.cc (revision 9462) @@ -106,5 +106,5 @@ return kTRUE; - UInt_t pattern = ~fEvtHeader->GetTriggerID(); + const UInt_t pattern = ~fEvtHeader->GetTriggerID(); // The the trigger pattern is currently written with inverted bits, @@ -133,4 +133,14 @@ fPattern->fUnprescaled = (pattern>>8) & 0xff; + // This is a workaround for the new scheme in which L1TPU (the signal + // comming directly from the L1 is connected, but the L1 (routed + // over L2 is disconnected) + if (!fRunHeader->IsMonteCarloRun() && fRunHeader->GetTelescopeNumber()==1 && + fRunHeader->GetRunNumber()>1006246) + { + fPattern->fPrescaled |= (pattern>> 6)&1; + fPattern->fUnprescaled |= (pattern>>14)&1; + } + return kTRUE; }