Index: /branches/Mars_MC/macros/rootlogon.C =================================================================== --- /branches/Mars_MC/macros/rootlogon.C (revision 17010) +++ /branches/Mars_MC/macros/rootlogon.C (revision 17011) @@ -190,5 +190,6 @@ dbg += " -std=c++0x "; - lib += " -lnova -lz"; + //lib += " -lnova -lz"; + lib +=" -lnova -L/cm/shared/apps/fact/libs/libnova-0.14.0/src/.libs -lz"; gSystem->SetLinkedLibs(lib); Index: /branches/Mars_MC/mjobs/MJSimulation.cc =================================================================== --- /branches/Mars_MC/mjobs/MJSimulation.cc (revision 17010) +++ /branches/Mars_MC/mjobs/MJSimulation.cc (revision 17011) @@ -260,4 +260,76 @@ write.AddContainer("IncidentAngle", "Events", kFALSE); } + +void MJSimulation::SetupHeaderKeys(MWriteFitsFile &write) const +{ + write.SetHeaderKey("TELESCOP", "FACT", ""); + write.SetHeaderKey("PACKAGE", "MARS Cheobs", ""); + write.SetHeaderKey("VERSION", "1.0", ""); + write.SetHeaderKey("CREATOR", "Ceres", ""); + write.SetHeaderKey("EXTREL", 1., ""); + write.SetHeaderKey("COMPILED", "1970-01-01T00:00:00", ""); + write.SetHeaderKey("REVISION", "0", ""); + write.SetHeaderKey("ORIGIN", "FACT", ""); + write.SetHeaderKey("DATE", "1970-01-01T00:00:00", ""); + write.SetHeaderKey("NIGHT", 19700101, ""); + write.SetHeaderKey("TIMESYS", "UTC", ""); + write.SetHeaderKey("TIMEUNIT", "d", ""); + write.SetHeaderKey("MJDREF", 40587, ""); + write.SetHeaderKey("BLDVER", 1, ""); + write.SetHeaderKey("RUNID", 1, ""); + write.SetHeaderKey("NBOARD", 40, ""); + write.SetHeaderKey("NPIX", 1440, ""); + write.SetHeaderKey("NROI", 150, ""); + write.SetHeaderKey("NROITM", 0, ""); + write.SetHeaderKey("TMSHIFT", 0, ""); + write.SetHeaderKey("CAMERA", "MGeomCamFACT", ""); + write.SetHeaderKey("DAQ", "DRS4", ""); + write.SetHeaderKey("ADCRANGE", 2000, "Dynamic range in mV"); + write.SetHeaderKey("ADC", 12, "Resolution in bits"); + write.SetHeaderKey("RUNTYPE", "data", ""); + write.SetHeaderKey("ID", 777, "Board 0: Board ID"); + write.SetHeaderKey("FMVER", 532, "Board 0: Firmware Version"); + write.SetHeaderKey("DNA", "0", ""); + write.SetHeaderKey("BOARD", 0, ""); + write.SetHeaderKey("PRESC", 40, ""); + write.SetHeaderKey("PHASE", 0, ""); + write.SetHeaderKey("DAC0", 26500, ""); + write.SetHeaderKey("DAC1", 0, ""); + write.SetHeaderKey("DAC2", 0, ""); + write.SetHeaderKey("DAC3", 0, ""); + write.SetHeaderKey("DAC4", 28800, ""); + write.SetHeaderKey("DAC5", 28800, ""); + write.SetHeaderKey("DAC6", 28800, ""); + write.SetHeaderKey("DAC7", 28800, ""); + write.SetHeaderKey("REFCLK", 2000.0768, ""); + write.SetHeaderKey("DRSCALIB", false, ""); + write.SetHeaderKey("TSTARTI", 0, ""); + write.SetHeaderKey("TSTARTF", 0., ""); + write.SetHeaderKey("TSTOPI", 0, ""); + write.SetHeaderKey("TSTOPF", 0., ""); + write.SetHeaderKey("DATE-OBS", "1970-01-01T00:00:00", ""); + write.SetHeaderKey("DATE-END", "1970-01-01T00:00:00", ""); + write.SetHeaderKey("NTRG", 0, ""); + write.SetHeaderKey("NTRGPED", 0, ""); + write.SetHeaderKey("NTRGLPE", 0, ""); + write.SetHeaderKey("NTRGTIM", 0, ""); + write.SetHeaderKey("NTRGLPI", 0, ""); + write.SetHeaderKey("NTRGEXT1", 0, ""); + write.SetHeaderKey("NTRGEXT2", 0, ""); + write.SetHeaderKey("NTRGMISC", 0, ""); +} + +void MJSimulation::SetupVetoColumns(MWriteFitsFile &write) const +{ + write.VetoColumn("MParameterD.fVal"); + write.VetoColumn("MRawEvtData.fLoGainPixId"); + write.VetoColumn("MRawEvtData.fLoGainFadcSamples"); + write.VetoColumn("MRawEvtData.fABFlags"); + write.VetoColumn("MRawEvtData.fIsSigned"); + write.VetoColumn("MRawEvtHeader.fNumTrigLvl2"); + write.VetoColumn("MRawEvtHeader.fTrigPattern"); + write.VetoColumn("MRawEvtHeader.fNumLoGainOn"); +} + Bool_t MJSimulation::Process(const MArgs &args, const MSequence &seq) { @@ -558,66 +630,6 @@ MTask &write3a = fWriteFitsFile ? static_cast(write3af) : static_cast(write3ar); - write3af.SetHeaderKey("DUMMY0", 3, "First dummy"); - write3af.SetHeaderKey("DUMMY1", 3.14159265358979323846, "Second dummy"); - write3af.SetHeaderKey("DUMMY2", true, "Third dummy"); - write3af.SetHeaderKey("DUMMY3", "one value", "Fourth dummy"); - - write3af.VetoColumn("MParameterD.fVal"); - write3af.VetoColumn("MCorsikaEvtHeader.fEvtNumber"); - write3af.VetoColumn("MCorsikaEvtHeader.fNumReuse"); - write3af.VetoColumn("MCorsikaEvtHeader.fTotalEnergy"); - write3af.VetoColumn("MCorsikaEvtHeader.fStartAltitude"); - write3af.VetoColumn("MCorsikaEvtHeader.fFirstTargetNum"); - write3af.VetoColumn("MCorsikaEvtHeader.fFirstInteractionHeight"); - write3af.VetoColumn("MCorsikaEvtHeader.fMomentumX"); - write3af.VetoColumn("MCorsikaEvtHeader.fMomentumY"); - write3af.VetoColumn("MCorsikaEvtHeader.fMomentumZ"); - write3af.VetoColumn("MCorsikaEvtHeader.fAz"); - write3af.VetoColumn("MCorsikaEvtHeader.fWeightedNumPhotons"); - write3af.VetoColumn("MCorsikaEvtHeader.fZd"); - write3af.VetoColumn("MCorsikaEvtHeader.fAd"); - write3af.VetoColumn("MCorsikaEvtHeader.fX"); - write3af.VetoColumn("MCorsikaEvtHeader.fY"); - write3af.VetoColumn("MCorsikaEvtHeader.fWeightNumPhotons"); - write3af.VetoColumn("MMcEvt.fEvtNumber"); - write3af.VetoColumn("MMcEvt.fThick0"); - write3af.VetoColumn("MMcEvt.fFirstTarget"); - write3af.VetoColumn("MMcEvt.fZFirstInteraction"); - write3af.VetoColumn("MMcEvt.fCoreD"); - write3af.VetoColumn("MMcEvt.fCoreX"); - write3af.VetoColumn("MMcEvt.fCoreY"); - write3af.VetoColumn("MMcEvt.fTimeFirst"); - write3af.VetoColumn("MMcEvt.fTimeLast"); - write3af.VetoColumn("MMcEvt.fLongiNmax"); - write3af.VetoColumn("MMcEvt.fLongit0"); - write3af.VetoColumn("MMcEvt.fLongitmax"); - write3af.VetoColumn("MMcEvt.fLongia"); - write3af.VetoColumn("MMcEvt.fLongib"); - write3af.VetoColumn("MMcEvt.fLongic"); - write3af.VetoColumn("MMcEvt.fLongichi2"); - write3af.VetoColumn("MMcEvt.fPhotIni"); - write3af.VetoColumn("MMcEvt.fPassPhotAtm"); - write3af.VetoColumn("MMcEvt.fPassPhotRef"); - write3af.VetoColumn("MMcEvt.fPassPhotCone"); - write3af.VetoColumn("MMcEvt.fPhotElfromShower"); - write3af.VetoColumn("MMcEvt.fPhotElinCamera"); - write3af.VetoColumn("MMcEvt.fElecCphFraction"); - write3af.VetoColumn("MMcEvt.fMuonCphFraction"); - write3af.VetoColumn("MMcEvt.fOtherCphFraction"); - write3af.VetoColumn("MMcEvt.fFadcTimeJitter"); - write3af.VetoColumn("MMcEvt.fEventReuse"); - - write3af.VetoColumn("MRawEvtData.fHiGainPixId"); - write3af.VetoColumn("MRawEvtData.fLoGainPixId"); - write3af.VetoColumn("MRawEvtData.fLoGainFadcSamples"); - write3af.VetoColumn("MRawEvtData.fABFlags"); - write3af.VetoColumn("MRawEvtData.fStartCells"); - write3af.VetoColumn("MRawEvtData.fNumBytesPerSample"); - write3af.VetoColumn("MRawEvtData.fIsSigned"); - write3af.VetoColumn("MRawEvtHeader.fDAQEvtNumber"); //EventNum ? - write3af.VetoColumn("MRawEvtHeader.fNumTrigLvl1"); - write3af.VetoColumn("MRawEvtHeader.fNumTrigLvl2"); - write3af.VetoColumn("MRawEvtHeader.fTrigPattern"); - write3af.VetoColumn("MRawEvtHeader.fNumLoGainOn"); + SetupHeaderKeys(write3af); + SetupVetoColumns(write3af); write3af.SetBytesPerSample("Data", 2); @@ -649,4 +661,5 @@ write3af.AddContainer("IntendedPulsePos", "RunHeaders"); write3af.AddContainer("MRawEvtData", "Events"); + write3af.AddContainer("MPhotonEvent","Events"); write3ar.AddContainer("ElectronicNoise", "RunHeaders"); Index: /branches/Mars_MC/mjobs/MJSimulation.h =================================================================== --- /branches/Mars_MC/mjobs/MJSimulation.h (revision 17010) +++ /branches/Mars_MC/mjobs/MJSimulation.h (revision 17011) @@ -34,4 +34,6 @@ void SetupCommonFileStructure(MWriteRootFile &write) const; void SetupCommonFileStructure(MWriteFitsFile& write) const; + void SetupHeaderKeys(MWriteFitsFile& write) const; + void SetupVetoColumns(MWriteFitsFile& write) const; public: Index: /branches/Mars_MC/mmc/MMcEvt.hxx =================================================================== --- /branches/Mars_MC/mmc/MMcEvt.hxx (revision 17010) +++ /branches/Mars_MC/mmc/MMcEvt.hxx (revision 17011) @@ -10,5 +10,5 @@ { private: - UInt_t fEvtNumber; + UInt_t fEvtNumber; // {fits: name=EventNum ; unit=uint32} Float_t fThick0; // [g/cm2] Float_t fFirstTarget; // [] Index: /branches/Mars_MC/mraw/MRawEvtData.h =================================================================== --- /branches/Mars_MC/mraw/MRawEvtData.h (revision 17010) +++ /branches/Mars_MC/mraw/MRawEvtData.h (revision 17011) @@ -36,5 +36,5 @@ MArrayB *fABFlags; //-> A information about the exact trigger position - MArrayS *fStartCells; // + MArrayS *fStartCells; // {fits: unit=uint16 ; name=StartCellData } UShort_t fNumBytesPerSample; Index: /branches/Mars_MC/msimcamera/MSimCamera.cc =================================================================== --- /branches/Mars_MC/msimcamera/MSimCamera.cc (revision 17010) +++ /branches/Mars_MC/msimcamera/MSimCamera.cc (revision 17011) @@ -135,4 +135,11 @@ return kFALSE; + fAccidentalPhotons = (MPedestalCam*)pList->FindObject("AccidentalPhotonRates","MPedestalCam"); + if(!fAccidentalPhotons) + { + *fLog << err << "AccidentalPhotonRates [MPedestalCam] not found... aborting." << endl; + return kFALSE; + } + MParSpline *pulse = (MParSpline*)pList->FindObject("PulseShape", "MParSpline"); if (!pulse) @@ -238,7 +245,45 @@ const Double_t rms = pix.GetPedestalRms(); + // FTemme: Implementation of AC-coupling: + // to calculate the value of the accoupling per slice I use the + // following equation: + // accouplingPerSlice = accidentalPhotonRate * (1 + crossTalkProb) + // * areaOfOnePulse / samplingRate; + // Therefore I need the following variables + Double_t accidentalPhotonRate = 0; // [MHz] + Float_t crossTalkProb = 0; // [1] + Double_t areaOfOnePulse = 0; // [ADC-Counts * s] + Double_t samplingRate = 0; // [slices * MHz] + + // The accidental photon rate is stored in GHz, so we have to multiply + // with 1E3 to get MHz: + MPedestalPix &accPhoPix = (*fAccidentalPhotons)[i]; + accidentalPhotonRate = accPhoPix.GetPedestal() * 1E3; + + // I don't know how to get the variable fCrosstalkProb from + // the class APD (see MAvalanchePhotoDiode.h), because there is no + // getter for the APD array(fAPDs) in MSimAPD. + // So I set the crossTalkProb hardcoded to the value 0.15, which is + // equal to the value of the apd of type 4 + crossTalkProb = 0.15; + + // To get the area of one Pulse, I only need to calculate the Integral + // of the Pulse Shape, which is stored in fSpline. Because the spline is + // normalized to a maximal amplitude of 1.0, I had to multiply it with + // the Default gain: + areaOfOnePulse = fSpline->Integral() * fDefaultGain; + + // The sampling rate I get from the RunHeader: + samplingRate = fRunHeader->GetFreqSampling(); + + Double_t accouplingPerSlice = accidentalPhotonRate * (1 + crossTalkProb) + * areaOfOnePulse / samplingRate; + + // The accoupling is substracted from the timeline by decreasing the + // mean of the gaussian noise which is added + if (!fBaselineGain) { - (*fCamera)[i].AddGaussianNoise(rms, val); + (*fCamera)[i].AddGaussianNoise(rms, val - accouplingPerSlice); continue; } @@ -250,5 +295,6 @@ // FIXME: How stable is the offset? // FIXME: Should we write a container AppliedGain for MSImTrigger? - (*fCamera)[i].AddGaussianNoise(rms*gain, val*gain); + + (*fCamera)[i].AddGaussianNoise(rms*gain, (val - accouplingPerSlice)*gain); } Index: /branches/Mars_MC/msimcamera/MSimCamera.h =================================================================== --- /branches/Mars_MC/msimcamera/MSimCamera.h (revision 17010) +++ /branches/Mars_MC/msimcamera/MSimCamera.h (revision 17011) @@ -24,9 +24,10 @@ MPedestalCam *fElectronicNoise; //! Electronic noise (baseline and rms) MPedestalCam *fGain; //! Electronic noise (baseline and rms) + MPedestalCam *fAccidentalPhotons;//! Accidental photon rates MAnalogChannels *fCamera; //! Output of the analog signals MMcEvt *fMcEvt; //! For information stored in MMcEvt - const MSpline3 *fSpline; // Pulse Shape + const MSpline3 *fSpline; // Pulse Shape Bool_t fBaselineGain; // Should the gain be applied to baseline and electronic noise? Index: /branches/Mars_MC/msimcamera/MSimReadout.cc =================================================================== --- /branches/Mars_MC/msimcamera/MSimReadout.cc (revision 17010) +++ /branches/Mars_MC/msimcamera/MSimReadout.cc (revision 17011) @@ -145,4 +145,5 @@ fData->InitRead(fRunHeader); fData->ResetPixels(); + fData->InitStartCells(); fData->SetIndices(); @@ -202,5 +203,5 @@ // Check if the whole requested signal can be digitized - if (trig+nslices>=nsamp) + if (trig+nslices>nsamp) { *fLog << err << "ERROR - Trigger position beyond valid analog signal range." << endl;