//#endif
// #include "AliFMDGeometryBuilder.h"
#include "AliFMDRawWriter.h" // ALIFMDRAWWRITER_H
+#include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
#include "AliTrackReference.h"
#include "AliFMDStripIndex.h"
+#include "AliFMDParameters.h"
+#include "AliFMDReconstructor.h"
//____________________________________________________________________
ClassImp(AliFMD)
//____________________________________________________________________
void
-AliFMD::AddDigitByFields(UShort_t detector,
- Char_t ring,
- UShort_t sector,
- UShort_t strip,
- UShort_t count1,
- Short_t count2,
- Short_t count3,
- Short_t count4)
+AliFMD::AddDigitByFields(UShort_t detector,
+ Char_t ring,
+ UShort_t sector,
+ UShort_t strip,
+ UShort_t count1,
+ Short_t count2,
+ Short_t count3,
+ Short_t count4,
+ const TArrayI& refs)
{
// add a real digit - as coming from data
//
TClonesArray& a = *(DigitsArray());
new (a[fNdigits++])
- AliFMDDigit(detector, ring, sector, strip, count1, count2, count3, count4);
- AliFMDDebug(15, ("Adding digit # %5d/%5d for FMD%d%c[%2d,%3d]=(%d,%d,%d,%d)",
+ AliFMDDigit(detector, ring, sector, strip,
+ count1, count2, count3, count4, refs);
+ AliFMDDebug(15, ("Adding digit # %5d/%5d for FMD%d%c[%2d,%3d]"
+ "=(%d,%d,%d,%d) with %d tracks",
fNdigits-1, a.GetEntriesFast(),
detector, ring, sector, strip,
- count1, count2, count3, count4));
+ count1, count2, count3, count4, refs.fN));
}
writer.Exec();
}
+//====================================================================
+//
+// Raw data reading
+//
+//__________________________________________________________________
+Bool_t
+AliFMD::Raw2SDigits(AliRawReader* reader)
+{
+ // Turn digits into raw data.
+ //
+ // This uses the class AliFMDRawWriter to do the job. Please refer
+ // to that class for more information.
+ AliFMDParameters::Instance()->Init();
+ MakeTree("S");
+ MakeBranch("S");
+
+ TClonesArray* sdigits = SDigits();
+ AliFMDReconstructor rec;
+
+ // The two boolean arguments
+ // Make sdigits instead of digits
+ // Subtract the pedestal off the signal
+ rec.Digitize(reader, sdigits);
+ //
+ // Bool_t ret = fmdReader.ReadAdcs(sdigits, kTRUE, kTRUE);
+ // sdigits->ls();
+ UShort_t ns = sdigits->GetEntriesFast();
+ for (UShort_t i = 0; i < ns; i++)
+ sdigits->At(i)->Print("pl");
+
+ AliFMDDebug(1, ("Got a total of %d SDigits", ns));
+
+ fLoader->TreeS()->Fill();
+ ResetSDigits();
+ fLoader->WriteSDigits("OVERWRITE");
+
+ return kTRUE;
+}
+
//====================================================================
//
@param count1 ADC count (a 10-bit word)
@param count2 ADC count (a 10-bit word), or -1 if not used
@param count3 ADC count (a 10-bit word), or -1 if not used */
- virtual void AddDigitByFields(UShort_t detector=0,
- Char_t ring='\0',
- UShort_t sector=0,
- UShort_t strip=0,
- UShort_t count1=0,
- Short_t count2=-1,
- Short_t count3=-1,
- Short_t count4=-1);
+ virtual void AddDigitByFields(UShort_t detector=0,
+ Char_t ring='\0',
+ UShort_t sector=0,
+ UShort_t strip=0,
+ UShort_t count1=0,
+ Short_t count2=-1,
+ Short_t count3=-1,
+ Short_t count4=-1,
+ const TArrayI& refs=TArrayI());
/** Add a digit to the Digit tree
@param digits
- digits[0] [UShort_t] Detector #
/** @{ */
/** @name Raw data */
- /** Turn digits into raw data. This uses the class AliFMDRawWriter
- to do the job. Please refer to that class for more
- information. */
- virtual void Digits2Raw();
+ /**
+ * Turn digits into raw data. This uses the class AliFMDRawWriter to
+ * do the job. Please refer to that class for more information.
+ */
+ virtual void Digits2Raw();
+ /**
+ * Convert raw data to sdigits
+ *
+ * @param reader Raw reader
+ *
+ * @return @c true on success
+ */
+ virtual Bool_t Raw2SDigits(AliRawReader* reader);
/** @}*/
/** @{ */
- /** @name Utility */
- /** Browse this object
- @param b Browser to show this object in */
+ /**
+ * @name Utility
+ */
+ /**
+ * Browse this object
+ *
+ * @param b Browser to show this object in
+ */
void Browse(TBrowser* b);
/** @}*/
protected:
}
//_____________________________________________________________________
-Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::ProcessHit( AliFMDHit* h, TParticle* p) {
+Bool_t
+AliFMDBackgroundCorrection::AliFMDInputBG::ProcessHit(AliFMDHit* h,
+ TParticle* /*p*/)
+{
if(!h)
return kTRUE;
- Bool_t retval = ProcessEvent(h->Detector(),h->Ring(),h->Sector(),h->Strip(),h->Track(),h->Q());
+ Bool_t retval = ProcessEvent(h->Detector(),
+ h->Ring(),
+ h->Sector()
+ ,h->Strip(),
+ h->Track(),
+ h->Q());
return retval;
}
//_____________________________________________________________________
-Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::ProcessTrackRef( AliTrackReference* tr, TParticle* p) {
-
+Bool_t
+AliFMDBackgroundCorrection::AliFMDInputBG::ProcessTrackRef(AliTrackReference* tr,
+ TParticle* p)
+{
if(!tr)
return kTRUE;
UShort_t det,sec,strip;
Char_t ring;
AliFMDStripIndex::Unpack(tr->UserId(),det,ring,sec,strip);
- Int_t nTrack = tr->GetTrack();
- TDatabasePDG* pdgDB = TDatabasePDG::Instance();
+ Int_t nTrack = tr->GetTrack();
+ TDatabasePDG* pdgDB = TDatabasePDG::Instance();
TParticlePDG* pdgPart = pdgDB->GetParticle(p->GetPdgCode());
- Float_t charge = (pdgPart ? pdgPart->Charge() : 0);
- Bool_t retval = ProcessEvent(det,ring,sec,strip,nTrack,charge);
+ Float_t charge = (pdgPart ? pdgPart->Charge() : 0);
+ Bool_t retval = ProcessEvent(det,ring,sec,strip,nTrack,charge);
return retval;
}
//_____________________________________________________________________
-Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::ProcessEvent(UShort_t det,
- Char_t ring,
- UShort_t sec,
- UShort_t strip,
- Int_t nTrack,
- Float_t charge)
+Bool_t
+AliFMDBackgroundCorrection::AliFMDInputBG::ProcessEvent(UShort_t det,
+ Char_t ring,
+ UShort_t sec,
+ UShort_t strip,
+ Int_t nTrack,
+ Float_t charge)
{
AliGenEventHeader* genHeader = fLoader->GetHeader()->GenEventHeader();
if(TMath::Abs(vertex.At(2)) > fZvtxCut)
return kTRUE;
- Double_t delta = 2*fZvtxCut/fNvtxBins;
+ Double_t delta = 2 * fZvtxCut / fNvtxBins;
Double_t vertexBinDouble = (vertex.At(2) + fZvtxCut) / delta;
- Int_t vertexBin = (Int_t)vertexBinDouble;
+ Int_t vertexBin = Int_t(vertexBinDouble);
- if(charge != 0 && (nTrack != fPrevTrack || det != fPrevDetector || ring != fPrevRing)) {
-
+ if(charge != 0 &&
+ (nTrack != fPrevTrack ||
+ det != fPrevDetector ||
+ ring != fPrevRing)) {
Double_t x,y,z;
AliFMDGeometry* fmdgeo = AliFMDGeometry::Instance();
fmdgeo->Detector2XYZ(det,ring,sec,strip,x,y,z);
TObjArray* vtxArray = (TObjArray*)detArray->At(iring);
TH2F* hHits = (TH2F*)vtxArray->At(vertexBin);
- Float_t phi = TMath::ATan2(y,x);
- if(phi<0)
- phi = phi+2*TMath::Pi();
- Float_t theta = TMath::ATan2(TMath::Sqrt(TMath::Power(x,2)+TMath::Power(y,2)),z+vertex.At(2));
- Float_t eta = -1*TMath::Log(TMath::Tan(0.5*theta));
+ Float_t phi = TMath::ATan2(y,x);
+ if(phi<0) phi = phi+2*TMath::Pi();
+ Float_t r = TMath::Sqrt(TMath::Power(x,2)+TMath::Power(y,2));
+ Float_t theta = TMath::ATan2(r,z+vertex.At(2));
+ Float_t eta = -1*TMath::Log(TMath::Tan(0.5*theta));
hHits->Fill(eta,phi);
fHits++;
fPrevDetector = det;
}
//_____________________________________________________________________
-Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Init() {
-
+Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Init()
+{
fPrimaryArray.SetOwner();
fPrimaryArray.SetName("FMD_primary");
fHitArray.SetOwner();
fHitArray.SetName("FMD_hits");
- for(Int_t det =1; det<=3;det++)
- {
- TObjArray* detArrayHits = new TObjArray();
- detArrayHits->SetName(Form("FMD%d",det));
- fHitArray.AddAtAndExpand(detArrayHits,det);
- Int_t nRings = (det==1 ? 1 : 2);
- for(Int_t ring = 0;ring<nRings;ring++)
- {
- Int_t nSec = (ring == 1 ? 40 : 20);
- Char_t ringChar = (ring == 0 ? 'I' : 'O');
- TObjArray* vtxArrayHits = new TObjArray();
- vtxArrayHits->SetName(Form("FMD%d%c",det,ringChar));
- detArrayHits->AddAtAndExpand(vtxArrayHits,ring);
- for(Int_t v=0; v<fNvtxBins;v++)
- {
-
- TH2F* hHits = new TH2F(Form("hHits_FMD%d%c_vtx%d",det,ringChar,v),
- Form("hHits_FMD%d%c_vtx%d",det,ringChar,v),
- fNbinsEta, -6,6, nSec, 0, 2*TMath::Pi());
- hHits->Sumw2();
- vtxArrayHits->AddAtAndExpand(hHits,v);
-
- }
- }
-
+ for(Int_t det =1; det<=3;det++) {
+ TObjArray* detArrayHits = new TObjArray();
+ detArrayHits->SetName(Form("FMD%d",det));
+ fHitArray.AddAtAndExpand(detArrayHits,det);
+ Int_t nRings = (det==1 ? 1 : 2);
+ for(Int_t ring = 0;ring<nRings;ring++) {
+ Int_t nSec = (ring == 1 ? 40 : 20);
+ Char_t ringChar = (ring == 0 ? 'I' : 'O');
+ TObjArray* vtxArrayHits = new TObjArray();
+ vtxArrayHits->SetName(Form("FMD%d%c",det,ringChar));
+ detArrayHits->AddAtAndExpand(vtxArrayHits,ring);
+ for(Int_t v=0; v<fNvtxBins;v++) {
+ TH2F* hHits = new TH2F(Form("hHits_FMD%d%c_vtx%d", det,ringChar,v),
+ Form("hHits_FMD%d%c_vtx%d", det,ringChar,v),
+ fNbinsEta, -6,6, nSec, 0, 2*TMath::Pi());
+ hHits->Sumw2();
+ vtxArrayHits->AddAtAndExpand(hHits,v);
+
+ }
}
+ }
-
-
AliFMDInput::Init();
-
return kTRUE;
}
//_____________________________________________________________________
-Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Begin(Int_t event ) {
+Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Begin(Int_t event )
+{
- Bool_t retVal = AliFMDInput::Begin(event);
-
- AliStack* partStack=fLoader->Stack();
-
- Int_t nTracks=partStack->GetNtrack();
+ Bool_t retVal = AliFMDInput::Begin(event);
+ AliStack* partStack = fLoader->Stack();
+ Int_t nTracks = partStack->GetNtrack();
AliGenEventHeader* genHeader = fLoader->GetHeader()->GenEventHeader();
TArrayF vertex;
genHeader->PrimaryVertex(vertex);
if(TMath::Abs(vertex.At(2)) > fZvtxCut)
return kTRUE;
- Double_t delta = 2*fZvtxCut/fNvtxBins;
+ Double_t delta = 2*fZvtxCut/fNvtxBins;
Double_t vertexBinDouble = (vertex.At(2) + fZvtxCut) / delta;
- Int_t vertexBin = (Int_t)vertexBinDouble;
-
- for(Int_t j=0;j<nTracks;j++)
- {
- TParticle* p = partStack->Particle(j);
- TDatabasePDG* pdgDB = TDatabasePDG::Instance();
- TParticlePDG* pdgPart = pdgDB->GetParticle(p->GetPdgCode());
- Float_t charge = (pdgPart ? pdgPart->Charge() : 0);
- Float_t phi = TMath::ATan2(p->Py(),p->Px());
-
- if(phi<0)
- phi = phi+2*TMath::Pi();
- if(p->Theta() == 0) continue;
- Float_t eta = -1*TMath::Log(TMath::Tan(0.5*p->Theta()));
-
- Bool_t primary = partStack->IsPhysicalPrimary(j);
- //(charge!=0)&&(TMath::Abs(p->Vx() - vertex.At(0))<0.01)&&(TMath::Abs(p->Vy() - vertex.At(1))<0.01)&&(TMath::Abs(p->Vz() - vertex.At(2))<0.01);
- if(primary && charge !=0) {
-
- fPrim++;
- TObjArray* innerArray = (TObjArray*)fPrimaryArray.At(0);
- TObjArray* outerArray = (TObjArray*)fPrimaryArray.At(1);
-
- TH2F* hPrimaryInner = (TH2F*)innerArray->At(vertexBin);
- TH2F* hPrimaryOuter = (TH2F*)outerArray->At(vertexBin);
- hPrimaryInner->Fill(eta,phi);
- hPrimaryOuter->Fill(eta,phi);
- }
- }
+ Int_t vertexBin = (Int_t)vertexBinDouble;
+
+ for(Int_t j=0;j<nTracks;j++) {
+ TParticle* p = partStack->Particle(j);
+ TDatabasePDG* pdgDB = TDatabasePDG::Instance();
+ TParticlePDG* pdgPart = pdgDB->GetParticle(p->GetPdgCode());
+ Float_t charge = (pdgPart ? pdgPart->Charge() : 0);
+ Float_t phi = TMath::ATan2(p->Py(),p->Px());
+
+ if(phi<0) phi = phi+2*TMath::Pi();
+ if(p->Theta() == 0) continue;
+ Float_t eta = -1*TMath::Log(TMath::Tan(0.5*p->Theta()));
+
+ Bool_t primary = partStack->IsPhysicalPrimary(j);
+ //(charge!=0)&&(TMath::Abs(p->Vx() - vertex.At(0))<0.01)&&(TMath::Abs(p->Vy() - vertex.At(1))<0.01)&&(TMath::Abs(p->Vz() - vertex.At(2))<0.01);
+ if(primary && charge !=0) {
+
+ fPrim++;
+ TObjArray* innerArray = (TObjArray*)fPrimaryArray.At(0);
+ TObjArray* outerArray = (TObjArray*)fPrimaryArray.At(1);
+
+ TH2F* hPrimaryInner = (TH2F*)innerArray->At(vertexBin);
+ TH2F* hPrimaryOuter = (TH2F*)outerArray->At(vertexBin);
+ hPrimaryInner->Fill(eta,phi);
+ hPrimaryOuter->Fill(eta,phi);
+ }
+ }
return retVal;
}
//////////////////////////////////////////////////////////////////////
#include "AliFMDBaseDigit.h" // ALIFMDDIGIT_H
+#include "AliFMDStripIndex.h"
#include "Riostream.h" // ROOT_Riostream
// #include <TString.h>
// #include <AliLog.h>
Char_t ring,
UShort_t sector,
UShort_t strip)
- : fDetector(detector),
+ : AliDigit(),
+ fDetector(detector),
fRing(ring),
fSector(sector),
fStrip(strip),
// strip Strip # (For inner/outer rings: 0-511/0-255)
}
+//____________________________________________________________________
+AliFMDBaseDigit::AliFMDBaseDigit(Int_t* tracks,
+ UShort_t detector,
+ Char_t ring,
+ UShort_t sector,
+ UShort_t strip)
+ : AliDigit(tracks),
+ fDetector(detector),
+ fRing(ring),
+ fSector(sector),
+ fStrip(strip),
+ fName("")
+{
+ //
+ // Creates a base data digit object
+ //
+ // Parameters
+ //
+ // tracks Array of 3 track labels
+ // detector Detector # (1, 2, or 3)
+ // ring Ring ID ('I' or 'O')
+ // sector Sector # (For inner/outer rings: 0-19/0-39)
+ // strip Strip # (For inner/outer rings: 0-511/0-255)
+}
+
//____________________________________________________________________
void
AliFMDBaseDigit::Print(Option_t* /* option*/) const
{
// Print digit to standard out
- cout << ClassName() << ": FMD" << GetName() << flush;
+ cout << ClassName() << ": " << GetName() << flush;
}
//____________________________________________________________________
AliFMDBaseDigit::Hash() const
{
// Calculate a hash value based on the detector coordinates.
+#if 1
+ return AliFMDStripIndex::Pack(fDetector, fRing, fSector, fStrip);
+#else
size_t ringi = (fRing == 'I' || fRing == 'i' ? 0 : 1);
return fStrip + fMaxStrips *
(fSector + fMaxSectors * (ringi + fMaxRings * (fDetector - 1)));
+#endif
}
return 1;
}
+//____________________________________________________________________
+void
+AliFMDBaseDigit::AddTrack(Int_t track)
+{
+ if (fTracks[0] == -1) fTracks[0] = track;
+ else if (fTracks[1] == -1) fTracks[1] = track;
+ else if (fTracks[2] == -1) fTracks[2] = track;
+ else
+ AliWarning(Form("While adding track label to %s for %s: "
+ "All 3 track labels used, cannot add reference to track %d",
+ ClassName(), GetName(), track));
+}
+
+//____________________________________________________________________
+UShort_t
+AliFMDBaseDigit::GetNTrack() const
+{
+ for (Int_t i = 3; i > 0; i--)
+ if (fTracks[i-1] != -1) return i;
+ return 0;
+}
+
//____________________________________________________________________
//
// AliFMDDigit - Normal (smeared) digit
// AliFMDSDigit - Summable (non-smeared) digit
//
-#ifndef ROOT_TObject
-# include <TObject.h>
+#ifndef ALIDIGIT_H
+# include <AliDigit.h>
#endif
#ifndef ROOT_TString
# include <TString.h>
#endif
//____________________________________________________________________
-/** @class AliFMDBaseDigit AliFMDDigit.h <FMD/AliFMDDigit.h>
- @brief base class for digits
- @ingroup FMD_base
+/**
+ * @class AliFMDBaseDigit AliFMDDigit.h <FMD/AliFMDDigit.h>
+ *
+ * @brief base class for digits
+ *
+ * @ingroup FMD_base
*/
-class AliFMDBaseDigit : public TObject
+class AliFMDBaseDigit : public AliDigit
{
public:
- /** CTOR */
+ /**
+ * CTOR
+ */
AliFMDBaseDigit();
- /** Constrctor
- @param detector Detector
- @param ring Ring
- @param sector Sector
- @param strip Strip */
+ /**
+ * Constrctor
+ *
+ * @param detector Detector
+ * @param ring Ring
+ * @param sector Sector
+ * @param strip Strip
+ */
AliFMDBaseDigit(UShort_t detector,
Char_t ring='\0',
UShort_t sector=0,
UShort_t strip=0);
- /** DTOR */
+ /**
+ * Constrctor
+ *
+ * @param tracks Array of 3 track indicies
+ * @param detector Detector
+ * @param ring Ring
+ * @param sector Sector
+ * @param strip Strip
+ */
+ AliFMDBaseDigit(Int_t* tracks,
+ UShort_t detector,
+ Char_t ring='\0',
+ UShort_t sector=0,
+ UShort_t strip=0);
+ /**
+ * DTOR
+ */
virtual ~AliFMDBaseDigit() {}
- /** @return Detector # */
+ /**
+ *
+ * @return Detector #
+ */
UShort_t Detector() const { return fDetector; }
- /** @return Ring ID */
+ /**
+ *
+ * @return Ring ID
+ */
Char_t Ring() const { return fRing; }
- /** @return sector # */
+ /**
+ *
+ * @return sector #
+ */
UShort_t Sector() const { return fSector; }
- /** @return strip # */
+ /**
+ *
+ * @return strip #
+ */
UShort_t Strip() const { return fStrip; }
- /** Print information
- @param opt Not used */
+ /**
+ * Print information
+ *
+ * @param opt Not used
+ */
virtual void Print(Option_t* opt="") const;
- /** @return Name */
+ /**
+ *
+ * @return Name
+ */
const char* GetName() const;
- /** @param rhs Other digit to compare to
- @return -1 if this is less than @a rhs, 0 if the refer to the
- same, and 1 if @a rhs is larger than this */
+ /**
+ * @param rhs Other digit to compare to
+ *
+ * @return -1 if this is less than @a rhs, 0 if the refer to the
+ * same, and 1 if @a rhs is larger than this
+ */
Int_t Compare(const TObject* o) const;
- /** @return Always true */
+ /**
+ *
+ * @return Always true
+ */
Bool_t IsSortable() const { return kTRUE; }
+
+ /**
+ * Add a track referenc
+ *
+ * @param trackno The track number
+ */
+ void AddTrack(Int_t trackno);
+
+ /**
+ * Get the number of track references (max 3)
+ *
+ *
+ * @return Number of valid track references.
+ */
+ UShort_t GetNTrack() const;
+
+ /**
+ * Set the count value
+ *
+ * @param s Sample number
+ * @param c Counts
+ */
+ virtual void SetCount(UShort_t s, Short_t c) = 0;
protected:
+ /**
+ * Calculate the hash value
+ *
+ *
+ * @return Hash value
+ */
ULong_t Hash() const;
UShort_t fDetector; // (Sub) Detector # (1,2, or 3)
Char_t fRing; // Ring ID ('I' or 'O')
UShort_t fSector; // Sector # (phi division)
UShort_t fStrip; // Strip # (radial division)
mutable TString fName; //! Name (cached, but not stored)
- ClassDef(AliFMDBaseDigit, 2) // Base class for FMD digits
+ ClassDef(AliFMDBaseDigit, 3) // Base class for FMD digits
};
#endif
AliFMDMap::kMaxSectors,
AliFMDMap::kMaxStrips),
fShapingTime(6),
- fStoreTrackRefs(kFALSE)
+ fStoreTrackRefs(kTRUE)
{
AliFMDDebug(1, ("Constructed"));
// Default ctor - don't use it
AliFMDMap::kMaxSectors,
AliFMDMap::kMaxStrips),
fShapingTime(6),
- fStoreTrackRefs(kFALSE)
+ fStoreTrackRefs(kTRUE)
{
// Normal CTOR
AliFMDDebug(1, ("Constructed"));
AliFMDMap::kMaxSectors,
AliFMDMap::kMaxStrips),
fShapingTime(6),
- fStoreTrackRefs(kFALSE)
+ fStoreTrackRefs(kTRUE)
{
// Normal CTOR
AliFMDDebug(1, (" Constructed"));
UShort_t strip,
Float_t edep,
Bool_t isPrimary,
- Int_t trackno)
+ Int_t nTrack,
+ Int_t* tracknos)
{
// Add edep contribution from (detector,ring,sector,strip) to cache
AliFMDParameters* param = AliFMDParameters::Instance();
// Sum energy deposition
entry.fEdep += edep;
- entry.fN += 1;
- if (isPrimary) entry.fNPrim += 1;
+ entry.fN += nTrack;
+ if (isPrimary) entry.fNPrim += nTrack;
if (fStoreTrackRefs) {
+ Int_t oldN = entry.fLabels.fN;
entry.fLabels.Set(entry.fN);
- entry.fLabels[entry.fN-1] = trackno;
+ for (Int_t i = 0; i < nTrack; i++)
+ entry.fLabels[oldN + i] = tracknos[i];
}
- AliFMDDebug(15, ("Adding contribution %f to FMD%d%c[%2d,%3d] (%f)",
+ AliFMDDebug(15, ("Adding contribution %f to FMD%d%c[%2d,%3d] (%f) track %d",
edep, detector, ring, sector, strip,
- entry.fEdep));
+ entry.fEdep, (nTrack > 0 ? tracknos[0] : -1)));
}
Short_t(counts[2]), Short_t(counts[3]),
ntot, nprim, labels);
AliFMDDebug(15, (" Adding digit in FMD%d%c[%2d,%3d]=%d",
- detector,ring,sector,strip,counts[0]));
+ detector,ring,sector,strip,counts[0]));
#if 0
// This checks if the digit created will give the `right'
// number of particles when reconstructed, using a naiive
Short_t count4,
UShort_t /* ntot */,
UShort_t /* nprim */,
- const TArrayI& /* refs */) const
+ const TArrayI& refs) const
{
// Add a digit or summable digit
fFMD->AddDigitByFields(detector, ring, sector, strip,
- count1, count2, count3, count4);
+ count1, count2, count3, count4, refs);
}
//____________________________________________________________________
UShort_t strip,
Float_t edep,
Bool_t isPrimary,
- Int_t trackno);
+ Int_t nTrackno,
+ Int_t* tracknos);
/** Add a digit to output */
virtual void AddDigit(UShort_t detector,
Char_t ring,
}
//____________________________________________________________________
-AliFMDDigit::AliFMDDigit(UShort_t detector,
- Char_t ring,
- UShort_t sector,
- UShort_t strip,
- UShort_t count1,
- Short_t count2,
- Short_t count3,
- Short_t count4)
+AliFMDDigit::AliFMDDigit(UShort_t detector,
+ Char_t ring,
+ UShort_t sector,
+ UShort_t strip,
+ UShort_t count1,
+ Short_t count2,
+ Short_t count3,
+ Short_t count4,
+ const TArrayI& refs)
: AliFMDBaseDigit(detector, ring, sector, strip),
fCount1(count1),
fCount2(count2),
// count1 ADC count (a 10-bit word)
// count2 ADC count (a 10-bit word) -1 if not used
// count3 ADC count (a 10-bit word) -1 if not used
+ for (Int_t i = 0; i < refs.fN; i++) AddTrack(refs.fArray[i]);
}
//____________________________________________________________________
//____________________________________________________________________
void
-AliFMDDigit::Print(Option_t* /* option*/) const
+AliFMDDigit::Print(Option_t* option) const
{
// Print digit to standard out
AliFMDBaseDigit::Print();
- cout << "\t"
- << fCount1 << " (" << fCount2 << "," << fCount3 << "," << fCount4
- << ") = " << Counts() << endl;
+ std::cout << "\t"
+ << std::setw(4) << fCount1 << " ("
+ << std::setw(4) << fCount2 << ","
+ << std::setw(4) << fCount3 << ","
+ << std::setw(4) << fCount4 << ") = "
+ << std::setw(4) << Counts() << std::flush;
+ TString opt(option);
+ if (opt.Contains("l", TString::kIgnoreCase)) {
+ std::cout << " ";
+ for (Int_t i = 0; i < 3; i++)
+ std::cout << (i == 0 ? "" : ",") << std::setw(5) << fTracks[i];
+ }
+ std::cout << std::endl;
}
//____________________________________________________________________
#ifndef ALIFMDBASEDIGIT_H
# include <AliFMDBaseDigit.h>
#endif
+#ifndef ROOT_TArrayI
+# include <TArrayI.h>
+#endif
+
//____________________________________________________________________
/** @class AliFMDDigit AliFMDDigit.h <FMD/AliFMDDigit.h>
public:
/** CTOR */
AliFMDDigit();
- /** Constrctor
- @param detector Detector
- @param ring Ring
- @param sector Sector
- @param strip Strip
- @param count ADC (first sample)
- @param count2 ADC (second sample, or -1 if not used)
- @param count3 ADC (third sample, or -1 if not used) */
- AliFMDDigit(UShort_t detector,
- Char_t ring='\0',
- UShort_t sector=0,
- UShort_t strip=0,
- UShort_t count=0,
- Short_t count2=-1,
- Short_t count3=-1,
- Short_t count4=-1);
- /** DTOR */
+ /**
+ * Constrctor
+ *
+ * @param detector Detector
+ * @param ring Ring
+ * @param sector Sector
+ * @param strip Strip
+ * @param count ADC (first sample)
+ * @param count2 ADC (second sample, or -1 if not used)
+ * @param count3 ADC (third sample, or -1 if not used)
+ * @param refs Track references
+ */
+ AliFMDDigit(UShort_t detector,
+ Char_t ring='\0',
+ UShort_t sector=0,
+ UShort_t strip=0,
+ UShort_t count=0,
+ Short_t count2=-1,
+ Short_t count3=-1,
+ Short_t count4=-1,
+ const TArrayI& refs=TArrayI());
+ /**
+ * DTOR
+ */
virtual ~AliFMDDigit() {}
- /** @param i # of sample to get
- @return sample # @a i */
+ /**
+ * @param i # of sample to get
+ *
+ * @return sample # @a i
+ */
Int_t Count(UShort_t i=0) const;
- /** @return ADC count (first sample) */
- UShort_t Count1() const { return fCount1; }
- /** @return ADC count (second sample, or -1 if not used) */
- Short_t Count2() const { return fCount2; }
- /** @return ADC count (third sample, or -1 if not used) */
- Short_t Count3() const { return fCount3; }
- /** @return ADC count (third sample, or -1 if not used) */
- Short_t Count4() const { return fCount4; }
- /** @return Canonical ADC counts */
- UShort_t Counts() const;
- /** Print info
- @param opt Not used */
+ /**
+ *
+ * @return ADC count (first sample)
+ */
+ UShort_t Count1() const { return fCount1; }
+ /**
+ *
+ * @return ADC count (second sample, or -1 if not used)
+ */
+ Short_t Count2() const { return fCount2; }
+ /**
+ *
+ * @return ADC count (third sample, or -1 if not used)
+ */
+ Short_t Count3() const { return fCount3; }
+ /**
+ *
+ * @return ADC count (third sample, or -1 if not used)
+ */
+ Short_t Count4() const { return fCount4; }
+ /**
+ *
+ * @return Canonical ADC counts
+ */
+ UShort_t Counts() const;
+ /**
+ * Print info
+ *
+ * @param opt Not used
+ */
void Print(Option_t* opt="") const;
- /** @return Title */
+ /**
+ *
+ * @return Title
+ */
const char* GetTitle() const;
- /** Set the count value
- @param s Sample number
- @param c Counts */
+ /**
+ * Set the count value
+ *
+ * @param s Sample number
+ * @param c Counts
+ */
void SetCount(UShort_t s, Short_t c);
protected:
UShort_t fCount1; // Digital signal
folderName.Data()));
return;
}
- runLoader->LoadgAlice();
+ if (!runLoader->GetAliRun()) runLoader->LoadgAlice();
runLoader->GetAliRun();
}
if (!gAlice) {
fFMD->SetTreeAddress();
// Sum contributions from the sdigits
+ AliFMDDebug(3, ("Will now sum contributions from SDigits"));
SumContributions(sdigitsBranch);
// Unload the sdigits
// Get the sdigit number `sdigit'
AliFMDSDigit* fmdSDigit =
static_cast<AliFMDSDigit*>(fmdSDigits->UncheckedAt(sdigit));
+
+ AliFMDDebug(5, ("Adding contribution of %d tracks",
+ fmdSDigit->GetNTrack()));
+ AliFMDDebug(15, ("Contrib from FMD%d%c[%2d,%3d] (%s) from track %d",
+ fmdSDigit->Detector(),
+ fmdSDigit->Ring(),
+ fmdSDigit->Sector(),
+ fmdSDigit->Strip(),
+ fmdSDigit->GetName(),
+ fmdSDigit->GetTrack(0)));
// Extract parameters
AddContribution(fmdSDigit->Detector(),
fmdSDigit->Strip(),
fmdSDigit->Edep(),
kTRUE,
- -1);
+ fmdSDigit->GetNTrack(),
+ fmdSDigit->GetTracks());
} // sdigit loop
} // event loop
AliFMDRing* r = fmd1->GetInner();
Double_t z = fmd1->GetInnerZ();
+ // `Top' or `Outside' master volume
TGeoVolume* fmd1TopVolume = new TGeoVolumeAssembly(Form(fgkFMDName,
fmd1->GetId(), 'T'));
fmd1TopVolume->SetTitle("FMD1 top half");
+
+ // `Bottom' or `Inside' master volume
TGeoVolume* fmd1BotVolume = new TGeoVolumeAssembly(Form(fgkFMDName,
fmd1->GetId(), 'B'));
fmd1BotVolume->SetTitle("FMD1 bottom half");
TGeoRotation* rot = new TGeoRotation("FMD1 rotatation");
rot->RotateZ(90);
TGeoMatrix* matrix = new TGeoCombiTrans("FMD1 trans", 0, 0, z, rot);
+
AliFMDDebug(5, ("Placing volumes %s and %s in ALIC at z=%f",
fmd1TopVolume->GetName(), fmd1BotVolume->GetName(), z));
top->AddNode(fmd1TopVolume, fmd1->GetId(), matrix);
top->AddNode(fmd1BotVolume, fmd1->GetId(), matrix);
+
+
+ // Survey points on V0A (screw holes for the FMD)
+ const Double_t icb[] = { +12.700, -21.997, 324.670 };
+ const Double_t ict[] = { +12.700, +21.997, 324.670 };
+ const Double_t ocb[] = { -12.700, -21.997, 324.670 };
+ const Double_t oct[] = { -12.700, +21.997, 324.670 };
+
+ TGeoTube* surveyShape = new TGeoTube("FMD1_survey_marker",
+ 0, .2, .001);
+
+ TGeoMatrix* outMat = matrix;
+#if 0
+ if (gGeoManager->cd("/ALIC_1/F1MT_1"))
+ outMat = gGeoManager->GetCurrentMatrix();
+ else
+ AliWarning("Couldn't cd to /ALIC_1/F1MT_1");
+#endif
+
+ Double_t loct[3], locb[3];
+ outMat->MasterToLocal(oct, loct);
+ outMat->MasterToLocal(ocb, locb);
+ TGeoVolume* vOct = new TGeoVolume("V0L_OCT", surveyShape, fPlastic);
+ TGeoVolume* vOcb = new TGeoVolume("V0L_OCB", surveyShape, fPlastic);
+
+ fmd1TopVolume->AddNode(vOct, 1, new TGeoTranslation(loct[0],loct[1],loct[2]));
+ fmd1TopVolume->AddNode(vOcb, 1, new TGeoTranslation(locb[0],locb[1],locb[2]));
+
+
+ TGeoMatrix* inMat = matrix;
+#if 0
+ if (gGeoManager->cd("/ALIC_1/F1MT_1"))
+ inMat = gGeoManager->GetCurrentMatrix();
+ else
+ AliWarning("Couldn't cd to /ALIC_1/F1MT_1");
+#endif
+
+ Double_t lict[3], licb[3];
+ inMat->MasterToLocal(ict, lict);
+ inMat->MasterToLocal(icb, licb);
+ TGeoVolume* vIct = new TGeoVolume("V0L_ICT", surveyShape, fPlastic);
+ TGeoVolume* vIcb = new TGeoVolume("V0L_ICB", surveyShape, fPlastic);
+ fmd1BotVolume->AddNode(vIct, 1, new TGeoTranslation(lict[0],lict[1],lict[2]));
+ fmd1BotVolume->AddNode(vIcb, 1, new TGeoTranslation(licb[0],licb[1],licb[2]));
+
return 0;
}
fmdHit->Strip(),
fmdHit->Edep(),
isPrimary,
- trackno);
+ 1,
+ &trackno);
} // hit loop
} // track loop
AliFMDDebug(5, ("Size of cache: %d bytes, read %d bytes",
{
// Add a digit or summable digit
if (fOutput == kDigits) {
+ AliFMDDebug(15,("Adding digit for FMD%d%c[%2d,%3d] = (%x,%x,%x,%x)",
+ detector, ring, sector, strip,
+ count1, count2, count3, count4));
AliFMDBaseDigitizer::AddDigit(detector, ring, sector, strip, 0,
count1, count2, count3, count4, 0, 0, refs);
return;
detector, ring, sector, strip));
return;
}
- AliFMDDebug(15, ("Adding digit for FMD%d%c[%2d,%3d] = (%x,%x,%x,%x) [%d/%d]",
+ AliFMDDebug(15, ("Adding sdigit for FMD%d%c[%2d,%3d] = "
+ "(%x,%x,%x,%x) [%d/%d] %d",
detector, ring, sector, strip,
- count1, count2, count3, count4, nprim, ntotal));
+ count1, count2, count3, count4, nprim, ntotal, refs.fN));
fFMD->AddSDigitByFields(detector, ring, sector, strip, edep,
count1, count2, count3, count4,
ntotal, nprim, refs);
fChainE(0),
fArrayE(0),
fArrayH(0),
+ fArrayTR(0),
fArrayD(0),
fArrayS(0),
fArrayR(0),
fChainE(0),
fArrayE(0),
fArrayH(0),
+ fArrayTR(0),
fArrayD(0),
fArrayS(0),
fArrayR(0),
{
// Get number of events
if (TESTBIT(fTreeMask, kRaw)) return fReader->GetNumberOfEvents();
+ if (fChainE) return fChainE->GetEntriesFast();
if (fTreeE) return fTreeE->GetEntries();
return -1;
}
AliWarning("Already initialized");
return fIsInit;
}
+ Info("Init","Initialising w/mask 0x%04x\n"
+ "\tHits: %d\n"
+ "\tKinematics: %d\n"
+ "\tDigits: %d\n"
+ "\tSDigits: %d\n"
+ "\tHeader: %d\n"
+ "\tRecPoints: %d\n"
+ "\tESD: %d\n"
+ "\tRaw: %d\n"
+ "\tGeometry: %d\n"
+ "\tTracks: %d\n"
+ "\tTracksRefs: %d",
+ fTreeMask,
+ TESTBIT(fTreeMask, kHits),
+ TESTBIT(fTreeMask, kKinematics),
+ TESTBIT(fTreeMask, kDigits),
+ TESTBIT(fTreeMask, kSDigits),
+ TESTBIT(fTreeMask, kHeader),
+ TESTBIT(fTreeMask, kRecPoints),
+ TESTBIT(fTreeMask, kESD),
+ TESTBIT(fTreeMask, kRaw),
+ TESTBIT(fTreeMask, kGeometry),
+ TESTBIT(fTreeMask, kTracks),
+ TESTBIT(fTreeMask, kTrackRefs));
// Get the run
- if (!TESTBIT(fTreeMask, kRaw) &&
- !TESTBIT(fTreeMask, kESD)) {
- if (fGAliceFile.IsNull()) fGAliceFile = "galice.root";
- // Get the loader
- fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read");
- if (!fLoader) {
- AliError(Form("Coulnd't read the file %s", fGAliceFile.Data()));
- return kFALSE;
+ if (TESTBIT(fTreeMask, kDigits) ||
+ TESTBIT(fTreeMask, kSDigits) ||
+ TESTBIT(fTreeMask, kKinematics) ||
+ TESTBIT(fTreeMask, kTrackRefs) ||
+ TESTBIT(fTreeMask, kTracks) ||
+ TESTBIT(fTreeMask, kHeader)) {
+ if (!gSystem->FindFile(".:/", fGAliceFile)) {
+ AliWarning(Form("Cannot find file %s in .:/", fGAliceFile.Data()));
}
-
- if (fLoader->LoadgAlice()) return kFALSE;
- fRun = fLoader->GetAliRun();
-
- // Get the FMD
- fFMD = static_cast<AliFMD*>(fRun->GetDetector("FMD"));
- if (!fFMD) {
- AliError("Failed to get detector FMD from loader");
- return kFALSE;
- }
-
- // Get the FMD loader
- fFMDLoader = fLoader->GetLoader("FMDLoader");
- if (!fFMDLoader) {
- AliError("Failed to get detector FMD loader from loader");
- return kFALSE;
- }
- if (fLoader->LoadHeader()) {
- AliError("Failed to get event header information from loader");
- return kFALSE;
+ else {
+ fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read");
+ if (!fLoader) {
+ AliError(Form("Coulnd't read the file %s", fGAliceFile.Data()));
+ return kFALSE;
+ }
+ AliInfo(Form("Opened GAlice file %s", fGAliceFile.Data()));
+
+ if (fLoader->LoadgAlice()) return kFALSE;
+
+ fRun = fLoader->GetAliRun();
+
+ // Get the FMD
+ fFMD = static_cast<AliFMD*>(fRun->GetDetector("FMD"));
+ if (!fFMD) {
+ AliError("Failed to get detector FMD from loader");
+ return kFALSE;
+ }
+
+ // Get the FMD loader
+ fFMDLoader = fLoader->GetLoader("FMDLoader");
+ if (!fFMDLoader) {
+ AliError("Failed to get detector FMD loader from loader");
+ return kFALSE;
+ }
+ if (fLoader->LoadHeader()) {
+ AliError("Failed to get event header information from loader");
+ return kFALSE;
+ }
+ fTreeE = fLoader->TreeE();
}
- fTreeE = fLoader->TreeE();
}
// Optionally, get the ESD files
// Possibly load FMD Sdigit information
if (TESTBIT(fTreeMask, kSDigits)) {
// AliInfo("Getting FMD summable digits");
- if (!fFMDLoader || fFMDLoader->LoadSDigits("READ")) return kFALSE;
+ if (!fFMDLoader || fFMDLoader->LoadSDigits("READ")) {
+ AliWarning("Failed to load SDigits!");
+ return kFALSE;
+ }
fTreeS = fFMDLoader->TreeS();
if (!fArrayS) fArrayS = fFMD->SDigits();
}
if (!ProcessTrackRefs()) return kFALSE;
if (TESTBIT(fTreeMask, kTracks))
if (!ProcessTracks()) return kFALSE;
- if (TESTBIT(fTreeMask, kDigits))
- if (!ProcessDigits()) return kFALSE;
if (TESTBIT(fTreeMask, kSDigits))
if (!ProcessSDigits()) return kFALSE;
+ if (TESTBIT(fTreeMask, kDigits))
+ if (!ProcessDigits()) return kFALSE;
if (TESTBIT(fTreeMask, kRaw))
if (!ProcessRawDigits()) return kFALSE;
if (TESTBIT(fTreeMask, kRecPoints))
Int_t nEv = fTreeS->GetEntries();
for (Int_t i = 0; i < nEv; i++) {
Int_t sdigitRead = fTreeS->GetEntry(i);
- if (sdigitRead <= 0) continue;
- Int_t nSdigit = fArrayS->GetEntries();
+ if (sdigitRead <= 0) {
+ AliInfo(Form("Read nothing from tree"));
+ continue;
+ }
+ Int_t nSdigit = fArrayS->GetEntriesFast();
AliFMDDebug(0, ("Got %5d digits for this event", nSdigit));
+ AliInfo(Form("Got %5d digits for this event", nSdigit));
if (nSdigit <= 0) continue;
for (Int_t j = 0; j < nSdigit; j++) {
AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(fArrayS->At(j));
};
inline Bool_t AliFMDInput::ProcessHit(AliFMDHit*,TParticle*) { return kTRUE; }
-inline Bool_t AliFMDInput::ProcessTrackRef(AliTrackReference* trackRef, TParticle* track) { return kTRUE; }
+inline Bool_t AliFMDInput::ProcessTrackRef(AliTrackReference*,
+ TParticle*) { return kTRUE; }
inline Bool_t AliFMDInput::ProcessTrack(Int_t,TParticle*,
AliFMDHit*) { return kTRUE; }
inline Bool_t AliFMDInput::ProcessDigit(AliFMDDigit*) { return kTRUE; }
#include "AliFMDDebug.h" // Better debug macros
#include "AliFMDParameters.h" // ALIFMDPARAMETERS_H
#include "AliFMDDigit.h" // ALIFMDDIGIT_H
+#include "AliFMDSDigit.h" // ALIFMDSDIGIT_H
#include "AliFMDRawStream.h" // ALIFMDRAWSTREAM_H
#include "AliRawReader.h" // ALIRAWREADER_H
#include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
: TTask("FMDRawReader", "Reader of Raw ADC values from the FMD"),
fTree(tree),
fReader(reader),
- fSampleRate(1),
+ // fSampleRate(1),
fData(0),
fNbytes(0),
fSeen()
{
// Default CTOR
+ for (Int_t i = 0; i < 3; i++) {
+ fSampleRate[i] = 0;
+ fZeroSuppress[i] = kFALSE;
+ fNoiseFactor[i] = 1;
+ }
}
//____________________________________________________________________
array->GetEntriesFast(), nWrite));
}
-
//____________________________________________________________________
Bool_t
-AliFMDRawReader::NextSignal(UShort_t& det, Char_t& rng,
- UShort_t& sec, UShort_t& str,
- Short_t& adc, Bool_t& zs,
- UShort_t& fac)
+AliFMDRawReader::NextSample(UShort_t& det, Char_t& rng, UShort_t& sec,
+ UShort_t& str, UShort_t& sam, UShort_t& rat,
+ Short_t& adc, Bool_t& zs, UShort_t& fac)
{
// Scan current event for next signal. It returns kFALSE when
// there's no more data in the event.
static AliFMDParameters* pars = 0;
static AliFMDAltroMapping* map = 0;
static Int_t ddl = -1;
- static UInt_t rate = 0;
+ // static UInt_t rate = 0;
static UShort_t tdet = 0;
static Char_t trng = '\0';
static UShort_t tsec = 0;
static Short_t tstr = 0;
static Short_t bstr = -1;
+ static Short_t tsam = -1;
+ static UInt_t trate = 0;
static Int_t hwaddr = -1;
static UShort_t stripMin = 0;
static UShort_t stripMax = 0; // 127;
// Reset variables
ddl = -1;
- rate = 0;
+ trate = 0;
tdet = 0;
trng = '\0';
tsec = 0;
tstr = 0;
+ tsam = -1;
hwaddr = -1;
}
do {
Int_t thisDDL = stream.GetDDLNumber();
AliFMDDebug(10, ("RCU @ DDL %d", thisDDL));
if (thisDDL != ddl) {
- ddl = thisDDL;
- zs = stream.GetZeroSupp();
- fac = stream.GetNPostsamples();
- tdet = map->DDL2Detector(ddl);
- rate = 0; // stream.GetNPresamples();
+ ddl = thisDDL;
+ fZeroSuppress[ddl] = zs = stream.GetZeroSupp();
+ fNoiseFactor[ddl] = fac = stream.GetNPostsamples();
+ fSampleRate[ddl] = trate = 0; // stream.GetNPresamples();
+ tdet = map->DDL2Detector(ddl);
AliFMDDebug(10, ("RCU @ DDL %d zero suppression: %s",ddl, zs?"yes":"no"));
AliFMDDebug(10, ("RCU @ DDL %d noise factor: %d", ddl,fac));
- AliFMDDebug(10, ("RCU @ DDL %d sample rate: %d", ddl, rate));
+ AliFMDDebug(10, ("RCU @ DDL %d sample rate: %d", ddl, trate));
}
Int_t thisAddr = stream.GetHWAddress();
AliFMDDebug(10, ("RCU @ DDL %d, Address 0x%03x", ddl, thisAddr));
stripMin = pars->GetMinStrip(tdet, trng, tsec, bstr);
stripMax = pars->GetMaxStrip(tdet, trng, tsec, bstr);
preSamp = pars->GetPreSamples(tdet, trng, tsec, bstr);
- if (rate == 0) rate = pars->GetSampleRate(tdet, trng, tsec, bstr);
+ if (trate == 0)
+ fSampleRate[ddl] = trate = pars->GetSampleRate(tdet, trng, tsec, bstr);
AliFMDDebug(10, ("RCU @ DDL %d, Address 0x%03x sample rate: %d",
- ddl, hwaddr, rate));
+ ddl, hwaddr, trate));
Int_t nChAddrMismatch = stream.GetNChAddrMismatch();
Int_t nChLenMismatch = stream.GetNChLengthMismatch();
Short_t strOff = 0;
UShort_t samp = 0;
- map->Timebin2Strip(tsec, t, preSamp, rate, strOff, samp);
+ map->Timebin2Strip(tsec, t, preSamp, trate, strOff, samp);
tstr = bstr + strOff;
+ tsam = samp;
AliFMDDebug(20, ("0x%04x/0x%03x/%04d maps to FMD%d%c[%2d,%3d]-%d",
ddl, hwaddr, t, tdet, trng, tsec, tstr, samp));
- Bool_t take = kFALSE;
- switch (rate) {
- case 1: take = kTRUE; break;
- case 2: if (samp == 1) take = kTRUE; break;
- case 3: if (samp == 1) take = kTRUE; break;
- case 4: if (samp == 2) take = kTRUE; break;
- default: if (samp == rate-2) take = kTRUE; break;
- }
- if (!take) continue;
-
-
// Local strip number to channel
Short_t l = (tstr > bstr ? tstr - bstr : bstr - tstr);
AliFMDDebug(10, ("Checking if strip %d in range [%d,%d]",
rng = trng;
sec = tsec;
str = tstr;
+ sam = tsam;
+ rat = trate;
// adc = stream.GetSignal();
break;
return kTRUE;
}
+//____________________________________________________________________
+Bool_t
+AliFMDRawReader::NextSignal(UShort_t& det, Char_t& rng,
+ UShort_t& sec, UShort_t& str,
+ Short_t& adc, Bool_t& zs,
+ UShort_t& fac)
+{
+
+ do {
+ UShort_t samp, rate;
+ if (!NextSample(det, rng, sec, str, samp, rate, adc, zs, fac))
+ return kFALSE;
+
+ Bool_t take = kFALSE;
+ switch (rate) {
+ case 1: take = kTRUE; break;
+ case 2: if (samp == 1) take = kTRUE; break;
+ case 3: if (samp == 1) take = kTRUE; break;
+ case 4: if (samp == 2) take = kTRUE; break;
+ default: if (samp == rate-2) take = kTRUE; break;
+ }
+ if (!take) continue;
+ break;
+ } while (true);
+ return kTRUE;
+}
+
+//____________________________________________________________________
+Bool_t
+AliFMDRawReader::SelectSample(UShort_t samp, UShort_t rate)
+{
+ Bool_t take = kFALSE;
+ switch (rate) {
+ case 1: take = kTRUE; break;
+ case 2: if (samp == 1) take = kTRUE; break;
+ case 3: if (samp == 1) take = kTRUE; break;
+ case 4: if (samp == 2) take = kTRUE; break;
+ default: if (samp == rate-2) take = kTRUE; break;
+ }
+
+ return take;
+}
+
#if 1
//____________________________________________________________________
Bool_t
Int_t oldddl = -1;
UInt_t ddl = 0;
- UInt_t rate = 0;
+ // UInt_t rate = 0;
UInt_t last = 0;
UInt_t hwaddr = 0;
// Data array is approx twice the size needed.
fNoiseFactor[ddl] = input.GetNPostsamples();
AliFMDDebug(20, ("RCU @ DDL %d noise factor: %d", ddl,fNoiseFactor[ddl]));
+ // WARNING: We store the noise factor in the 2nd baseline
+ // filters excluded post samples, since we'll never use that
+ // mode.
+ fSampleRate[ddl] = input.GetNPretriggerSamples();
+ AliFMDDebug(20, ("RCU @ DDL %d Sample rate: %d", ddl,fNoiseFactor[ddl]));
+
Int_t nChAddrMismatch = input.GetNChAddrMismatch();
Int_t nChLenMismatch = input.GetNChLengthMismatch();
if (nChAddrMismatch != 0)
"hardware address 0x%03x", ddl, hwaddr));
continue;
}
- AliFMDDebug(1, ("Board: 0x%02x, Altro: 0x%x, Channel: 0x%x, Length: %4d",
+ AliFMDDebug(5, ("Board: 0x%02x, Altro: 0x%x, Channel: 0x%x, Length: %4d",
board, chip, channel, last));
stripMin = pars->GetMinStrip(det, ring, sec, strbase);
stripMax = pars->GetMaxStrip(det, ring, sec, strbase);
preSamp = pars->GetPreSamples(det, ring, sec, strbase);
- // WARNING: We use the number of pre-samples to store the
- // oversampling rate in.
- rate = input.GetNPretriggerSamples();
- if (rate == 0) rate = pars->GetSampleRate(det, ring, sec, strbase);
+ if (fSampleRate[ddl] == 0)
+ fSampleRate[ddl] = pars->GetSampleRate(det, ring, sec, strbase);
// Loop over the `timebins', and make the digits
for (size_t i = 0; i < last; i++) {
AliFMDDebug(15, ("0x%04x/0x%03x/%04d %4d", ddl, hwaddr, i, data[i]));
Short_t stroff = 0;
- map->Timebin2Strip(sec, i, preSamp, rate, stroff, samp);
+ map->Timebin2Strip(sec, i, preSamp, fSampleRate[ddl], stroff, samp);
Short_t str = strbase + stroff;
AliFMDDebug(10, ("0x%04x/0x%03x/%04d maps to FMD%d%c[%2d,%3d]-%d",
continue;
}
- Short_t lstrip = (i - preSamp) / rate + stripMin;
+ Short_t lstrip = (i - preSamp) / fSampleRate[ddl] + stripMin;
AliFMDDebug(15, ("Checking if strip %d (%d) in range [%d,%d]",
lstrip, str, stripMin, stripMax));
data[i] = 0; // Reset cache
continue;
}
+ // Possibly do pedestal subtraction of signal
+ Int_t counts = data[i];
+
// Check the cache of indicies
Int_t idx = fSeen(det, ring, sec, str);
det, ring, sec, str, samp, i));
new ((*array)[idx]) AliFMDDigit(det, ring, sec, str);
}
- AliFMDDigit* digit = static_cast<AliFMDDigit*>(array->At(idx));
+ AliFMDBaseDigit* digit = static_cast<AliFMDBaseDigit*>(array->At(idx));
AliFMDDebug(10,
- ("Setting from FMD%d%c[%2d,%3d]-%d from timebin %4d = %4d",
- det, ring, sec, str, samp, i, data[i]));
- digit->SetCount(samp, data[i]);
+ ("Setting FMD%d%c[%2d,%3d]-%d from timebin %4d=%4d (%4d)",
+ det, ring, sec, str, samp, i, counts, data[i]));
+ digit->SetCount(samp, counts);
data[i] = 0; // Reset cache
}
}
class AliFMDRawReader : public TTask
{
public:
- /** CTOR
- @param reader Raw reader
- @param array Output tree */
+ /**
+ * CTOR
+ *
+ * @param reader Raw reader
+ * @param array Output tree
+ */
AliFMDRawReader(AliRawReader* reader, TTree* array);
- /** DTOR */
+ /**
+ * DTOR
+ */
virtual ~AliFMDRawReader() {}
- /** Read in, and store in output tree
- @param option Not used */
+ /**
+ * Read in, and store in output tree
+ *
+ * @param option Not used
+ */
virtual void Exec(Option_t* option="");
- /** Read ADC's into a TClonesArray of AliFMDDigit objects.
- @param array Array to read into
- @return @c true on success */
+ /**
+ * Read ADC's into a TClonesArray of AliFMDDigit objects.
+ *
+ * @param array Array to read into
+ * @param summable Create SDigits rather than digits
+ * @param pedSub Whether to do pedestal subtraction.
+ * @param noiseFactor If we do pedestal subtraction, then this is
+ * the number we use to suppress remenants of the noise.
+ *
+ * @return @c true on success
+ */
virtual Bool_t ReadAdcs(TClonesArray* array);
/**
* Read SOD event into passed objects.
*/
UShort_t NoiseFactor(UShort_t ddl) const { return fNoiseFactor[ddl]; }
+ /**
+ * Get the next signal
+ *
+ * @param det On return, the detector
+ * @param rng On return, the ring
+ * @param sec On return, the sector
+ * @param str On return, the strip
+ * @param sam On return, the sample
+ * @param rat On return, the sample rate
+ * @param adc On return, the ADC value
+ * @param zs On return, whether zero-supp. is enabled
+ * @param fac On return, the usd noise factor
+ *
+ * @return true if valid data is returned
+ */
+ Bool_t NextSample(UShort_t& det, Char_t& rng, UShort_t& sec, UShort_t& str,
+ UShort_t& sam, UShort_t& rat, Short_t& adc,
+ Bool_t& zs, UShort_t& fac);
+ /**
+ * Get the next signal
+ *
+ * @param det On return, the detector
+ * @param rng On return, the ring
+ * @param sec On return, the sector
+ * @param str On return, the strip
+ * @param adc On return, the ADC value
+ * @param zs On return, whether zero-supp. is enabled
+ * @param fac On return, the usd noise factor
+ *
+ * @return true if valid data is returned
+ */
Bool_t NextSignal(UShort_t& det, Char_t& rng,
UShort_t& sec, UShort_t& str,
Short_t& adc, Bool_t& zs,
UShort_t& fac);
+ /**
+ * Whether to keep a sample based on the rate used.
+ *
+ * @param samp Sample number
+ * @param rate Over sampling rate
+ *
+ * @return Whether to keep the sample or not
+ */
+ static Bool_t SelectSample(UShort_t samp, UShort_t rate);
protected:
AliFMDRawReader(const AliFMDRawReader& o)
: TTask(o),
fTree(0),
fReader(0),
- fSampleRate(0),
fData(0),
fNbytes(0),
fSeen()
Int_t GetHalfringIndex(UShort_t det, Char_t ring, UShort_t board);
TTree* fTree; //! Pointer to tree to read into
AliRawReader* fReader; //! Pointer to raw reader
- UShort_t fSampleRate; // The sample rate (if 0, inferred from data)
+ UShort_t fSampleRate[3]; // The sample rate (if 0, inferred from data)
UChar_t* fData;
ULong_t fNbytes;
Bool_t fZeroSuppress[3];
#include "AliFMDParameters.h" // ALIFMDPARAMETERS_H
#include "AliFMDAltroMapping.h" // ALIFMDALTROMAPPING_H
#include "AliFMDDigit.h" // ALIFMDDIGIT_H
+#include "AliFMDSDigit.h" // ALIFMDDIGIT_H
#include "AliFMDReconstructor.h" // ALIFMDRECONSTRUCTOR_H
#include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
#include "AliFMDRecPoint.h" // ALIFMDMULTNAIIVE_H
TTree* digitsTree) const
{
// Convert Raw digits to AliFMDDigit's in a tree
- AliFMDDebug(2, ("Reading raw data into digits tree"));
+ AliFMDDebug(1, ("Reading raw data into digits tree"));
AliFMDRawReader rawRead(reader, digitsTree);
// rawRead.SetSampleRate(fFMD->GetSampleRate());
rawRead.Exec();
}
}
+//____________________________________________________________________
+void
+AliFMDReconstructor::Digitize(AliRawReader* reader, TClonesArray* sdigits) const
+{
+ // Reconstruct directly from raw data (no intermediate output on
+ // digit tree or rec point tree).
+ //
+ // Parameters:
+ // reader Raw event reader
+ // ctree Not used.
+ AliFMDRawReader rawReader(reader, 0);
+
+ UShort_t det, sec, str, sam, rat, fac;
+ Short_t adc, oldDet = -1;
+ Bool_t zs;
+ Char_t rng;
+
+ while (rawReader.NextSample(det, rng, sec, str, sam, rat, adc, zs, fac)) {
+ if (!rawReader.SelectSample(sam, rat)) continue;
+ if (det != oldDet) {
+ fZS[det-1] = zs;
+ fZSFactor[det-1] = fac;
+ oldDet = det;
+ }
+ DigitizeSignal(sdigits, det, rng, sec, str, sam, adc);
+ }
+}
+
//____________________________________________________________________
void
AliFMDReconstructor::Reconstruct(TTree* digitsTree,
AliFMDDebug(5, ("Getting entry 0 from digit branch"));
digitBranch->GetEntry(0);
- AliFMDDebug(5, ("Processing digits"));
+ AliFMDDebug(1, ("Processing digits"));
ProcessDigits(digits);
Int_t written = clusterTree->Fill();
}
+//____________________________________________________________________
+void
+AliFMDReconstructor::DigitizeSignal(TClonesArray* sdigits,
+ UShort_t det,
+ Char_t rng,
+ UShort_t sec,
+ UShort_t str,
+ UShort_t /* sam */,
+ Short_t adc) const
+{
+ // Process the signal from a single strip
+ //
+ // Parameters:
+ // det Detector ID
+ // rng Ring ID
+ // sec Sector ID
+ // rng Strip ID
+ // adc ADC counts
+ //
+ AliFMDParameters* param = AliFMDParameters::Instance();
+ // Check that the strip is not marked as dead
+ if (param->IsDead(det, rng, sec, str)) {
+ AliFMDDebug(10, ("FMD%d%c[%2d,%3d] is dead", det, rng, sec, str));
+ return;
+ }
+
+ // Substract pedestal.
+ UShort_t counts = SubtractPedestal(det, rng, sec, str, adc);
+ if(counts == USHRT_MAX || counts == 0) return;
+
+ // Gain match digits.
+ Double_t edep = Adc2Energy(det, rng, sec, str, counts);
+ // Get rid of nonsense energy
+ if(edep < 0) return;
+
+ Int_t n = sdigits->GetEntriesFast();
+ // AliFMDSDigit* sdigit =
+ new ((*sdigits)[n])
+ AliFMDSDigit(det, rng, sec, str, edep, counts, counts, counts, counts);
+ // sdigit->SetCount(sam, counts);
+}
+
+//____________________________________________________________________
+UShort_t
+AliFMDReconstructor::SubtractPedestal(UShort_t det,
+ Char_t rng,
+ UShort_t sec,
+ UShort_t str,
+ UShort_t adc,
+ Float_t noiseFactor,
+ Bool_t zsEnabled,
+ UShort_t zsNoiseFactor) const
+{
+ AliFMDParameters* param = AliFMDParameters::Instance();
+ Float_t ped = (zsEnabled ? 0 :
+ param->GetPedestal(det, rng, sec, str));
+ Float_t noise = param->GetPedestalWidth(det, rng, sec, str);
+ if(ped < 0 || noise < 0) {
+ AliWarningClass(Form("Invalid pedestal (%f) or noise (%f) "
+ "for FMD%d%c[%02d,%03d]",
+ ped, noise, det, rng, sec, str));
+ return USHRT_MAX;
+ }
+ AliDebugClass(15, Form("Subtracting pedestal for FMD%d%c[%2d,%3d]=%4d "
+ "(%s w/factor %d, noise factor %f, "
+ "pedestal %8.2f+/-%8.2f)",
+ det, rng, sec, str, adc,
+ (zsEnabled ? "zs'ed" : "straight"),
+ zsNoiseFactor, noiseFactor, ped, noise));
+
+ Int_t counts = adc + Int_t(zsEnabled ? zsNoiseFactor * noise : - ped);
+ counts = TMath::Max(Int_t(counts), 0);
+ // Calculate the noise factor for suppressing remenants of the noise
+ // peak. If we have done on-line zero suppression, we only check
+ // for noise signals that are larger than the suppressed noise. If
+ // the noise factor used on line is larger than the factor used
+ // here, we do not do this check at all.
+ //
+ // For example:
+ // Online factor | Read factor | Result
+ // ---------------+--------------+-------------------------------
+ // 2 | 3 | Check if signal > 1 * noise
+ // 3 | 3 | Check if signal > 0
+ // 3 | 2 | Check if signal > 0
+ //
+ // In this way, we make sure that we do not suppress away too much
+ // data, and that the read-factor is the most stringent cut.
+ Float_t nf = TMath::Max(0.F, noiseFactor - (zsEnabled ? zsNoiseFactor : 0));
+ if (counts < noise * nf) counts = 0;
+ if (counts > 0) AliDebugClass(15, "Got a hit strip");
+
+ return counts;
+}
//____________________________________________________________________
// Return:
// Pedestal subtracted signal or USHRT_MAX in case of problems
//
- AliFMDParameters* param = AliFMDParameters::Instance();
- Bool_t zs = fZS[det-1];
- UShort_t fac = fZSFactor[det-1];
- Float_t ped = (zs ? 0 :
- param->GetPedestal(det, rng, sec, str));
- Float_t noise = param->GetPedestalWidth(det, rng, sec, str);
- if(ped < 0 || noise < 0) {
- AliWarning(Form("Invalid pedestal (%f) or noise (%f) "
- "for FMD%d%c[%02d,%03d]", ped, noise, det, rng, sec, str));
- return USHRT_MAX;
- }
-
- AliFMDDebug(5, ("Subtracting pedestal %f from signal %d", ped, adc));
- // if (digit->Count3() > 0) adc = digit->Count3();
- // else if (digit->Count2() > 0) adc = digit->Count2();
- // else adc = digit->Count1();
- Int_t counts = adc + Int_t(zs ? fac * noise : - ped);
- counts = TMath::Max(Int_t(counts), 0);
- if (counts < noise * fNoiseFactor) counts = 0;
- if (counts > 0) AliFMDDebug(15, ("Got a hit strip"));
+ UShort_t counts = SubtractPedestal(det, rng, sec, str, adc,
+ fNoiseFactor, fZS[det-1],
+ fZSFactor[det-1]);
if (fDiagStep1) fDiagStep1->Fill(adc, counts);
- return UShort_t(counts);
+ return counts;
}
//____________________________________________________________________
Char_t rng,
UShort_t sec,
UShort_t str,
- Float_t eta,
UShort_t count) const
{
// Converts number of ADC counts to energy deposited.
// rng Ring ID
// sec Sector ID
// rng Strip ID
- // eta Psuedo-rapidity
// counts Number of ADC counts over pedestal
// Return
// The energy deposited in a single strip, or -1 in case of problems
Double_t edep = ((count * param->GetEdepMip())
/ (gain * param->GetDACPerMIP()));
+ return edep;
+}
+
+//____________________________________________________________________
+Float_t
+AliFMDReconstructor::Adc2Energy(UShort_t det,
+ Char_t rng,
+ UShort_t sec,
+ UShort_t str,
+ Float_t eta,
+ UShort_t count) const
+{
+ // Converts number of ADC counts to energy deposited.
+ // Note, that this member function can be overloaded by derived
+ // classes to do strip-specific look-ups in databases or the like,
+ // to find the proper gain for a strip.
+ //
+ // In the first simple version, we calculate the energy deposited as
+ //
+ // EnergyDeposited = cos(theta) * gain * count
+ //
+ // where
+ //
+ // Pre_amp_MIP_Range
+ // gain = ----------------- * Energy_deposited_per_MIP
+ // ADC_channel_size
+ //
+ // is constant and the same for all strips.
+ //
+ // For the production we use the conversion measured in the NBI lab.
+ // The total conversion is then:
+ //
+ // gain = ADC / DAC
+ //
+ // EdepMip * count
+ // => energy = ----------------
+ // gain * DACPerADC
+ //
+ // Parameters:
+ // det Detector ID
+ // rng Ring ID
+ // sec Sector ID
+ // rng Strip ID
+ // eta Psuedo-rapidity
+ // counts Number of ADC counts over pedestal
+ // Return
+ // The energy deposited in a single strip, or -1 in case of problems
+ //
+ Double_t edep = Adc2Energy(det, rng, sec, str, count);
+
if (fDiagStep2) fDiagStep2->Fill(count, edep);
if (fAngleCorrect) {
Double_t theta = 2 * TMath::ATan(TMath::Exp(-eta));
*/
virtual void FillESD(AliRawReader*, TTree* clusterTree,
AliESDEvent* esd) const;
+
+ /**
+ * Create SDigits from raw data
+ *
+ * @param reader The raw reader
+ * @param sdigits Array to fill with AliFMDSDigit objects.
+ */
+ virtual void Digitize(AliRawReader* reader,
+ TClonesArray* sdigits) const;
+
/**
* Not used
*/
UShort_t sec,
UShort_t str,
Short_t adc) const;
+ /**
+ * Process the signal from a single strip.
+ *
+ * @param sdigits Array to fill
+ * @param det Detector number
+ * @param rng Ring identifier
+ * @param sec Sector number
+ * @param str Strip number
+ * @param sam Sample number
+ * @param adc Number of ADC counts for this strip
+ */
+ virtual void DigitizeSignal(TClonesArray* sdigits,
+ UShort_t det,
+ Char_t rng,
+ UShort_t sec,
+ UShort_t str,
+ UShort_t sam,
+ Short_t adc) const;
+ /**
+ * Subtract the pedestal off the ADC counts.
+ *
+ * @param det Detector number
+ * @param rng Ring identifier
+ * @param sec Sector number
+ * @param str Strip number
+ * @param adc ADC counts
+ * @param noiseFactor If pedestal substracted pedestal is less then
+ * this times the noise, then consider this to be 0.
+ * @param zsEnabled Whether zero-suppression is on.
+ * @param zsNoiseFactor Noise factor used in on-line pedestal
+ * subtraction.
+ *
+ * @return The pedestal subtracted ADC counts (possibly 0), or @c
+ * USHRT_MAX in case of problems.
+ */
+ virtual UShort_t SubtractPedestal(UShort_t det,
+ Char_t rng,
+ UShort_t sec,
+ UShort_t str,
+ UShort_t adc,
+ Float_t noiseFactor,
+ Bool_t zsEnabled,
+ UShort_t zsNoiseFactor) const;
/**
* Substract pedestals from raw ADC in @a digit
*
*
* @return Energy deposited @f$ E_i@f$
*/
+ virtual Float_t Adc2Energy(UShort_t det,
+ Char_t rng,
+ UShort_t sec,
+ UShort_t str,
+ UShort_t count) const;
+ /**
+ * Converts number of ADC counts to energy deposited. This is
+ * done by
+ * @f[
+ * E_i = A_i g_i
+ * @f]
+ * where @f$ A_i@f$ is the pedestal subtracted ADC counts, and @f$
+ * g_i@f$ is the gain for the @f$ i^{\mbox{th}}@f$ strip.
+ *
+ * @param det Detector number
+ * @param rng Ring identifier
+ * @param sec Sector number
+ * @param str Strip number
+ * @param eta Psuedo-rapidity of digit.
+ * @param count Pedestal subtracted ADC counts
+ *
+ * @return Energy deposited @f$ E_i@f$
+ */
virtual Float_t Adc2Energy(UShort_t det,
Char_t rng,
UShort_t sec,
fCount2(-1),
fCount3(-1),
fCount4(-1),
- fLabels(0)
+ fNParticles(0),
+ fNPrimaries(0)
+ // , fLabels(0)
{
// cTOR
}
fCount3(count3),
fCount4(count4),
fNParticles(npart),
- fNPrimaries(nprim),
- fLabels(refs)
+ fNPrimaries(nprim)
+ // , fLabels(refs)
{
//
// Creates a real data digit object
// count1 ADC count (a 10-bit word)
// count2 ADC count (a 10-bit word) -1 if not used
// count3 ADC count (a 10-bit word) -1 if not used
+ for (Int_t i = 0; i < refs.fN; i++) AddTrack(refs.fArray[i]);
}
//____________________________________________________________________
{
// Print digit to standard out
AliFMDBaseDigit::Print();
- std::cout << "\t" << fEdep << " -> "
- << fCount1 << " (" << fCount2 << "," << fCount3 << ","
- << fCount4 << ") = " << Counts() << std::flush;
+ std::cout << "\t"
+ << std::setw(10) << fEdep << " -> "
+ << std::setw(4) << fCount1 << " ("
+ << std::setw(4) << fCount2 << ","
+ << std::setw(4) << fCount3 << ","
+ << std::setw(4) << fCount4 << ") = "
+ << std::setw(4) << Counts() << std::flush;
+
TString opt(option);
if (opt.Contains("p", TString::kIgnoreCase))
- std::cout << " [" << fNPrimaries << "/" << fNParticles << "]"
+ std::cout << " ["
+ << std::setw(2) << fNPrimaries << "/"
+ << std::setw(2) << fNParticles << "]"
<< std::flush;
if (opt.Contains("l", TString::kIgnoreCase)) {
std::cout << " ";
+ for (Int_t i = 0; i < GetNTrack(); i++)
+ std::cout << (i == 0 ? "" : ",") << std::setw(5) << fTracks[i];
+#if 0
for (Int_t i = 0; i < fLabels.fN; i++)
std::cout << (i == 0 ? "" : ",") << fLabels.fArray[i];
+#endif
}
std::cout << std::endl;
}
class AliFMDSDigit : public AliFMDBaseDigit
{
public:
- /** CTOR */
+ /**
+ * CTOR
+ */
AliFMDSDigit();
- /** Constrctor
- @param detector Detector
- @param ring Ring
- @param sector Sector
- @param strip Strip
- @param edep Energy deposited
- @param count ADC (first sample)
- @param count2 ADC (second sample, or -1 if not used)
- @param count3 ADC (third sample, or -1 if not used) */
+ /**
+ * Constrctor
+ *
+ * @param detector Detector
+ * @param ring Ring
+ * @param sector Sector
+ * @param strip Strip
+ * @param edep Energy deposited
+ * @param count ADC (first sample)
+ * @param count2 ADC (second sample, or -1 if not used)
+ * @param count3 ADC (third sample, or -1 if not used)
+ */
AliFMDSDigit(UShort_t detector,
Char_t ring='\0',
UShort_t sector=0,
UShort_t npart=0,
UShort_t nprim=0,
const TArrayI& refs=TArrayI());
- /** DTOR */
+ /**
+ * DTOR
+ */
virtual ~AliFMDSDigit() {}
- /** @return ADC count (first sample) */
- UShort_t Count1() const { return fCount1; }
- /** @return ADC count (second sample, or -1 if not used) */
- Short_t Count2() const { return fCount2; }
- /** @return ADC count (third sample, or -1 if not used) */
- Short_t Count3() const { return fCount3; }
- /** @return ADC count (third sample, or -1 if not used) */
- Short_t Count4() const { return fCount4; }
- /** @return Canonical ADC counts */
- UShort_t Counts() const;
- /** @return Energy deposited */
- Float_t Edep() const { return fEdep; }
- /** @return Number of particles that hit this strip */
+ /**
+ *
+ * @return ADC count (first sample)
+ */
+ UShort_t Count1() const { return fCount1; }
+ /**
+ *
+ * @return ADC count (second sample, or -1 if not used)
+ */
+ Short_t Count2() const { return fCount2; }
+ /**
+ *
+ * @return ADC count (third sample, or -1 if not used)
+ */
+ Short_t Count3() const { return fCount3; }
+ /**
+ *
+ * @return ADC count (third sample, or -1 if not used)
+ */
+ Short_t Count4() const { return fCount4; }
+ /**
+ *
+ * @return Canonical ADC counts
+ */
+ UShort_t Counts() const;
+ /**
+ *
+ * @return Energy deposited
+ */
+ Float_t Edep() const { return fEdep; }
+ /**
+ *
+ * @return Number of particles that hit this strip
+ */
UShort_t NParticles() const { return fNParticles; }
- /** @return Number of primary particles that hit this strip */
+ /**
+ *
+ * @return Number of primary particles that hit this strip
+ */
UShort_t NPrimaries() const { return fNPrimaries; }
+#if 0
/** @return the track labels */
const TArrayI& TrackLabels() const { return fLabels; }
- /** Print info
- @param opt Not used */
+#endif
+ /**
+ * Print info
+ *
+ * @param opt Not used
+ */
void Print(Option_t* opt="") const;
+ /**
+ * Set the count value
+ *
+ * @param s Sample number
+ * @param c Counts
+ */
+ void SetCount(UShort_t s, Short_t c);
protected:
Float_t fEdep; // Energy deposited
UShort_t fCount1; // Digital signal
Short_t fCount4; // Digital signal (-1 if not used)
UShort_t fNParticles; // Total number of particles that hit this strip
UShort_t fNPrimaries; // Number of primary particles that his this strip
+#if 0
TArrayI fLabels; // MC-truth track labels
- ClassDef(AliFMDSDigit,4) // Summable FMD digit
+#endif
+ ClassDef(AliFMDSDigit,5) // Summable FMD digit
};
inline UShort_t
return fCount1;
}
+inline void
+AliFMDSDigit::SetCount(UShort_t i, Short_t c)
+{
+ switch (i) {
+ case 0: fCount1 = c; break;
+ case 1: fCount2 = c; break;
+ case 2: fCount3 = c; break;
+ case 3: fCount4 = c; break;
+ }
+}
#endif
//____________________________________________________________________
* provided "as is" without express or implied warranty. *
**************************************************************************/
-//Struct to encode a strip address into one integer
-//developed by Christian Holm Christensen (cholm@nbi.dk).
-// The functions are static
-// to ensure applicability from anywhere. This is needed to smoothly store
-//strip addresses in track references.
+// Struct to encode a strip address into one integer
+// developed by Christian Holm Christensen (cholm@nbi.dk).
+//
+// The functions are static to ensure applicability from
+// anywhere. This is needed to smoothly store strip addresses in track
+// references.
+//
// Added by Hans H. Dalsgaard (hans.dalsgaard@cern.ch)
-struct AliFMDStripIndex
+class AliFMDStripIndex
{
+public:
+ /**
+ * Constructor
+ *
+ */
+ AliFMDStripIndex() {}
+ /**
+ * Destructor
+ *
+ */
+ virtual ~AliFMDStripIndex() {}
+ /**
+ * Pack an identifier from detector coordinates
+ *
+ * @param det Detector
+ * @param rng Ring
+ * @param sec Sector
+ * @param str Strip
+ *
+ * @return Packed identifier
+ */
static UInt_t Pack(UShort_t det, Char_t rng, UShort_t sec, UShort_t str)
{
UInt_t irg = (rng == 'I' || rng == 'i' ? 0 : 1);
((det & 0x003) << 17));
return id;
}
+ /**
+ * Unpack an identifier to detector coordinates
+ *
+ * @param id Identifier to unpack
+ * @param det On return, the detector
+ * @param rng On return, the ring
+ * @param sec On return, the sector
+ * @param str On return, the strip
+ */
static void Unpack(UInt_t id,
- UShort_t& det, Char_t& rng, UShort_t& sec, UShort_t& str)
+ UShort_t& det, Char_t& rng, UShort_t& sec, UShort_t& str)
{
str = ((id >> 0) & 0x1FF);
sec = ((id >> 9) & 0x03F);
det = ((id >> 17) & 0x003);
}
- ClassDef(AliFMDStripIndex,0)
+ ClassDef(AliFMDStripIndex,1)
};
#endif
+//
+// Local Variables:
+// mode: C++
+// End:
+//
#include <TGeoPhysicalNode.h>
#include "AliFMDGeometry.h"
+//____________________________________________________________________
+Double_t
+AliFMDSurveyToAlignObjs::GetUnitFactor() const
+{
+ // Returns the conversion factor from the measured values to
+ // centimeters.
+ TString units(fSurveyObj->GetUnits());
+ if (units.CompareTo("mm", TString::kIgnoreCase) == 0) return .1;
+ else if (units.CompareTo("cm", TString::kIgnoreCase) == 0) return 1.;
+ else if (units.CompareTo("m", TString::kIgnoreCase) == 0) return 100.;
+ return 1;
+}
+
//____________________________________________________________________
Bool_t
-AliFMDSurveyToAlignObjs::GetFMD2Plane(Double_t* rot,
- Double_t* trans)
+AliFMDSurveyToAlignObjs::GetPoint(const char* name,
+ TVector3& point,
+ TVector3& error) const
{
+ // Get named point. On return, point will contain the point
+ // coordinates in centimeters, and error will contain the
+ // meassurement errors in centimeters too. If no point is found,
+ // returns false, otherwise true.
+ Double_t unit = GetUnitFactor();
+ TObject* obj = fSurveyPoints->FindObject(name);
+ if (!obj) return kFALSE;
+
+ AliSurveyPoint* p = static_cast<AliSurveyPoint*>(obj);
+ point.SetXYZ(unit * p->GetX(),
+ unit * p->GetY(),
+ unit * p->GetZ());
+ error.SetXYZ(unit * p->GetPrecisionX(),
+ unit * p->GetPrecisionY(),
+ unit * p->GetPrecisionZ());
+ return kTRUE;
+}
- // The possile survey points
- const char* names[] = { "FMD2_ITOP", "FMD2_OTOP",
- "FMD2_IBOTM", "FMD2_OBOTM",
- "FMD2_IBOT", "FMD2_OBOT",
- 0 };
- const char** name = names;
- Int_t i = 0;
- TGraph2DErrors g;
+//____________________________________________________________________
+Bool_t
+AliFMDSurveyToAlignObjs::CalculatePlane(const TVector3& a,
+ const TVector3& b,
+ const TVector3& c,
+ Double_t* trans,
+ Double_t* rot) const
+{
+ // Vector a->b, b->c, and normal to plane defined by these two
+ // vectors.
+ TVector3 ab(b-a), bc(c-b);
- Double_t unit = 1.;
- TString units(fSurveyObj->GetUnits());
- if (units.CompareTo("mm", TString::kIgnoreCase) == 0) unit = .1;
- else if (units.CompareTo("cm", TString::kIgnoreCase) == 0) unit = 1.;
- else if (units.CompareTo("m", TString::kIgnoreCase) == 0) unit = 100.;
+ // Normal vector to the plane of the fiducial marks obtained
+ // as cross product of the two vectors on the plane d0^d1
+ TVector3 nn(ab.Cross(bc));
+ if (nn.Mag() < 1e-8) {
+ Info("DoIt", "Normal vector is null vector");
+ return kFALSE;
+ }
+ // We express the plane in Hessian normal form.
+ //
+ // n x = -p,
+ //
+ // where n is the normalised normal vector given by
+ //
+ // n_x = a / l, n_y = b / l, n_z = c / l, p = d / l
+ //
+ // with l = sqrt(a^2+b^2+c^2) and a, b, c, and d are from the
+ // normal plane equation
+ //
+ // ax + by + cz + d = 0
+ //
+ // Normalize
+ TVector3 n(nn.Unit());
+ Double_t p = - (n * a);
+
+ // The center of the square with the fiducial marks as the
+ // corners. The mid-point of one diagonal - md. Used to get the
+ // center of the surveyd box.
+ TVector3 md(a + c);
+ md *= 1/2.;
+
+ // The center of the box.
+ TVector3 orig(md);
+ trans[0] = orig[0];
+ trans[1] = orig[1];
+ trans[2] = orig[2];
+
+ // Normalize the spanning vectors
+ TVector3 uab(ab.Unit());
+ TVector3 ubc(bc.Unit());
+
+ for (size_t i = 0; i < 3; i++) {
+ rot[i * 3 + 0] = uab[i];
+ rot[i * 3 + 1] = ubc[i];
+ rot[i * 3 + 2] = n[i];
+ }
+ return kTRUE;
+}
+
+
+//____________________________________________________________________
+Bool_t
+AliFMDSurveyToAlignObjs::GetFMD1Plane(Double_t* rot, Double_t* trans) const
+{
+
+ // The possile survey points
+ TVector3 icb, ict, ocb, oct, dummy;
+ Int_t missing = 0;
+ if (!GetPoint("V0L_ICB", icb, dummy)) missing++;
+ if (!GetPoint("V0L_ICT", icb, dummy)) missing++;
+ if (!GetPoint("V0L_OCB", icb, dummy)) missing++;
+ if (!GetPoint("V0L_OCT", icb, dummy)) missing++;
+
+ // Check that we have enough points
+ if (missing > 1) {
+ Error("GetFMD1Plane", "Only got %d survey points - no good for FMD1 plane",
+ 4-missing);
+ return kFALSE;
+ }
+
+ if (!CalculatePlane(icb, ict, oct, rot, trans)) return kFALSE;
+
+ return kTRUE;
+}
+
+//____________________________________________________________________
+Bool_t
+AliFMDSurveyToAlignObjs::DoFMD1()
+{
+ return kTRUE;
+}
+
+//____________________________________________________________________
+Bool_t
+AliFMDSurveyToAlignObjs::GetFMD2Plane(Double_t* rot, Double_t* trans) const
+{
+
+ // The possile survey points
+ const char* names[] = { "FMD2_ITOP", "FMD2_OTOP",
+ "FMD2_IBOTM", "FMD2_OBOTM",
+ "FMD2_IBOT", "FMD2_OBOT",
+ 0 };
+ const char** name = names;
+ Int_t i = 0;
+ TGraph2DErrors g;
// Loop and fill graph
while (*name) {
- TObject* obj = fSurveyPoints->FindObject(*name);
- name++;
- if (!obj) continue;
-
- AliSurveyPoint* p = static_cast<AliSurveyPoint*>(obj);
- Double_t x = unit * p->GetX();
- Double_t y = unit * p->GetY();
- Double_t z = unit * p->GetZ();
- Double_t ex = unit * p->GetPrecisionX();
- Double_t ey = unit * p->GetPrecisionY();
- Double_t ez = unit * p->GetPrecisionZ();
-
- g.SetPoint(i, x, y, z);
- g.SetPointError(i, ex, ey, ez);
+ TVector3 p, e;
+ if (!GetPoint(*name++, p, e)) continue;
+
+ g.SetPoint(i, p.X(), p.Y(), p.Z());
+ g.SetPointError(i, e.X(), e.Y(), e.Z());
i++;
}
TVector3 ua(a.Unit());
TVector3 ub(b.Unit());
Double_t angAb = ua.Angle(ub);
- PrintVector("ua: ", ua);
- PrintVector("ub: ", ub);
- std::cout << "Angle: " << angAb * 180 / TMath::Pi() << std::endl;
+ // PrintVector("ua: ", ua);
+ // PrintVector("ub: ", ub);
+ // std::cout << "Angle: " << angAb * 180 / TMath::Pi() << std::endl;
for (size_t i = 0; i < 3; i++) {
rot[i * 3 + 0] = ua[i];
void Run();
Bool_t CreateAlignObjs() { return kTRUE; }
protected:
+ Bool_t DoFMD1();
+ Bool_t GetFMD1Plane(Double_t* rot, Double_t* trans) const;
Bool_t DoFMD2();
- Bool_t GetFMD2Plane(Double_t* rot, Double_t* trans);
-
+ Bool_t GetFMD2Plane(Double_t* rot, Double_t* trans) const;
+
+ Double_t GetUnitFactor() const;
+ Bool_t GetPoint(const char* name, TVector3& p, TVector3& e) const;
+ Bool_t CalculatePlane(const TVector3& a,
+ const TVector3& b,
+ const TVector3& c,
+ Double_t* trans,
+ Double_t* rot) const;
static void PrintVector(const char* text, const Double_t* v);
static void PrintVector(const char* text, const TVector3& v);
static void PrintRotation(const char* text, const Double_t* rot);
// Functions
Float_t EtaToTheta(Float_t eta);
EG_t LookupEG(const Char_t* name);
-AliGenerator* GeneratorFactory(EG_t eg, Rad_t rad, TString& comment);
+AliGenerator* GeneratorFactory(EG_t eg, Rad_t rad);
AliGenHijing* HijingStandard();
void ProcessEnvironmentVars(EG_t& eg, Int_t& seed);
Int_t seed = 12345; //Set 0 to use the current time
MC_t mc = kGEANT3TGEO;
- //____________________________________________________________________
- // Comment line
- static TString comment;
-
//____________________________________________________________________
// Get settings from environment variables
ProcessEnvironmentVars(eg, seed);
//__________________________________________________________________
// Generator Configuration
- AliGenerator* gener = GeneratorFactory(eg, rad, comment);
+ AliGenerator* gener = GeneratorFactory(eg, rad);
gener->SetOrigin(0, 0, 0); // vertex position
gener->SetSigma(0, 0, 5.3); // Sigma in (X,Y,Z) (cm) on IP position
gener->SetCutVertexZ(1.); // Truncate at 1 sigma
gener->SetTrackingFlag(1);
gener->Init();
- //__________________________________________________________________
- //
- // Comments
- //
- switch (mag) {
- case k2kG: comment = comment.Append(" | L3 field 0.2 T"); break;
- case k4kG: comment = comment.Append(" | L3 field 0.4 T"); break;
- case k5kG: comment = comment.Append(" | L3 field 0.5 T"); break;
- }
-
- switch (rad) {
- case kGluonRadiation:
- comment = comment.Append(" | Gluon Radiation On"); break;
- default:
- comment = comment.Append(" | Gluon Radiation Off"); break;
- }
-
- switch(geo) {
- case kHoles: comment = comment.Append(" | Holes for PHOS/HMPID"); break;
- default: comment = comment.Append(" | No holes for PHOS/HMPID"); break;
- }
-
- std::cout << "\n\n Comment: " << comment << "\n" << std::endl;
-
//__________________________________________________________________
// Field (L3 0.4 T)
AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., mag);
//____________________________________________________________________
AliGenerator*
-GeneratorFactory(EG_t eg, Rad_t rad, TString& comment)
+GeneratorFactory(EG_t eg, Rad_t rad)
{
Int_t isw = 3;
if (rad == kNoGluonRadiation) isw = 0;
switch (eg) {
case test50:
{
- comment = comment.Append(":HIJINGparam test 50 particles");
AliGenHIJINGpara *gener = new AliGenHIJINGpara(50);
gener->SetMomentumRange(0, 999999.);
gener->SetPhiRange(0., 360.);
break;
case kParam_8000:
{
- comment = comment.Append(":HIJINGparam N=8000");
AliGenHIJINGpara *gener = new AliGenHIJINGpara(86030);
gener->SetMomentumRange(0, 999999.);
gener->SetPhiRange(0., 360.);
break;
case kParam_4000:
{
- comment = comment.Append("HIJINGparam N=4000");
AliGenHIJINGpara *gener = new AliGenHIJINGpara(43015);
gener->SetMomentumRange(0, 999999.);
gener->SetPhiRange(0., 360.);
break;
case kParam_2000:
{
- comment = comment.Append("HIJINGparam N=2000");
AliGenHIJINGpara *gener = new AliGenHIJINGpara(21507);
gener->SetMomentumRange(0, 999999.);
gener->SetPhiRange(0., 360.);
break;
case kParam_fmd:
{
- comment = comment.Append("HIJINGparam N=100");
AliGenHIJINGpara *gener = new AliGenHIJINGpara(500);
gener->SetMomentumRange(0, 999999.);
gener->SetPhiRange(0., 360.);
//
case kHijing_cent1:
{
- comment = comment.Append("HIJING cent1");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
break;
case kHijing_cent2:
{
- comment = comment.Append("HIJING cent2");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 2.);
//
case kHijing_per1:
{
- comment = comment.Append("HIJING per1");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(5., 8.6);
break;
case kHijing_per2:
{
- comment = comment.Append("HIJING per2");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(8.6, 11.2);
break;
case kHijing_per3:
{
- comment = comment.Append("HIJING per3");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(11.2, 13.2);
break;
case kHijing_per4:
{
- comment = comment.Append("HIJING per4");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(13.2, 15.);
break;
case kHijing_per5:
{
- comment = comment.Append("HIJING per5");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(15., 100.);
//
case kHijing_jj25:
{
- comment = comment.Append("HIJING Jet 25 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_jj50:
{
- comment = comment.Append("HIJING Jet 50 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_jj75:
{
- comment = comment.Append("HIJING Jet 75 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_jj100:
{
- comment = comment.Append("HIJING Jet 100 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_jj200:
{
- comment = comment.Append("HIJING Jet 200 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
//
case kHijing_gj25:
{
- comment = comment.Append("HIJING Gamma 25 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_gj50:
{
- comment = comment.Append("HIJING Gamma 50 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_gj75:
{
- comment = comment.Append("HIJING Gamma 75 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_gj100:
{
- comment = comment.Append("HIJING Gamma 100 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
case kHijing_gj200:
{
- comment = comment.Append("HIJING Gamma 200 GeV");
AliGenHijing *gener = HijingStandard();
// impact parameter range
gener->SetImpactParameterRange(0., 5.);
break;
case kHijing_pA:
{
- comment = comment.Append("HIJING pA");
AliGenCocktail *gener = new AliGenCocktail();
break;
case kPythia6:
{
- comment = comment.Append(":Pythia p-p @ 14 TeV");
AliGenPythia *gener = new AliGenPythia(-1);
gener->SetMomentumRange(0,999999);
gener->SetThetaRange(0., 180.);
break;
case kPythia6Jets20_24:
{
- comment = comment.Append(":Pythia jets 20-24 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets24_29:
{
- comment = comment.Append(":Pythia jets 24-29 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets29_35:
{
- comment = comment.Append(":Pythia jets 29-35 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets35_42:
{
- comment = comment.Append(":Pythia jets 35-42 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets42_50:
{
- comment = comment.Append(":Pythia jets 42-50 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets50_60:
{
- comment = comment.Append(":Pythia jets 50-60 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets60_72:
{
- comment = comment.Append(":Pythia jets 60-72 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets72_86:
{
- comment = comment.Append(":Pythia jets 72-86 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets86_104:
{
- comment = comment.Append(":Pythia jets 86-104 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets104_125:
{
- comment = comment.Append(":Pythia jets 105-125 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets125_150:
{
- comment = comment.Append(":Pythia jets 125-150 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kPythia6Jets150_180:
{
- comment = comment.Append(":Pythia jets 150-180 GeV @ 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);// Centre of mass energy
gener->SetProcess(kPyJets);// Process type
break;
case kD0PbPb5500:
{
- comment = comment.Append(" D0 in Pb-Pb at 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(10);
gener->SetProcess(kPyD0PbPbMNR);
gener->SetStrucFunc(kCTEQ4L);
break;
case kCharmSemiElPbPb5500:
{
- comment = comment.Append(" Charm in Pb-Pb at 5.5 TeV");
AliGenPythia * gener = new AliGenPythia(10);
gener->SetProcess(kPyCharmPbPbMNR);
gener->SetStrucFunc(kCTEQ4L);
break;
case kBeautySemiElPbPb5500:
{
- comment = comment.Append(" Beauty in Pb-Pb at 5.5 TeV");
AliGenPythia *gener = new AliGenPythia(10);
gener->SetProcess(kPyBeautyPbPbMNR);
gener->SetStrucFunc(kCTEQ4L);
break;
case kCocktailTRD:
{
- comment = comment.Append(" Cocktail for TRD at 5.5 TeV");
AliGenCocktail *gener = new AliGenCocktail();
AliGenParam *jpsi = new AliGenParam(10,
break;
case kPyJJ:
{
- comment = comment.Append(" Jet-jet at 5.5 TeV");
AliGenPythia *gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);
gener->SetProcess(kPyJets);
break;
case kPyGJ:
{
- comment = comment.Append(" Gamma-jet at 5.5 TeV");
AliGenPythia *gener = new AliGenPythia(-1);
gener->SetEnergyCMS(5500.);
gener->SetProcess(kPyDirectGamma);
break;
case kMuonCocktailCent1:
{
- comment = comment.Append(" Muon Cocktail Cent1");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailPer1:
{
- comment = comment.Append(" Muon Cocktail Per1");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailPer4:
{
- comment = comment.Append(" Muon Cocktail Per4");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailCent1HighPt:
{
- comment = comment.Append(" Muon Cocktail HighPt Cent1");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailPer1HighPt :
{
- comment = comment.Append(" Muon Cocktail HighPt Per1");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailPer4HighPt:
{
- comment = comment.Append(" Muon Cocktail HighPt Per4");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailCent1Single:
{
- comment = comment.Append(" Muon Cocktail Single Cent1");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailPer1Single :
{
- comment = comment.Append(" Muon Cocktail Single Per1");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kMuonCocktailPer4Single:
{
- comment = comment.Append(" Muon Cocktail Single Per4");
AliGenMUONCocktail * gener = new AliGenMUONCocktail();
gener->SetPtRange(1.0,100.); // Transverse momentum range
gener->SetPhiRange(0.,360.); // Azimuthal angle range
break;
case kFMD1Flat:
{
- comment = comment.Append(" Flat in FMD1 range");
AliGenBox* gener = new AliGenBox(2000);
gener->SetPart(kPiPlus);
gener->SetMomentumRange(3,4);
break;
case kFMD2Flat:
{
- comment = comment.Append(" Flat in FMD2 range");
AliGenBox* gener = new AliGenBox(10);
gener->SetPart(kPiPlus);
gener->SetMomentumRange(3,4);
break;
case kFMD3Flat:
{
- comment = comment.Append(" Flat in FMD3 range");
AliGenBox* gener = new AliGenBox(2000);
gener->SetPart(kPiPlus);
gener->SetMomentumRange(3,4);
break;
case kFMDFlat:
{
- comment = comment.Append(" Flat in FMD range");
AliGenCocktail* gener = new AliGenCocktail();
gener->SetMomentumRange(3,4);
gener->SetPhiRange(0, 360);
Reconstruct()
{
// Debug the FMD.
- // AliLog::SetModuleDebugLevel("FMD", 1);
+ AliLog::SetModuleDebugLevel("FMD", 1);
// To reconstruct raw data from FDR-I, please enable below lines:
// AliFMDParameters::Instance()->UseRcuTrailer(false);
// AliFMDParameters::Instance()->UseCompleteHeader(false);
- TFile* magF = TFile::Open("mag.root", "READ");
- AliMagF* mag = static_cast<AliMagF*>(magF->Get("mag"));
- if (!mag) return;
- AliTracker::SetFieldMap(mag, true);
+ // TFile* magF = TFile::Open("mag.root", "READ");
+ // AliMagF* mag = static_cast<AliMagF*>(magF->Get("mag"));
+ // if (!mag) return;
+ // AliTracker::SetFieldMap(mag, true);
AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT");
AliCDBManager::Instance()->SetRun(0);
rec.SetRunReconstruction("FMD");
rec.SetRunTracking("");
rec.SetFillESD("FMD");
- rec.SetRunQA("");
- rec.SetInput(".");
+ rec.SetRunQA(":");
+ rec.SetInput("raw.root");
// rec.SetRecoParam("TOF", new AliTOFRecoParam());
rec.Run();
sim.SetConfigFile("$(ALICE_ROOT)/FMD/Config.C");
sim.SetMakeSDigits("FMD");
sim.SetMakeDigits("FMD");
- sim.SetWriteRawData("FMD");
+ sim.SetWriteRawData("FMD", "raw.root");
// sim.SetMakeDigitsFromHits("FMD");
- sim.SetRunQA(":");
+ sim.SetRunQA("FMD:ALL");
TStopwatch w;
w.Start();
sim.Run(n);
--- /dev/null
+> Title:
+
+Measurement of targets on FMD2
+
+> Date:
+
+12.08.2008
+
+> Subdetector:
+
+FMD1
+
+> Report URL:
+
+https://edms.cern.ch/document/943928
+
+> Version:
+
+1
+
+> General Observations:
+
+The V0A detector has been installed and first measured on 18.07.2008.
+Afterwards it has been moved close to final position and measured on
+25.07.2008. The final alignment of the V0A detector has been done in
+three iterations (called sessions) on 12.08.2008.
+
+> Coordinate System:
+
+ALICEPH
+
+> Units:
+
+mm
+
+> Nr Columns:
+
+7
+
+> Column Names:
+
+Point Name,XPH,YPH,ZPH,Point Type,Target Used,Precision(mm)
+
+
+> Data:
+V0L_ICB 127.73 -220.41 3246.58 M Y 1
+V0L_ICT 127.37 219.48 3247.24 M Y 1
+V0L_OCB -126.99 -220.40 3247.27 M Y 1
+V0L_OCT -126.99 219.51 3247.49 M Y 1
--- /dev/null
+> Title:
+
+Measurement of targets on FMD2
+
+> Date:
+
+05.06.2008
+
+> Subdetector:
+
+FMD2
+
+> Report URL:
+
+https://edms.cern.ch/document/976326
+
+> Version:
+
+1
+
+> General Observations:
+
+During the installation period of detectors in the Space Frame,
+several targets on the FMD2 have been measured regularly in order to
+detect relative movements of the TPC.
+
+
+> Coordinate System:
+
+ALICEPH
+
+> Units:
+
+mm
+
+> Nr Columns:
+
+7
+
+> Column Names:
+
+Point Name,XPH,YPH,ZPH,Point Type,Target Used,Precision(mm)
+
+
+> Data:
+
+FMD2_IBOT 192.7 -257.4 892.4 M Y 1
+FMD2_OBOT -194.5 -257.2 894.3 M Y 1
+FMD2_OBOTM -271.7 -122.0 892.1 M Y 1
+FMD2_ITOP 246.9 100.4 892.5 M Y 1
+FMD2_OTOP -232.4 113.6 891.7 M Y 1
--- /dev/null
+//____________________________________________________________________
+//
+// $Id: DrawBothDigits.C 20907 2007-09-25 08:44:03Z cholm $
+//
+// Script that contains a class to draw eloss from hits, versus ADC
+// counts from digits, using the AliFMDInputHits class in the util library.
+//
+// It draws the energy loss versus the p/(mq^2). It can be overlayed
+// with the Bethe-Bloc curve to show how the simulation behaves
+// relative to the expected.
+//
+// Use the script `Compile.C' to compile this class using ACLic.
+//
+#include <TH2F.h>
+#include <AliFMDDigit.h>
+#include <AliFMDSDigit.h>
+#include <AliFMDInput.h>
+#include <AliFMDEdepMap.h>
+#include <iostream>
+#include <TStyle.h>
+#include <TArrayF.h>
+#include <AliLog.h>
+
+/** @class DrawBothDigits
+ @brief Draw hit energy loss versus digit ADC
+ @code
+ Root> .L Compile.C
+ Root> Compile("DrawBothDigits.C")
+ Root> DrawBothDigits c
+ Root> c.Run();
+ @endcode
+ @ingroup FMD_script
+ */
+class DrawBothDigits : public AliFMDInput
+{
+private:
+ TH2F* fTrackNos; // Histogram
+ AliFMDEdepMap fCache;
+public:
+ //__________________________________________________________________
+ DrawBothDigits(Int_t max=300)
+ : AliFMDInput("galice.root")
+ {
+ AddLoad(kDigits);
+ AddLoad(kSDigits);
+ fTrackNos = new TH2F("trackNos", "Track numbers",
+ max+1, -1.5, max-.5, max+1, -1.5, max-.5);
+ fTrackNos->SetXTitle("Digit track");
+ fTrackNos->SetYTitle("SDigit track");
+ }
+ //__________________________________________________________________
+ Bool_t Begin(Int_t evno)
+ {
+ fCache.Reset();
+ return AliFMDInput::Begin(evno);
+ }
+ //__________________________________________________________________
+ Bool_t ProcessSDigit(AliFMDSDigit* sdigit)
+ {
+ if (!sdigit) return kTRUE;
+ AliFMDEdepHitPair& entry = fCache(sdigit->Detector(),
+ sdigit->Ring(),
+ sdigit->Sector(),
+ sdigit->Strip());
+ entry.fLabels.Set(sdigit->GetNTrack());
+ Info("ProcessSDigit", "Got %d SDigit tracks", sdigit->GetNTrack());
+ for (size_t i = 0; i < sdigit->GetNTrack(); i++)
+ entry.fLabels.fArray[i] = sdigit->GetTrack(i);
+ return kTRUE;
+ }
+ //__________________________________________________________________
+ Bool_t ProcessDigit(AliFMDDigit* digit)
+ {
+ if (!digit) return kTRUE;
+ AliFMDEdepHitPair& entry = fCache(digit->Detector(),
+ digit->Ring(),
+ digit->Sector(),
+ digit->Strip());
+ TArrayI& stracks = entry.fLabels;
+ Info("ProcessDigit", "Got %d SDigit tracks, and %d Digit tracks",
+ stracks.fN, digit->GetNTrack());
+ for (Int_t i = 0; (i < stracks.fN || i < digit->GetNTrack()); i++) {
+ Int_t strack = (i < stracks.fN ? stracks.fArray[i] : -1);
+ Int_t dtrack = (i < digit->GetNTrack() ? digit->GetTrack(i) : -1);
+ fTrackNos->Fill(strack, dtrack);
+ }
+ return kTRUE;
+ }
+ //__________________________________________________________________
+ Bool_t Finish()
+ {
+ gStyle->SetPalette(1);
+ gStyle->SetOptTitle(0);
+ gStyle->SetCanvasColor(0);
+ gStyle->SetCanvasBorderSize(0);
+ gStyle->SetPadColor(0);
+ gStyle->SetPadBorderSize(0);
+ fTrackNos->SetStats(kFALSE);
+ // fTrackNos->Scale(1. / fAdc->GetEntries());
+ fTrackNos->Draw("colz");
+ return kTRUE;
+ }
+
+ ClassDef(DrawBothDigits,0);
+};
+
+//____________________________________________________________________
+//
+// EOF
+//
public:
//__________________________________________________________________
DrawDigits(Int_t m=1100, Double_t amin=-0.5, Double_t amax=1023.5)
+ : AliFMDInput("galice.root")
{
AddLoad(kDigits);
fAdc = new TH1D("adc", "ADC", m, amin, amax);
{
if (!digit) return kTRUE;
fAdc->Fill(digit->Counts());
+ digit->Print("l");
return kTRUE;
}
//__________________________________________________________________
return AliFMDInput::Begin(ev);
}
//__________________________________________________________________
- Bool_t ProcessESD(UShort_t /* det */,
- Char_t /* ring */,
- UShort_t /* sec */,
- UShort_t /* strip */,
+ Bool_t ProcessESD(UShort_t det,
+ Char_t rng,
+ UShort_t sec,
+ UShort_t str,
Float_t /* eta */,
Float_t mult)
{
// Cache the energy loss
- if (mult/fCorr > 0.01) fMult->Fill(mult/fCorr);
+ if (mult > 0)
+ Info("ProcessESD", "FMD%d%c[%2d,%3d]=%f", det, rng, sec, str, mult);
+ if (mult/fCorr > 0.001) fMult->Fill(mult/fCorr);
return kTRUE;
}
//__________________________________________________________________
TF1* FitPeak(Int_t n, TH1D* hist, Double_t min, Double_t& max)
{
+ if (TMath::Abs(max-min) < .25) return 0;
+ std::cout << "Fit peack in range " << min << " to " << max << std::endl;
TF1* l = new TF1(Form("l%02d", n), "landau", min, max);
hist->GetXaxis()->SetRangeUser(0, 4);
hist->Fit(l, "0Q", "", min, max);
//__________________________________________________________________
Bool_t Finish()
{
+ Info("Finish", "Will draw results");
gStyle->SetPalette(1);
gStyle->SetOptTitle(0);
gStyle->SetOptFit(1111111);
gStyle->SetPadColor(0);
gStyle->SetPadBorderSize(0);
+ if (fMult->GetEntries() <= 0) return kFALSE;
+
TCanvas* c = new TCanvas("c", "C");
c->cd();
- c->SetLogy();
+ // c->SetLogy();
fMult->GetXaxis()->SetRangeUser(0,4);
fMult->Scale(1. / fMult->GetEntries());
fMult->SetStats(kFALSE);
fMult->SetFillColor(2);
fMult->SetFillStyle(3001);
fMult->Draw("hist e");
-
+
+ return kTRUE;
+
Double_t mean, rms;
MaxInRange(fMult, 0.2, mean, rms);
Double_t x1 = mean-rms/2; // .75; // .8; // .65 / fCorr;
MaxInRange(fMult, x2, mean, rms);
Double_t x3 = mean + rms;
TF1* l2 = FitPeak(2, fMult, x2, x3);
- Double_t diff = l2->GetParameter(1)-l1->GetParameter(1);
- TF1* f = new TF1("user", "landau(0)+landau(3)", x1, x3);
-
+ TF1* f = 0;
+ Double_t diff = 0;
+ if (l2) {
+ diff = l2->GetParameter(1)-l1->GetParameter(1);
+ f = new TF1("user", "landau(0)+landau(3)", x1, x3);
+ f->SetParNames("A_{1}", "Mpv_{1}", "#sigma_{1}",
+ "A_{2}", "Mpv_{2}", "#sigma_{2}",
+ "A_{3}", "Mpv_{3}", "#sigma_{3}");
+ f->SetParameters(l1->GetParameter(0),
+ l1->GetParameter(1),
+ l1->GetParameter(2),
+ l2 ? l2->GetParameter(0) : 0,
+ l2 ? l2->GetParameter(1) : 0,
+ l2 ? l2->GetParameter(2) : 0,
+ l2->GetParameter(0)/10,
+ l2->GetParameter(1) + diff,
+ l2->GetParameter(2));
+ }
+ else {
+ x3 = x2;
+ f = new TF1("usr", "landau", x1, x3);
+ }
+
+ std::cout << "Range: " << x1 << "-" << x3 << std::endl;
+
fMult->GetXaxis()->SetRangeUser(0, 4);
- f->SetParNames("A_{1}", "Mpv_{1}", "#sigma_{1}",
- "A_{2}", "Mpv_{2}", "#sigma_{2}",
- "A_{3}", "Mpv_{3}", "#sigma_{3}");
- f->SetParameters(l1->GetParameter(0),
- l1->GetParameter(1),
- l1->GetParameter(2),
- l2->GetParameter(0),
- l2->GetParameter(1),
- l2->GetParameter(2),
- l2->GetParameter(0)/10,
- l2->GetParameter(1) + diff,
- l2->GetParameter(2));
- fMult->Fit(f, "0Q", "", x1, x3);
- fMult->Fit(f, "ME0", "E1", x1, x3);
+ fMult->Fit(f, "0QR", "", x1, x3);
+ fMult->Fit(f, "ME0R", "E1", x1, x3);
fMult->DrawClone("same hist");
l1->SetLineColor(3);
l1->SetLineWidth(2);
l1->SetRange(0, 4);
l1->Draw("same");
- l2->SetLineColor(4);
- l2->SetLineWidth(2);
- l2->SetRange(0, 4);
- l2->Draw("same");
+ if (l2) {
+ l2->SetLineColor(4);
+ l2->SetLineWidth(2);
+ l2->SetRange(0, 4);
+ l2->Draw("same");
+ }
f->SetLineWidth(2);
f->SetRange(0, 4);
f->Draw("same");
TLegend* l = new TLegend(0.6, 0.6, .89, .89);
l->AddEntry(l1, "1 particle Landau", "l");
- l->AddEntry(l2, "2 particle Landau", "l");
+ if (l2) l->AddEntry(l2, "2 particle Landau", "l");
l->AddEntry(f, "1+2 particle Landau", "l");
l->SetFillColor(0);
l->Draw("same");
-#if 1
+#if 0
c = new TCanvas("c2", "Landaus");
- c->SetLogy();
+ // c->SetLogy();
fMult->DrawClone("axis");
f->Draw("same");
Double_t* p1 = f->GetParameters();
std::cout << "No track" << std::endl;
return kFALSE;
}
- if (!p->IsPrimary()) return kTRUE;
+ // if (!p->IsPrimary()) return kTRUE;
if (hit->IsStop()) return kTRUE;
Float_t x = hit->P();
y /= q * q;
fElossVsPMQ->Fill(x, y);
fEloss->Fill(y);
-
+ fNHits++;
return kTRUE;
}
//__________________________________________________________________
c = new TCanvas("eloss", "Energy loss per unit material");
// c->SetLogx();
- c->SetLogy();
- fEloss->Scale(1. / fEloss->GetEntries());
- fEloss->GetXaxis()->SetRangeUser(1, 10);
- fEloss->Fit("landau", "", "", 1, 10);
- TF1* land = fEloss->GetFunction("landau");
- land->SetLineWidth(2);
- Double_t max = fEloss->GetMaximum();
- TGraph* resp = GetResp();
- Double_t* x = resp->GetX();
- Double_t* y = resp->GetY();
- TGraph* g = new TGraph(resp->GetN());
- TGraph* co = GetCorr();
- std::cout << "Correction factor: " << co->Eval(fBetaGammaMip) << std::endl;
- Double_t xs = fRho; // * 1.19; // /
- for (Int_t i = 0; i < g->GetN(); i++)
- g->SetPoint(i, x[i] * xs, y[i] * max);
- g->Draw("C same");
-
+ if (fEloss->GetEntries() != 0) {
+ c->SetLogy();
+ fEloss->Scale(1. / fEloss->GetEntries());
+ fEloss->GetXaxis()->SetRangeUser(1, 10);
+ fEloss->Fit("landau", "", "", 1, 10);
+ TF1* land = fEloss->GetFunction("landau");
+ land->SetLineWidth(2);
+ Double_t max = fEloss->GetMaximum();
+ TGraph* resp = GetResp();
+ Double_t* x = resp->GetX();
+ Double_t* y = resp->GetY();
+ TGraph* g = new TGraph(resp->GetN());
+ TGraph* co = GetCorr();
+ std::cout << "Correction factor: " << co->Eval(fBetaGammaMip)
+ << std::endl;
+ Double_t xs = fRho; // * 1.19; // /
+ for (Int_t i = 0; i < g->GetN(); i++)
+ g->SetPoint(i, x[i] * xs, y[i] * max);
+ g->Draw("C same");
+ }
+
l = new TLegend(.6, .6, .89, .89);
l->AddEntry(fEloss, fEloss->GetTitle(), "lf");
l->AddEntry(land, "Landau fit", "l");
//__________________________________________________________________
Bool_t ProcessHit(AliFMDHit* hit, TParticle*)
{
+ // Info("ProcessHit", "Processing hit %p", hit);
// Cache the energy loss
if (!hit) return kFALSE;
UShort_t det = hit->Detector();
}
fMap(det, rng, sec, str).fEdep += hit->Edep();
fMap(det, rng, sec, str).fN++;
+ // hit->Print();
return kTRUE;
}
//__________________________________________________________________
Bool_t ProcessDigit(AliFMDDigit* digit)
{
+ // Info("ProcessDigit", "Processing digit %p", digit);
if (!digit) return kFALSE;
UShort_t det = digit->Detector();
Char_t rng = digit->Ring();
return kFALSE;
}
fElossVsAdc->Fill(fMap(det, rng, sec, str).fEdep, digit->Counts());
+ if (fMap(det, rng, sec, str).fEdep > 0)
+ Info("ProcessDigit", "FMD%d%c[%2d,%3d] Edep=%8.5f -> ADC=0x%03x",
+ det, rng, sec, str,
+ fMap(det, rng, sec, str).fEdep, digit->Counts());
return kTRUE;
}
//__________________________________________________________________
--- /dev/null
+//____________________________________________________________________
+//
+// $Id: DrawHitsSDigits.C 20907 2007-09-25 08:44:03Z cholm $
+//
+// Script that contains a class to draw eloss from hits, versus ADC
+// counts from sdigits, using the AliFMDInputHits class in the util library.
+//
+// It draws the energy loss versus the p/(mq^2). It can be overlayed
+// with the Bethe-Bloc curve to show how the simulation behaves
+// relative to the expected.
+//
+// Use the script `Compile.C' to compile this class using ACLic.
+//
+#include <TH2D.h>
+#include <AliFMDHit.h>
+#include <AliFMDSDigit.h>
+#include <AliFMDInput.h>
+#include <AliFMDEdepMap.h>
+#include <iostream>
+#include <TStyle.h>
+#include <TArrayF.h>
+#include <AliLog.h>
+#include <TMath.h>
+
+/** @class DrawHitsSDigits
+ @brief Draw hit energy loss versus sdigit ADC
+ @code
+ Root> .L Compile.C
+ Root> Compile("DrawHitsSDigits.C")
+ Root> DrawHitsSDigits c
+ Root> c.Run();
+ @endcode
+ @ingroup FMD_script
+ */
+class DrawHitsSDigits : public AliFMDInput
+{
+private:
+ TH2D* fElossVsAdc; // Histogram
+ AliFMDEdepMap fMap;
+public:
+ //__________________________________________________________________
+ DrawHitsSDigits(Int_t n=900, Double_t emin=1e-3, Double_t emax=10,
+ Int_t m=1100, Double_t amin=-0.5, Double_t amax=1099.5)
+ {
+ AddLoad(kSDigits);
+ AddLoad(kHits);
+ TArrayF eloss(MakeLogScale(n, emin, emax));
+ TArrayF adcs(m+1);
+ adcs[0] = amin;
+ for (Int_t i = 1; i < m+1; i++) adcs[i] = adcs[i-1] + (amax-amin)/m;
+ fElossVsAdc = new TH2D("bad", "#Delta E vs. ADC",
+ eloss.fN-1, eloss.fArray, adcs.fN-1, adcs.fArray);
+ fElossVsAdc->SetXTitle("#Delta E/#Delta x [MeV/cm]");
+ fElossVsAdc->SetYTitle("ADC value");
+ }
+ //__________________________________________________________________
+ /** Begining of event
+ @param ev Event number
+ @return @c false on error */
+ Bool_t Begin(Int_t ev)
+ {
+ fMap.Reset();
+ return AliFMDInput::Begin(ev);
+ }
+ //__________________________________________________________________
+ Bool_t ProcessHit(AliFMDHit* hit, TParticle*)
+ {
+ Info("ProcessHit", "Processing hit %p", hit);
+ // Cache the energy loss
+ if (!hit) return kFALSE;
+ UShort_t det = hit->Detector();
+ Char_t rng = hit->Ring();
+ UShort_t sec = hit->Sector();
+ UShort_t str = hit->Strip();
+ if (str > 511) {
+ AliWarning(Form("Bad strip number %d in hit", str));
+ return kTRUE;
+ }
+ fMap(det, rng, sec, str).fEdep += hit->Edep();
+ fMap(det, rng, sec, str).fN++;
+ hit->Print();
+ return kTRUE;
+ }
+ //__________________________________________________________________
+ Bool_t ProcessSDigit(AliFMDSDigit* sdigit)
+ {
+ Info("ProcessSDigit", "Processing sdigit %p", sdigit);
+ if (!sdigit) return kFALSE;
+ UShort_t det = sdigit->Detector();
+ Char_t rng = sdigit->Ring();
+ UShort_t sec = sdigit->Sector();
+ UShort_t str = sdigit->Strip();
+ if (str > 511) {
+ AliWarning(Form("Bad strip number %d in sdigit", str));
+ return kFALSE;
+ }
+ fElossVsAdc->Fill(fMap(det, rng, sec, str).fEdep, sdigit->Counts());
+ sdigit->Print();
+ return kTRUE;
+ }
+ //__________________________________________________________________
+ Bool_t Finish()
+ {
+ gStyle->SetPalette(1);
+ gStyle->SetOptTitle(0);
+ gStyle->SetCanvasColor(0);
+ gStyle->SetCanvasBorderSize(0);
+ gStyle->SetPadColor(0);
+ gStyle->SetPadBorderSize(0);
+ fElossVsAdc->SetStats(kFALSE);
+ fElossVsAdc->Draw("COLZ");
+ return kTRUE;
+ }
+
+ ClassDef(DrawHitsSDigits,0);
+};
+
+//____________________________________________________________________
+//
+// EOF
+//
//__________________________________________________________________
DrawSDigits(Int_t m=1100, Double_t amin=-0.5, Double_t amax=1023.5)
+ : AliFMDInput("galice.root")
{
AddLoad(kSDigits);
AddLoad(kGeometry);
{
if (!digit) return kTRUE;
fAdc->Fill(digit->Counts());
+ digit->Print("lp");
if (digit->NParticles() == 0) return kTRUE;
Double_t ratio = digit->NPrimaries() / digit->NParticles();
if (phi < 0) phi += 360;
fPrimRatio[primIdx]->Fill(eta, phi, ratio);
-
return kTRUE;
}
//__________________________________________________________________
--- /dev/null
+void
+TestRaw2SDigits()
+{
+ AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT");
+ AliCDBManager::Instance()->SetRun(0);
+ AliLog::SetModuleDebugLevel("FMD", 2);
+ AliSimulation sim;
+ sim.SetConfigFile("$ALICE_ROOT/FMD/Config.C");
+ sim.ConvertRaw2SDigits("raw.root", 0);
+}