/*
$Log$
+Revision 1.12 1999/11/01 20:41:51 fca
+Added protections against using the wrong version of FRAME
+
Revision 1.11 1999/09/29 09:24:34 fca
Introduction of the Copyright and cvs Log
// //
// Transition Radiation Detector //
// This class contains the basic functions for the Transition Radiation //
-// detector. Functions specific to one particular geometry are //
-// contained in the derived classes //
+// Detector, as well as the geometry. //
+// Functions specific to one particular geometry are contained in the //
+// derived classes. //
// //
//Begin_Html
/*
// //
///////////////////////////////////////////////////////////////////////////////
+#include <stdlib.h>
+
#include <TMath.h>
#include <TNode.h>
#include <TPGON.h>
#include "AliTRD.h"
#include "AliRun.h"
-
#include "AliConst.h"
ClassImp(AliTRD)
fGasMix = 0;
fHits = 0;
fDigits = 0;
+ fHole = 0;
+
+ fClusters = 0;
+ fNclusters = 0;
// The chamber dimensions
for (Int_t iplan = 0; iplan < kNplan; iplan++) {
- fClengthI[iplan] = 0.;
- fClengthM[iplan] = 0.;
- fClengthO[iplan] = 0.;
+ fClengthI[iplan] = 0.;
+ fClengthM1[iplan] = 0.;
+ fClengthM2[iplan] = 0.;
+ fClengthO1[iplan] = 0.;
+ fClengthO2[iplan] = 0.;
+ fClengthO3[iplan] = 0.;
+ fCwidth[iplan] = 0.;
+ }
+
+ for (Int_t iplan = 0; iplan < kNplan; iplan++) {
+ for (Int_t icham = 0; icham < kNcham; icham++) {
+ for (Int_t isect = 0; isect < kNsect; isect++) {
+ fRowMax[iplan][icham][isect] = 0;
+ }
+ }
+ fColMax[iplan] = 0;
}
+ fTimeMax = 0;
+
+ fRowPadSize = 0;
+ fColPadSize = 0;
+ fTimeBinSize = 0;
}
// Standard constructor for the TRD
//
+ // Check that FRAME is there otherwise we have no place where to
+ // put TRD
+ AliModule* FRAME=gAlice->GetModule("FRAME");
+ if (!FRAME) {
+ Error("Ctor","TRD needs FRAME to be present\n");
+ exit(1);
+ }
+
+ // Define the TRD geometry according to the FRAME geometry
+ if (FRAME->IsVersion() == 0)
+ // With hole
+ fHole = 1;
+ else
+ // Without hole
+ fHole = 0;
+
// Allocate the hit array
- fHits = new TClonesArray("AliTRDhit", 405);
+ fHits = new TClonesArray("AliTRDhit" , 405);
// Allocate the digits array
- fDigits = new TClonesArray("AliTRDdigit",10000);
+ fDigits = new TClonesArray("AliTRDdigit" ,10000);
+
+ // Allocate the cluster array
+ fClusters = new TClonesArray("AliTRDcluster", 400);
+ fNclusters = 0;
fIshunt = 0;
fGasMix = 0;
// The chamber dimensions
for (Int_t iplan = 0; iplan < kNplan; iplan++) {
- fClengthI[iplan] = 0.;
- fClengthM[iplan] = 0.;
- fClengthO[iplan] = 0.;
- fCwidth[iplan] = 0.;
+ fClengthI[iplan] = 0.;
+ fClengthM1[iplan] = 0.;
+ fClengthM2[iplan] = 0.;
+ fClengthO1[iplan] = 0.;
+ fClengthO2[iplan] = 0.;
+ fClengthO3[iplan] = 0.;
+ fCwidth[iplan] = 0.;
}
+ for (Int_t iplan = 0; iplan < kNplan; iplan++) {
+ for (Int_t icham = 0; icham < kNcham; icham++) {
+ for (Int_t isect = 0; isect < kNsect; isect++) {
+ fRowMax[iplan][icham][isect] = 0;
+ }
+ }
+ fColMax[iplan] = 0;
+ }
+ fTimeMax = 0;
+
+ fRowPadSize = 0;
+ fColPadSize = 0;
+ fTimeBinSize = 0;
+
SetMarkerColor(kWhite);
}
delete fHits;
delete fDigits;
+ delete fClusters;
+
+}
+
+//_____________________________________________________________________________
+void AliTRD::AddCluster(Int_t *tracks, Int_t *clusters, Float_t *position)
+{
+ //
+ // Add a cluster for the TRD
+ //
+
+ TClonesArray &lclusters = *fClusters;
+ new(lclusters[fNclusters++]) AliTRDcluster(tracks,clusters,position);
}
// Author: Christoph Blume (C.Blume@gsi.de) 20/07/99
//
// The volumes:
- // TRD (Air) --- The TRD mother volume for one sector.
+ // TRD1-3 (Air) --- The TRD mother volumes for one sector.
// To be placed into the spaceframe.
//
// UAFI(/M/O) (Al) --- The aluminum frame of the inner(/middle/outer) chambers (readout)
Float_t par_dum[3];
Float_t par_trd[npar_trd];
Float_t par_cha[npar_cha];
- Int_t iplan;
Float_t xpos, ypos, zpos;
- Int_t *idtmed = fIdtmed->GetArray()-1299;
+ Int_t *idtmed = fIdtmed->GetArray() - 1299;
// The length of the inner chambers
- for (iplan = 0; iplan < kNplan; iplan++) fClengthI[iplan] = 110.0;
+ for (Int_t iplan = 0; iplan < kNplan; iplan++)
+ fClengthI[iplan] = 110.0;
+
// The length of the middle chambers
- fClengthM[0] = 123.5;
- fClengthM[1] = 131.0;
- fClengthM[2] = 138.5;
- fClengthM[3] = 146.0;
- fClengthM[4] = 153.0;
- fClengthM[5] = 160.5;
+ fClengthM1[0] = 123.5;
+ fClengthM1[1] = 131.0;
+ fClengthM1[2] = 138.5;
+ fClengthM1[3] = 146.0;
+ fClengthM1[4] = 153.0;
+ fClengthM1[5] = 160.5;
+
+ fClengthM2[0] = 123.5 - 7.0;
+ fClengthM2[1] = 131.0 - 7.0;
+ fClengthM2[2] = 138.5 - 7.0;
+ fClengthM2[3] = 146.0 - 7.0;
+ fClengthM2[4] = 153.0 - 7.0;
+ fClengthM2[5] = 160.4 - 7.0;
+
// The length of the outer chambers
- fClengthO[0] = 123.5;
- fClengthO[1] = 131.0;
- fClengthO[2] = 134.5;
- fClengthO[3] = 142.0;
- fClengthO[4] = 142.0;
- fClengthO[5] = 134.5;
+ fClengthO1[0] = 123.5;
+ fClengthO1[1] = 131.0;
+ fClengthO1[2] = 134.5;
+ fClengthO1[3] = 142.0;
+ fClengthO1[4] = 142.0;
+ fClengthO1[5] = 134.5;
+
+ fClengthO2[0] = 123.5;
+ fClengthO2[1] = 131.0;
+ fClengthO2[2] = 134.5;
+ fClengthO2[3] = 142.0;
+ fClengthO2[4] = 142.0;
+ fClengthO2[5] = 134.5;
+
+ fClengthO3[0] = 86.5;
+ fClengthO3[1] = 101.5;
+ fClengthO3[2] = 112.5;
+ fClengthO3[3] = 127.5;
+ fClengthO3[4] = 134.5;
+ fClengthO3[5] = 134.5;
// The width of the chambers
- fCwidth[0] = 99.6;
- fCwidth[1] = 104.1;
- fCwidth[2] = 108.5;
- fCwidth[3] = 112.9;
- fCwidth[4] = 117.4;
- fCwidth[5] = 121.8;
-
- // The TRD mother volume for one sector (Air) (dimensions identical to BTR1-3)
+ fCwidth[0] = 99.6;
+ fCwidth[1] = 104.1;
+ fCwidth[2] = 108.5;
+ fCwidth[3] = 112.9;
+ fCwidth[4] = 117.4;
+ fCwidth[5] = 121.8;
+
+ // The TRD mother volume for one sector (Air) (dimensions identical to BTR1)
par_trd[0] = kSwidth1/2.;
par_trd[1] = kSwidth2/2.;
- par_trd[2] = kSlength/2.;
+ par_trd[2] = kSlenTR1/2.;
par_trd[3] = kSheight/2.;
- gMC->Gsvolu("TRD ","TRD1",idtmed[1302-1],par_trd,npar_trd);
+ gMC->Gsvolu("TRD1","TRD1",idtmed[1302-1],par_trd,npar_trd);
+
+ // The TRD mother volume for one sector (Air) (dimensions identical to BTR2 + BTR3).
+ // Only used for the geometry with holes.
+ if (fHole) {
+
+ par_trd[0] = kSwidth1/2.;
+ par_trd[1] = kSwidth2/2.;
+ par_trd[2] = kSlenTR2/2.;
+ par_trd[3] = kSheight/2.;
+ gMC->Gsvolu("TRD2","TRD1",idtmed[1302-1],par_trd,npar_trd);
+
+ par_trd[0] = kSwidth1/2.;
+ par_trd[1] = kSwidth2/2.;
+ par_trd[2] = kSlenTR3/2.;
+ par_trd[3] = kSheight/2.;
+ gMC->Gsvolu("TRD3","TRD1",idtmed[1302-1],par_trd,npar_trd);
+
+ }
// The aluminum frames - readout + electronics (Al)
// The inner chambers
gMC->Gspos("UL11",1,"UAIO",xpos,ypos,zpos,0,"ONLY");
// Position the chambers in the TRD mother volume
- for (iplan = 1; iplan <= kNplan; iplan++) {
+ for (Int_t iplan = 1; iplan <= kNplan; iplan++) {
// The inner chambers ---------------------------------------------------------------
// the aluminum frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2.;
par_cha[0] = fCwidth[iplan-1]/2.;
par_cha[1] = fClengthI[iplan-1]/2.;
par_cha[2] = kCaframe/2.;
xpos = 0.;
ypos = 0.;
zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UAFI",iplan ,"TRD ",xpos,ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UAFI",iplan ,"TRD1",xpos,ypos,zpos,0,"MANY",par_cha,npar_cha);
// the inner part of the aluminum frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2. - kCathick;
par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
par_cha[1] = fClengthI[iplan-1]/2. - kCathick;
par_cha[2] = kCaframe/2.;
xpos = 0.;
ypos = 0.;
zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UAII",iplan ,"TRD ",xpos,ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UAII",iplan ,"TRD1",xpos,ypos,zpos,0,"ONLY",par_cha,npar_cha);
// the carbon frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2.;
par_cha[0] = fCwidth[iplan-1]/2.;
par_cha[1] = fClengthI[iplan-1]/2.;
par_cha[2] = kCcframe/2.;
xpos = 0.;
ypos = 0.;
zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UCFI",iplan ,"TRD ",xpos,ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UCFI",iplan ,"TRD1",xpos,ypos,zpos,0,"MANY",par_cha,npar_cha);
// the inner part of the carbon frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2. - kCcthick;
par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
par_cha[1] = fClengthI[iplan-1]/2. - kCcthick;
par_cha[2] = kCcframe/2.;
xpos = 0.;
ypos = 0.;
zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UCII",iplan ,"TRD ",xpos,ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UCII",iplan ,"TRD1",xpos,ypos,zpos,0,"ONLY",par_cha,npar_cha);
// The middle chambers --------------------------------------------------------------
// the aluminum frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2.;
par_cha[0] = fCwidth[iplan-1]/2.;
- par_cha[1] = fClengthM[iplan-1]/2.;
+ par_cha[1] = fClengthM1[iplan-1]/2.;
par_cha[2] = kCaframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1]/2.;
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1]/2.;
zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UAFM",iplan ,"TRD ",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
- gMC->Gsposp("UAFM",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UAFM",iplan ,"TRD1",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UAFM",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
// the inner part of the aluminum frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2. - kCathick;
- par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
- par_cha[1] = fClengthM[iplan-1]/2. - kCathick;
+ par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
+ par_cha[1] = fClengthM1[iplan-1]/2. - kCathick;
par_cha[2] = kCaframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1]/2.;
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1]/2.;
zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UAIM",iplan ,"TRD ",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
- gMC->Gsposp("UAIM",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UAIM",iplan ,"TRD1",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UAIM",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
// the carbon frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2.;
par_cha[0] = fCwidth[iplan-1]/2.;
- par_cha[1] = fClengthM[iplan-1]/2.;
+ par_cha[1] = fClengthM1[iplan-1]/2.;
par_cha[2] = kCcframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1]/2.;
- zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UCFM",iplan, "TRD ",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
- gMC->Gsposp("UCFM",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1]/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCFM",iplan ,"TRD1",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UCFM",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
// the inner part of the carbon frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2. - kCcthick;
- par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
- par_cha[1] = fClengthM[iplan-1]/2. - kCcthick;
+ par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
+ par_cha[1] = fClengthM1[iplan-1]/2. - kCcthick;
par_cha[2] = kCcframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1]/2.;
- zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UCIM",iplan ,"TRD ",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
- gMC->Gsposp("UCIM",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1]/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCIM",iplan ,"TRD1",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UCIM",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ // Only for the geometry with holes
+ if (fHole) {
+
+ // the aluminum frame
+ par_cha[0] = fCwidth[iplan-1]/2.;
+ par_cha[1] = fClengthM2[iplan-1]/2.;
+ par_cha[2] = kCaframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UAFM",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+
+ // the inner part of the aluminum frame
+ par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
+ par_cha[1] = fClengthM2[iplan-1]/2. - kCathick;
+ par_cha[2] = kCaframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UAIM",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ // the carbon frame
+ par_cha[0] = fCwidth[iplan-1]/2.;
+ par_cha[1] = fClengthM2[iplan-1]/2.;
+ par_cha[2] = kCcframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCFM",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+
+ // the inner part of the carbon frame
+ par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
+ par_cha[1] = fClengthM2[iplan-1]/2. - kCcthick;
+ par_cha[2] = kCcframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCIM",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ }
// The outer chambers ---------------------------------------------------------------
// the aluminum frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2.;
par_cha[0] = fCwidth[iplan-1]/2.;
- par_cha[1] = fClengthO[iplan-1]/2.;
+ par_cha[1] = fClengthO1[iplan-1]/2.;
par_cha[2] = kCaframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1] + fClengthO[iplan-1]/2.;
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1] + fClengthO1[iplan-1]/2.;
zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UAFO",iplan ,"TRD ",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
- gMC->Gsposp("UAFO",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UAFO",iplan ,"TRD1",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UAFO",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
// the inner part of the aluminum frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2. - kCathick;
- par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
- par_cha[1] = fClengthO[iplan-1]/2. - kCathick;
+ par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
+ par_cha[1] = fClengthO1[iplan-1]/2. - kCathick;
par_cha[2] = kCaframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1] + fClengthO[iplan-1]/2.;
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1] + fClengthO1[iplan-1]/2.;
zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UAIO",iplan ,"TRD ",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
- gMC->Gsposp("UAIO",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UAIO",iplan ,"TRD1",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UAIO",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
// the carbon frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2.;
par_cha[0] = fCwidth[iplan-1]/2.;
- par_cha[1] = fClengthO[iplan-1]/2.;
+ par_cha[1] = fClengthO1[iplan-1]/2.;
par_cha[2] = kCcframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1] + fClengthO[iplan-1]/2.;
- zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UCFO",iplan, "TRD ",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
- gMC->Gsposp("UCFO",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1] + fClengthO1[iplan-1]/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCFO",iplan, "TRD1",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+ gMC->Gsposp("UCFO",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"MANY",par_cha,npar_cha);
// the inner part of the carbon frame
- //par_cha[0] = kSwidth1/2. + (iplan-1) * kCwidcha/2. - kCcthick;
- par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
- par_cha[1] = fClengthO[iplan-1]/2. - kCcthick;
+ par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
+ par_cha[1] = fClengthO1[iplan-1]/2. - kCcthick;
par_cha[2] = kCcframe/2.;
xpos = 0.;
- ypos = fClengthI[iplan-1]/2. + fClengthM[iplan-1] + fClengthO[iplan-1]/2.;
- zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
- gMC->Gsposp("UCIO",iplan ,"TRD ",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
- gMC->Gsposp("UCIO",iplan+kNplan,"TRD ",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ ypos = fClengthI[iplan-1]/2. + fClengthM1[iplan-1] + fClengthO1[iplan-1]/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCIO",iplan ,"TRD1",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+ gMC->Gsposp("UCIO",iplan+ kNplan,"TRD1",xpos,-ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ // Only for the geometry with holes
+ if (fHole) {
+
+ // the aluminum frame
+ par_cha[0] = fCwidth[iplan-1]/2.;
+ par_cha[1] = fClengthO2[iplan-1]/2.;
+ par_cha[2] = kCaframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1] + fClengthO2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UAFO",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+
+ // the inner part of the aluminum frame
+ par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
+ par_cha[1] = fClengthO2[iplan-1]/2. - kCathick;
+ par_cha[2] = kCaframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1] + fClengthO2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UAIO",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ // the carbon frame
+ par_cha[0] = fCwidth[iplan-1]/2.;
+ par_cha[1] = fClengthO2[iplan-1]/2.;
+ par_cha[2] = kCcframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1] + fClengthO2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCFO",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+
+ // the inner part of the carbon frame
+ par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
+ par_cha[1] = fClengthO2[iplan-1]/2. - kCcthick;
+ par_cha[2] = kCcframe/2.;
+ xpos = 0.;
+ ypos = fClengthM2[iplan-1] + fClengthO2[iplan-1]/2. - kSlenTR2/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCIO",iplan+2*kNplan,"TRD2",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ // the aluminum frame
+ par_cha[0] = fCwidth[iplan-1]/2.;
+ par_cha[1] = fClengthO3[iplan-1]/2.;
+ par_cha[2] = kCaframe/2.;
+ xpos = 0.;
+ ypos = fClengthO3[iplan-1]/2. - kSlenTR3/2.;
+ zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UAFO",iplan+4*kNplan,"TRD3",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+
+ // the inner part of the aluminum frame
+ par_cha[0] = fCwidth[iplan-1]/2. - kCathick;
+ par_cha[1] = fClengthO3[iplan-1]/2. - kCathick;
+ par_cha[2] = kCaframe/2.;
+ xpos = 0.;
+ ypos = fClengthO3[iplan-1]/2. - kSlenTR3/2.;
+ zpos = kCheight - kCaframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UAIO",iplan+4*kNplan,"TRD3",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ // the carbon frame
+ par_cha[0] = fCwidth[iplan-1]/2.;
+ par_cha[1] = fClengthO3[iplan-1]/2.;
+ par_cha[2] = kCcframe/2.;
+ xpos = 0.;
+ ypos = fClengthO3[iplan-1]/2. - kSlenTR3/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCFO",iplan+4*kNplan,"TRD3",xpos, ypos,zpos,0,"MANY",par_cha,npar_cha);
+
+ // the inner part of the carbon frame
+ par_cha[0] = fCwidth[iplan-1]/2. - kCcthick;
+ par_cha[1] = fClengthO3[iplan-1]/2. - kCcthick;
+ par_cha[2] = kCcframe/2.;
+ xpos = 0.;
+ ypos = fClengthO3[iplan-1]/2. - kSlenTR3/2.;
+ zpos = kCcframe/2. - kSheight/2. + (iplan-1) * (kCheight + kCspace);
+ gMC->Gsposp("UCIO",iplan+4*kNplan,"TRD3",xpos, ypos,zpos,0,"ONLY",par_cha,npar_cha);
+
+ }
}
+ if (fHole) {
+ xpos = 0.;
+ ypos = 0.;
+ zpos = 0.;
+ gMC->Gspos("TRD1",1,"BTR1",xpos,ypos,zpos,0,"ONLY");
+ gMC->Gspos("TRD2",1,"BTR2",xpos,ypos,zpos,0,"ONLY");
+ gMC->Gspos("TRD3",1,"BTR3",xpos,ypos,zpos,0,"ONLY");
+ }
+ else {
+ xpos = 0.;
+ ypos = 0.;
+ zpos = 0.;
+ gMC->Gspos("TRD1",1,"BTR1",xpos,ypos,zpos,0,"ONLY");
+ gMC->Gspos("TRD1",2,"BTR2",xpos,ypos,zpos,0,"ONLY");
+ gMC->Gspos("TRD1",3,"BTR3",xpos,ypos,zpos,0,"ONLY");
+ }
+
}
//_____________________________________________________________________________
gMC->Gsatt("ALIC","SEEN", 0);
// Set the volumes visible
- gMC->Gsatt("B032","SEEN", 0);
- gMC->Gsatt("B028","SEEN", 0);
- gMC->Gsatt("B029","SEEN", 0);
- gMC->Gsatt("B030","SEEN", 0);
- gMC->Gsatt("BTR1","SEEN", 0);
- gMC->Gsatt("BTR2","SEEN", 0);
- gMC->Gsatt("BTR3","SEEN", 0);
- gMC->Gsatt("TRD" ,"SEEN", 0);
+ if (fHole) {
+ gMC->Gsatt("B071","SEEN", 0);
+ gMC->Gsatt("B074","SEEN", 0);
+ gMC->Gsatt("B075","SEEN", 0);
+ gMC->Gsatt("B077","SEEN", 0);
+ gMC->Gsatt("BTR1","SEEN", 0);
+ gMC->Gsatt("BTR2","SEEN", 0);
+ gMC->Gsatt("BTR3","SEEN", 0);
+ gMC->Gsatt("TRD1","SEEN", 0);
+ gMC->Gsatt("TRD2","SEEN", 0);
+ gMC->Gsatt("TRD3","SEEN", 0);
+ }
+ else {
+ gMC->Gsatt("B071","SEEN", 0);
+ gMC->Gsatt("B074","SEEN", 0);
+ gMC->Gsatt("B075","SEEN", 0);
+ gMC->Gsatt("B077","SEEN", 0);
+ gMC->Gsatt("BTR1","SEEN", 0);
+ gMC->Gsatt("BTR2","SEEN", 0);
+ gMC->Gsatt("BTR3","SEEN", 0);
+ gMC->Gsatt("TRD1","SEEN", 0);
+ }
gMC->Gsatt("UCII","SEEN", 0);
gMC->Gsatt("UCIM","SEEN", 0);
gMC->Gsatt("UCIO","SEEN", 0);
// Initialise the TRD detector after the geometry has been created
//
+ Int_t i;
+
+ printf("\n");
+ for(i=0;i<35;i++) printf("*");
+ printf(" TRD_INIT ");
+ for(i=0;i<35;i++) printf("*");
+ printf("\n");
+
// Here the TRD initialisation code (if any!)
if (fGasMix == 1)
printf(" Gas Mixture: 90%% Xe + 10%% CO2\n");
else
printf(" Gas Mixture: 97%% Xe + 3%% Isobutane\n");
+ if (fHole)
+ printf(" Geometry with holes\n");
+ else
+ printf(" Full geometry\n");
+
+ // The default pad dimensions
+ if (!(fRowPadSize)) fRowPadSize = 4.5;
+ if (!(fColPadSize)) fColPadSize = 1.0;
+ if (!(fTimeBinSize)) fTimeBinSize = 0.1;
+
+ // The maximum number of pads
+ // and the position of pad 0,0,0
+ //
+ // chambers seen from the top:
+ // +----------------------------+
+ // | |
+ // | | ^
+ // | | rphi|
+ // | | |
+ // |0 | |
+ // +----------------------------+ +------>
+ // z
+ // chambers seen from the side: ^
+ // +----------------------------+ time|
+ // | | |
+ // |0 | |
+ // +----------------------------+ +------>
+ // z
+ //
+ for (Int_t iplan = 0; iplan < kNplan; iplan++) {
+
+ // The pad row (z-direction)
+ for (Int_t isect = 0; isect < kNsect; isect++) {
+ Float_t clengthI = fClengthI[iplan];
+ Float_t clengthM = fClengthM1[iplan];
+ Float_t clengthO = fClengthO1[iplan];
+ if (fHole) {
+ switch (isect) {
+ case 12:
+ case 13:
+ case 14:
+ case 15:
+ case 16:
+ clengthM = fClengthM2[iplan];
+ clengthO = fClengthO2[iplan];
+ break;
+ case 4:
+ case 5:
+ case 6:
+ clengthO = fClengthO3[iplan];
+ break;
+ };
+ }
+ fRowMax[iplan][0][isect] = 1 + TMath::Nint((clengthO - 2. * kCcthick)
+ / fRowPadSize - 0.5);
+ fRowMax[iplan][1][isect] = 1 + TMath::Nint((clengthM - 2. * kCcthick)
+ / fRowPadSize - 0.5);
+ fRowMax[iplan][2][isect] = 1 + TMath::Nint((clengthI - 2. * kCcthick)
+ / fRowPadSize - 0.5);
+ fRowMax[iplan][3][isect] = 1 + TMath::Nint((clengthM - 2. * kCcthick)
+ / fRowPadSize - 0.5);
+ fRowMax[iplan][4][isect] = 1 + TMath::Nint((clengthO - 2. * kCcthick)
+ / fRowPadSize - 0.5);
+ fRow0[iplan][0][isect] = -clengthI/2. - clengthM - clengthO + kCcthick;
+ fRow0[iplan][1][isect] = -clengthI/2. - clengthM + kCcthick;
+ fRow0[iplan][2][isect] = -clengthI/2. + kCcthick;
+ fRow0[iplan][3][isect] = clengthI/2. + kCcthick;
+ fRow0[iplan][4][isect] = clengthI/2. + clengthM + kCcthick;
+ }
+
+ // The pad column (rphi-direction)
+ fColMax[iplan] = 1 + TMath::Nint((fCwidth[iplan] - 2. * kCcthick)
+ / fColPadSize - 0.5);
+ fCol0[iplan] = -fCwidth[iplan]/2. + kCcthick;
+
+ }
+
+ // The time bucket
+ fTimeMax = 1 + TMath::Nint(kDrThick / fTimeBinSize - 0.5);
+ for (Int_t iplan = 0; iplan < kNplan; iplan++) {
+ fTime0[iplan] = kRmin + kCcframe/2. + kDrZpos - 0.5 * kDrThick
+ + iplan * (kCheight + kCspace);
+ }
+
+}
+
+//_____________________________________________________________________________
+void AliTRD::MakeBranch(Option_t* option)
+{
+ //
+ // Create Tree branches for the TRD digits and cluster.
+ //
+
+ Int_t buffersize = 4000;
+ Char_t branchname[15];
+
+ AliDetector::MakeBranch(option);
+
+ Char_t *D = strstr(option,"D");
+ sprintf(branchname,"%s",GetName());
+ if (fDigits && gAlice->TreeD() && D) {
+ gAlice->TreeD()->Branch(branchname,&fDigits, buffersize);
+ printf("* AliTRD::MakeBranch * Making Branch %s for digits in TreeD\n",branchname);
+ }
+
+ sprintf(branchname,"%scluster",GetName());
+ if (fClusters && gAlice->TreeD() && D) {
+ gAlice->TreeD()->Branch(branchname,&fClusters,buffersize);
+ printf("* AliTRD::MakeBranch * Making Branch %s for cluster in TreeD\n",branchname);
+ }
+
+}
+
+//_____________________________________________________________________________
+void AliTRD::SetTreeAddress()
+{
+ //
+ // Set the branch addresses for the trees.
+ //
+
+ Char_t branchname[15];
+
+ AliDetector::SetTreeAddress();
+
+ TBranch *branch;
+ TTree *treeD = gAlice->TreeD();
+
+ if (treeD) {
+ sprintf(branchname,"%scluster",GetName());
+ if (fClusters) {
+ branch = treeD->GetBranch(branchname);
+ if (branch) branch->SetAddress(&fClusters);
+ }
+ }
+
}
//_____________________________________________________________________________
}
+//______________________________________________________________________________
+void AliTRD::Streamer(TBuffer &R__b)
+{
+ // Stream an object of class AliTRD.
+
+ if (R__b.IsReading()) {
+ Version_t R__v = R__b.ReadVersion(); if (R__v) { }
+ AliDetector::Streamer(R__b);
+ R__b >> fGasMix;
+ R__b.ReadStaticArray(fClengthI);
+ R__b.ReadStaticArray(fClengthM1);
+ R__b.ReadStaticArray(fClengthM2);
+ R__b.ReadStaticArray(fClengthO1);
+ R__b.ReadStaticArray(fClengthO2);
+ R__b.ReadStaticArray(fClengthO3);
+ R__b.ReadStaticArray(fCwidth);
+ R__b.ReadStaticArray((int*)fRowMax);
+ R__b.ReadStaticArray(fColMax);
+ R__b >> fTimeMax;
+ R__b.ReadStaticArray((float*)fRow0);
+ R__b.ReadStaticArray(fCol0);
+ R__b.ReadStaticArray(fTime0);
+ R__b >> fRowPadSize;
+ R__b >> fColPadSize;
+ R__b >> fTimeBinSize;
+ R__b >> fHole;
+ // Stream the pointers but not the TClonesArray
+ R__b >> fClusters; // diff
+ //R__b >> fNclusters;
+ } else {
+ R__b.WriteVersion(AliTRD::IsA());
+ AliDetector::Streamer(R__b);
+ R__b << fGasMix;
+ R__b.WriteArray(fClengthI, 6);
+ R__b.WriteArray(fClengthM1, 6);
+ R__b.WriteArray(fClengthM2, 6);
+ R__b.WriteArray(fClengthO1, 6);
+ R__b.WriteArray(fClengthO2, 6);
+ R__b.WriteArray(fClengthO3, 6);
+ R__b.WriteArray(fCwidth, 6);
+ R__b.WriteArray((int*)fRowMax, 540);
+ R__b.WriteArray(fColMax, 6);
+ R__b << fTimeMax;
+ R__b.WriteArray((float*)fRow0, 540);
+ R__b.WriteArray(fCol0, 6);
+ R__b.WriteArray(fTime0, 6);
+ R__b << fRowPadSize;
+ R__b << fColPadSize;
+ R__b << fTimeBinSize;
+ R__b << fHole;
+ // Stream the pointers but not the TClonesArrays
+ R__b << fClusters; // diff
+ //R__b << fNclusters;
+ }
+
+}
+
ClassImp(AliTRDhit)
//_____________________________________________________________________________
-AliTRDhit::AliTRDhit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits):
- AliHit(shunt, track)
+AliTRDhit::AliTRDhit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits)
+ :AliHit(shunt, track)
{
//
// Create a TRD hit
fSignal = digits[6];
}
+
+ClassImp(AliTRDcluster)
+
+//_____________________________________________________________________________
+AliTRDcluster::AliTRDcluster(Int_t *tracks, Int_t *cluster, Float_t* position)
+ :TObject()
+{
+ //
+ // Create a TRD cluster
+ //
+
+ fSector = cluster[0];
+ fChamber = cluster[1];
+ fPlane = cluster[2];
+
+ fTimeSlice = cluster[3];
+ fEnergy = cluster[4];
+
+ fX = position[0];
+ fY = position[1];
+ fZ = position[2];
+
+ fTracks[0] = tracks[0];
+ fTracks[1] = tracks[1];
+ fTracks[2] = tracks[2];
+
+}
+
class AliTRD : public AliDetector {
protected:
- Int_t fGasMix; // Gas mixture. 0: Xe/Isobutane 1: Xe/CO2
+ Int_t fGasMix; // Gas mixture. 0: Xe/Isobutane 1: Xe/CO2
- Float_t fClengthI[kNplan]; // Length of the inner chambers
- Float_t fClengthM[kNplan]; // Length of the middle chambers
- Float_t fClengthO[kNplan]; // Length of the outer chambers
- Float_t fCwidth[kNplan]; // Width of the chambers
+ Float_t fClengthI[kNplan]; // Length of the inner chambers
+ Float_t fClengthM1[kNplan]; // Length of the middle chambers
+ Float_t fClengthM2[kNplan]; //
+ Float_t fClengthO1[kNplan]; // Length of the outer chambers
+ Float_t fClengthO2[kNplan]; //
+ Float_t fClengthO3[kNplan]; //
+ Float_t fCwidth[kNplan]; // Width of the chambers
+
+ Int_t fRowMax[kNplan][kNcham][kNsect]; // Number of pad-rows
+ Int_t fColMax[kNplan]; // Number of pad-columns
+ Int_t fTimeMax; // Number of time buckets
+
+ Float_t fRow0[kNplan][kNcham][kNsect]; // Row-position of pad 0
+ Float_t fCol0[kNplan]; // Column-position of pad 0
+ Float_t fTime0[kNplan]; // Time-position of pad 0
+
+ Float_t fRowPadSize; // Pad size in z-direction
+ Float_t fColPadSize; // Pad size in rphi-direction
+ Float_t fTimeBinSize; // Size of the time buckets
+
+ Int_t fHole; // Geometry without (0) / with (1) hole
+
+ TClonesArray *fClusters; // List of clusters
+ Int_t fNclusters; // Number of clusters
public:
AliTRD();
AliTRD(const char *name, const char *title);
- virtual ~AliTRD();
- virtual void AddHit(Int_t, Int_t*, Float_t*);
- virtual void AddDigit(Int_t*, Int_t*);
- virtual void BuildGeometry();
- virtual void CreateGeometry();
- virtual void CreateMaterials();
- virtual void DrawModule();
- Int_t DistancetoPrimitive(Int_t px, Int_t py);
- virtual void Init();
- virtual Int_t IsVersion() const = 0;
- virtual void StepManager() = 0;
- virtual void SetGasMix(Int_t imix = 0);
- virtual void SetHits(Int_t ) {};
- virtual void SetSensPlane(Int_t) {};
- virtual void SetSensChamber(Int_t) {};
- virtual void SetSensSector(Int_t ) {};
+ virtual ~AliTRD();
+ virtual void AddHit(Int_t, Int_t*, Float_t*);
+ virtual void AddDigit(Int_t*, Int_t*);
+ virtual void AddCluster(Int_t*, Int_t*, Float_t*);
+ virtual void BuildGeometry();
+ virtual void CreateGeometry();
+ virtual void CreateMaterials();
+ virtual void DrawModule();
+ Int_t DistancetoPrimitive(Int_t px, Int_t py);
+ TClonesArray *Cluster() { return fClusters; };
+ virtual void Init();
+ virtual Int_t IsVersion() const = 0;
+ virtual void MakeBranch(Option_t* option);
+ virtual void StepManager() = 0;
+ virtual void SetTreeAddress();
+
+ virtual void SetHits(Int_t ) {};
+ virtual void SetSensPlane(Int_t) {};
+ virtual void SetSensChamber(Int_t) {};
+ virtual void SetSensSector(Int_t ) {};
+
+ virtual void SetGasMix(Int_t imix = 0);
+
+ virtual void SetRowPadSize(Float_t size) { fRowPadSize = size; };
+ virtual void SetColPadSize(Float_t size) { fColPadSize = size; };
+ virtual void SetTimeBinSize(Float_t size) { fTimeBinSize = size; };
+
+ virtual Float_t GetRowPadSize() { return fRowPadSize; };
+ virtual Float_t GetColPadSize() { return fColPadSize; };
+ virtual Float_t GetTimeBinSize() { return fTimeBinSize; };
+
+ virtual Int_t GetRowMax(Int_t p, Int_t c, Int_t s) { return fRowMax[p-1][c-1][s-1]; };
+ virtual Int_t GetColMax(Int_t p) { return fColMax[p-1]; };
+ virtual Int_t GetTimeMax() { return fTimeMax; };
+
+ virtual Int_t Hole() { return fHole; };
ClassDef(AliTRD,1) // Transition Radiation Detector base class
Int_t fSector; // TRD sector number
Int_t fChamber; // TRD chamber number
Int_t fPlane; // TRD plane number
- Float_t fQ; // Charge created by a hit (geometry 2)
+ Float_t fQ; // Charge created by a hit (slow simulator only)
public:
AliTRDhit() {}
};
+//_____________________________________________________________________________
+class AliTRDcluster : public TObject {
+
+public:
+ Int_t fSector; // TRD sector number
+ Int_t fChamber; // TRD chamber number
+ Int_t fPlane; // TRD plane number
+
+ Int_t fTimeSlice; // Timeslice in chamber where cluster has been found
+ Int_t fEnergy; // Charge sum of this cluster
+
+ Float_t fX; // X coord in ALICE reference frame
+ Float_t fY; // Y coord in ALICE reference frame
+ Float_t fZ; // Z coord in ALICE reference frame
+
+ Int_t fTracks[3]; // Track information
+
+public:
+ AliTRDcluster() {};
+ AliTRDcluster(Int_t *tracks, Int_t *cluster, Float_t *pos);
+ virtual ~AliTRDcluster() {};
+
+ inline virtual Int_t *GetTracks() { return &fTracks[0]; }
+
+ ClassDef(AliTRDcluster,1) // Cluster for Transition Radiation Detector
+
+};
#endif
const Float_t kZmax1 = 378.35; // z-dimensions of the TRD
const Float_t kZmax2 = 302.0;
-const Float_t kSheight = 74.0; // Height of the TRD-volume in spaceframe
-const Float_t kSwidth1 = 103.674; // Lower width of the TRD-volume in spaceframe
-const Float_t kSwidth2 = 129.768; // Upper width of the TRD-volume in spaceframe
-const Float_t kSlength = 751.0; // Length of the TRD-volume in spaceframe
+const Float_t kSheight = 74.0; // Height of the TRD-volume in spaceframe (BTR1-3)
+const Float_t kSwidth1 = 99.613; // Lower width of the TRD-volume in spaceframe (BTR1-3)
+const Float_t kSwidth2 = 125.707; // Upper width of the TRD-volume in spaceframe (BTR1-3)
+const Float_t kSlenTR1 = 751.0; // Length of the TRD-volume in spaceframe (BTR1)
+const Float_t kSlenTR2 = 313.5; // Length of the TRD-volume in spaceframe (BTR2)
+const Float_t kSlenTR3 = 159.5; // Length of the TRD-volume in spaceframe (BTR3)
const Float_t kCheight = 11.0; // Height of the chambers
const Float_t kCspace = 1.6; // Vertical spacing of the chambers
virtual Int_t GetChamber() { return fChamber; };
virtual Int_t GetPlane() { return fPlane; };
- ClassDef(AliTRDmatrix,1)
+ ClassDef(AliTRDmatrix,1) // The TRD detector matrix for one readout chamber
};
virtual Float_t GetSignal() { return fSignal; };
virtual Int_t GetTrack(Int_t i) { return fTrack[i]; };
- ClassDef(AliTRDpixel,1)
+ ClassDef(AliTRDpixel,1) // Information for one detector pixel
};
/*
$Log$
+Revision 1.11 1999/11/01 20:41:51 fca
+Added protections against using the wrong version of FRAME
+
Revision 1.10 1999/09/29 09:24:35 fca
Introduction of the Copyright and cvs Log
///////////////////////////////////////////////////////////////////////////////
// //
-// Transition Radiation Detector version 0 -- coarse simulation //
+// Transition Radiation Detector version 0 -- fast simulator //
// //
//Begin_Html
/*
-<img src="picts/AliTRDv0Class.gif">
+<img src="picts/AliTRDfullClass.gif">
*/
//End_Html
// //
// //
///////////////////////////////////////////////////////////////////////////////
-#include <stdlib.h>
-
#include <TMath.h>
#include <TRandom.h>
#include <TVector.h>
fIdSens = 0;
fHitsOn = 0;
- fIdSpace1 = 0;
- fIdSpace2 = 0;
- fIdSpace3 = 0;
-
fIdChamber1 = 0;
fIdChamber2 = 0;
fIdChamber3 = 0;
+ fRphiSigma = 0;
+ fRphiDist = 0;
+
}
//_____________________________________________________________________________
// Create the GEANT geometry for the Transition Radiation Detector - Version 0
// This version covers the full azimuth.
//
- // Author: Christoph Blume (C.Blume@gsi.de) 20/07/99
- //
-
- Float_t xpos, ypos, zpos;
// Check that FRAME is there otherwise we have no place where to put the TRD
AliModule* FRAME = gAlice->GetModule("FRAME");
// Define the chambers
AliTRD::CreateGeometry();
- // Position the the TRD-sectors in all TRD-volumes in the spaceframe
- xpos = 0.;
- ypos = 0.;
- zpos = 0.;
- gMC->Gspos("TRD ",1,"BTR1",xpos,ypos,zpos,0,"ONLY");
- gMC->Gspos("TRD ",2,"BTR2",xpos,ypos,zpos,0,"ONLY");
- gMC->Gspos("TRD ",3,"BTR3",xpos,ypos,zpos,0,"ONLY");
-
}
//_____________________________________________________________________________
}
+//_____________________________________________________________________________
+void AliTRDv0::Hits2Clusters()
+{
+ // A simple cluster generator. It takes the hits from the
+ // fast simulator (one hit per plane) and transforms them
+ // into cluster, by applying position smearing and merging
+ // of nearby cluster. The searing is done uniformly in z-direction
+ // over the length of a readout pad. In rphi-direction a Gaussian
+ // smearing is applied with a sigma given by fRphiSigma.
+ // Clusters are considered as overlapping when they are closer in
+ // rphi-direction than the value defined in fRphiDist.
+ // Use the macro fastClusterCreate.C to create the cluster.
+
+ printf("AliTRDv0::Hits2Clusters -- Start creating cluster\n");
+
+ Int_t nBytes = 0;
+
+ AliTRDhit *TRDhit;
+
+ // Get the pointer to the hit tree
+ TTree *HitTree = gAlice->TreeH();
+ // Get the pointer to the reconstruction tree
+ TTree *ClusterTree = gAlice->TreeD();
+
+ TObjArray *Chamber = new TObjArray();
+
+ // Get the number of entries in the hit tree
+ // (Number of primary particles creating a hit somewhere)
+ Int_t nTrack = (Int_t) HitTree->GetEntries();
+
+ // Loop through all the chambers
+ for (Int_t icham = 0; icham < kNcham; icham++) {
+ for (Int_t iplan = 0; iplan < kNplan; iplan++) {
+ for (Int_t isect = 0; isect < kNsect; isect++) {
+
+ // Loop through all entries in the tree
+ for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) {
+
+ gAlice->ResetHits();
+ nBytes += HitTree->GetEvent(iTrack);
+
+ // Get the number of hits in the TRD created by this particle
+ Int_t nHit = fHits->GetEntriesFast();
+
+ // Loop through the TRD hits
+ for (Int_t iHit = 0; iHit < nHit; iHit++) {
+
+ if (!(TRDhit = (AliTRDhit *) fHits->UncheckedAt(iHit)))
+ continue;
+
+ Float_t x = TRDhit->fX;
+ Float_t y = TRDhit->fY;
+ Float_t z = TRDhit->fZ;
+ Int_t track = TRDhit->fTrack;
+ Int_t plane = TRDhit->fPlane;
+ Int_t sector = TRDhit->fSector;
+ Int_t chamber = TRDhit->fChamber;
+
+ if ((sector != isect+1) ||
+ (plane != iplan+1) ||
+ (chamber != icham+1))
+ continue;
+
+ // Rotate the sectors on top of each other
+ Float_t phi = 2.0 * kPI / (Float_t) kNsect
+ * ((Float_t) sector - 0.5);
+ Float_t xRot = -x * TMath::Cos(phi) + y * TMath::Sin(phi);
+ Float_t yRot = x * TMath::Sin(phi) + y * TMath::Cos(phi);
+ Float_t zRot = z;
+
+ // Add this cluster to the temporary cluster-array for this chamber
+ Int_t tracks[3];
+ tracks[0] = track;
+ Int_t clusters[5];
+ clusters[0] = sector;
+ clusters[1] = chamber;
+ clusters[2] = plane;
+ clusters[3] = 0;
+ clusters[4] = 0;
+ Float_t position[3];
+ position[0] = zRot;
+ position[1] = yRot;
+ position[2] = xRot;
+ AliTRDcluster *Cluster = new AliTRDcluster(tracks,clusters,position);
+ Chamber->Add(Cluster);
+
+ }
+
+ }
+
+ // Loop through the temporary cluster-array
+ for (Int_t iClus1 = 0; iClus1 < Chamber->GetEntries(); iClus1++) {
+
+ AliTRDcluster *Cluster1 = (AliTRDcluster *) Chamber->UncheckedAt(iClus1);
+ Float_t x1 = Cluster1->fX;
+ Float_t y1 = Cluster1->fY;
+ Float_t z1 = Cluster1->fZ;
+
+ if (!(z1)) continue; // Skip marked cluster
+
+ const Int_t nSave = 2;
+ Int_t idxSave[nSave];
+ Int_t iSave = 0;
+
+ Int_t tracks[3];
+ tracks[0] = Cluster1->fTracks[0];
+
+ // Check the other cluster to see, whether there are close ones
+ for (Int_t iClus2 = iClus1 + 1; iClus2 < Chamber->GetEntries(); iClus2++) {
+ AliTRDcluster *Cluster2 = (AliTRDcluster *) Chamber->UncheckedAt(iClus2);
+ Float_t x2 = Cluster2->fX;
+ Float_t y2 = Cluster2->fY;
+ if ((TMath::Abs(x1 - x2) < fRowPadSize) ||
+ (TMath::Abs(y1 - y2) < fRphiDist)) {
+ if (iSave == nSave) {
+ printf("AliTRDv0::Hits2Clusters -- Boundary error: iSave = %d, nSave = %d\n"
+ ,iSave,nSave);
+ }
+ else {
+ idxSave[iSave] = iClus2;
+ tracks[iSave+1] = Cluster2->fTracks[0];
+ }
+ iSave++;
+ }
+ }
+
+ // Merge close cluster
+ Float_t yMerge = y1;
+ Float_t xMerge = x1;
+ if (iSave) {
+ for (Int_t iMerge = 0; iMerge < iSave; iMerge++) {
+ AliTRDcluster *Cluster2 = (AliTRDcluster *) Chamber->UncheckedAt(idxSave[iMerge]);
+ xMerge += Cluster2->fX;
+ yMerge += Cluster2->fY;
+ Cluster2->fZ = 0; // Mark merged cluster
+ }
+ xMerge /= (iSave + 1);
+ yMerge /= (iSave + 1);
+ }
+
+ // The position smearing in z-direction (uniform over pad width)
+ Int_t row = (Int_t) ((xMerge - fRow0[iplan][icham][isect]) / fRowPadSize);
+ Float_t xSmear = (row + gRandom->Rndm()) * fRowPadSize
+ + fRow0[iplan][icham][isect];
+
+ // The position smearing in rphi-direction (Gaussian)
+ Float_t ySmear = 0;
+ do
+ ySmear = gRandom->Gaus(yMerge,fRphiSigma);
+ while ((ySmear < fCol0[iplan]) ||
+ (ySmear > fCol0[iplan] + fColMax[iplan] * fColPadSize));
+
+ // Time direction stays unchanged
+ Float_t zSmear = z1;
+
+ Int_t clusters[5];
+ clusters[0] = Cluster1->fSector;
+ clusters[1] = Cluster1->fChamber;
+ clusters[2] = Cluster1->fPlane;
+ clusters[3] = 0;
+ clusters[4] = 0;
+ Float_t position[3];
+ // Rotate the sectors back into their real position
+ Float_t phi = 2*kPI / kNsect * ((Float_t) Cluster1->fSector - 0.5);
+ position[0] = -zSmear * TMath::Cos(phi) + ySmear * TMath::Sin(phi);
+ position[1] = zSmear * TMath::Sin(phi) + ySmear * TMath::Cos(phi);
+ position[2] = xSmear;
+
+ // Add the smeared cluster to the output array
+ AddCluster(tracks,clusters,position);
+
+ }
+
+ // Clear the temporary cluster-array and delete the cluster
+ Chamber->Delete();
+
+ }
+ }
+ }
+
+ printf("AliTRDv0::Hits2Clusters -- Found %d cluster\n",fClusters->GetEntries());
+ printf("AliTRDv0::Hits2Clusters -- Fill the cluster tree\n");
+ ClusterTree->Fill();
+
+}
+
//_____________________________________________________________________________
void AliTRDv0::Init()
{
AliTRD::Init();
- printf("**************************************"
- " TRD "
- "**************************************\n");
- printf("\n Version 0 of TRD initialing, "
- "symmetric TRD\n\n");
-
- //
- // Check that FRAME is there otherwise we have no place where to
- // put TRD
- AliModule* FRAME=gAlice->GetModule("FRAME");
- if(!FRAME) {
- Error("Ctor","TRD needs FRAME to be present\n");
- exit(1);
- } else
- if(FRAME->IsVersion()!=1) {
- Error("Ctor","FRAME version 1 needed with this version of TRD\n");
- exit(1);
- }
-
- for (Int_t i = 0; i < 80; i++) printf("*");
- printf("\n");
-
// Identifier of the sensitive volume (amplification region)
fIdSens = gMC->VolId("UL06");
- // Identifier of the TRD-spaceframe volumina
- fIdSpace1 = gMC->VolId("B028");
- fIdSpace2 = gMC->VolId("B029");
- fIdSpace3 = gMC->VolId("B030");
-
// Identifier of the TRD-driftchambers
fIdChamber1 = gMC->VolId("UCIO");
fIdChamber2 = gMC->VolId("UCIM");
fIdChamber3 = gMC->VolId("UCII");
- printf("**************************************"
- " TRD "
- "**************************************\n");
+ // Parameter for Hits2Cluster
+
+ // Position resolution in rphi-direction
+ fRphiSigma = 0.02;
+ // Minimum distance of non-overlapping cluster
+ fRphiDist = 1.0;
+
+ printf(" Fast simulator\n");
+ for (Int_t i = 0; i < 80; i++) printf("*");
+ printf("\n");
+
}
//_____________________________________________________________________________
Int_t vol[3];
Int_t iIdSens, icSens;
- Int_t iIdSpace, icSpace;
Int_t iIdChamber, icChamber;
- Int_t secMap1[10] = { 3, 7, 8, 9, 10, 11, 2, 1, 18, 17 };
- Int_t secMap2[ 5] = { 16, 15, 14, 13, 12 };
- Int_t secMap3[ 3] = { 5, 6, 4 };
-
Float_t hits[4];
TLorentzVector p;
// No charge created
hits[3] = 0;
- iIdSpace = gMC->CurrentVolOffID(4,icSpace );
- iIdChamber = gMC->CurrentVolOffID(1,icChamber);
-
// The sector number
- if (iIdSpace == fIdSpace1)
- vol[0] = secMap1[icSpace-1];
- else if (iIdSpace == fIdSpace2)
- vol[0] = secMap2[icSpace-1];
- else if (iIdSpace == fIdSpace3)
- vol[0] = secMap3[icSpace-1];
+ Float_t phi = hits[1] != 0 ? TMath::Atan2(hits[0],hits[1]) : (hits[0] > 0 ? 180. : 0.);
+ vol[0] = ((Int_t) (phi / 20)) + 1;
// The chamber number
// 1: outer left
// 3: inner
// 4: middle right
// 5: outer right
+ iIdChamber = gMC->CurrentVolOffID(1,icChamber);
if (iIdChamber == fIdChamber1)
vol[1] = (hits[2] < 0 ? 1 : 5);
else if (iIdChamber == fIdChamber2)
public:
AliTRDv0() {}
AliTRDv0(const char *name, const char *title);
- virtual ~AliTRDv0() {}
- virtual void CreateGeometry();
- virtual void CreateMaterials();
- virtual Int_t IsVersion() const { return 0; };
- virtual void SetHits(Int_t ihit = 1) { fHitsOn = ihit; };
- virtual void StepManager();
- virtual void Init();
+ virtual ~AliTRDv0() {}
+ virtual void CreateGeometry();
+ virtual void CreateMaterials();
+ virtual Int_t IsVersion() const { return 0; };
+ virtual void Hits2Clusters();
+ virtual void SetHits(Int_t ihit = 1) { fHitsOn = ihit; };
+ virtual void StepManager();
+ virtual void Init();
+ virtual void SetRphiSigma(Float_t sigma) { fRphiSigma = sigma; };
+ virtual void SetRphiDist(Float_t dist) { fRphiDist = dist; };
+
+ virtual Float_t GetRphiSigma() { return fRphiSigma; };
+ virtual Float_t GetRphiDist() { return fRphiDist; };
+
protected:
Int_t fIdSens; // Sensitive volume identifier
- Int_t fIdSpace1; // Spaceframe volume identifier
- Int_t fIdSpace2; //
- Int_t fIdSpace3; //
-
Int_t fIdChamber1; // Driftchamber volume identifier
Int_t fIdChamber2; //
Int_t fIdChamber3; //
Int_t fHitsOn; // Used to switch hits on
- ClassDef(AliTRDv0,1) // Transition Radiation Detector version 0
+ Float_t fRphiSigma; // Gaussian position smearing in rphi-direction
+ Float_t fRphiDist; // Maximum distnace for non-overlapping cluster
+
+ ClassDef(AliTRDv0,1) // Transition Radiation Detector version 0 (fast simulator)
};
/*
$Log$
+Revision 1.11 1999/11/01 20:41:51 fca
+Added protections against using the wrong version of FRAME
+
Revision 1.10 1999/09/29 09:24:35 fca
Introduction of the Copyright and cvs Log
///////////////////////////////////////////////////////////////////////////////
// //
-// Transition Radiation Detector version 1 -- coarse simulation //
-// This version has two detector arms, leaving the space in front of the //
-// HMPID and PHOS empty //
+// Transition Radiation Detector version 2 -- slow simulator //
// //
//Begin_Html
/*
-<img src="picts/AliTRDv1Class.gif">
+<img src="picts/AliTRDfullClass.gif">
*/
//End_Html
// //
// //
///////////////////////////////////////////////////////////////////////////////
-#include <stdlib.h>
-
#include <TMath.h>
-#include <TRandom.h>
#include <TVector.h>
+#include <TRandom.h>
#include "AliTRDv1.h"
+#include "AliTRDmatrix.h"
#include "AliRun.h"
#include "AliMC.h"
#include "AliConst.h"
-
+
ClassImp(AliTRDv1)
//_____________________________________________________________________________
:AliTRD(name, title)
{
//
- // Standard constructor for the Transition Radiation Detector version 1
+ // Standard constructor for Transition Radiation Detector version 2
//
- fIdSens = 0;
- fHitsOn = 0;
+ fIdSens = 0;
+
+ fIdChamber1 = 0;
+ fIdChamber2 = 0;
+ fIdChamber3 = 0;
+
+ fSensSelect = 0;
+ fSensPlane = 0;
+ fSensChamber = 0;
+ fSensSector = 0;
- fIdSpace1 = 0;
- fIdSpace2 = 0;
- fIdSpace3 = 0;
+ fGasGain = 0;
+ fNoise = 0;
+ fChipGain = 0;
+ fADCoutRange = 0;
+ fADCinRange = 0;
+ fADCthreshold = 0;
+
+ fDiffusionT = 0;
+ fDiffusionL = 0;
+
+ fClusMaxThresh = 0;
+ fClusSigThresh = 0;
+ fClusMethod = 0;
+
+ fDeltaE = NULL;
+
+ SetBufferSize(128000);
+
+}
+
+//_____________________________________________________________________________
+AliTRDv1::~AliTRDv1()
+{
- fIdChamber1 = 0;
- fIdChamber2 = 0;
- fIdChamber3 = 0;
+ if (fDeltaE) delete fDeltaE;
}
void AliTRDv1::CreateGeometry()
{
//
- // Create the GEANT geometry for the Transition Radiation Detector - Version 1
- // This version covers only part of the azimuth.
- //
- // Author: Christoph Blume (C.Blume@gsi.de) 20/07/99
+ // Create the GEANT geometry for the Transition Radiation Detector - Version 2
+ // This version covers the full azimuth.
//
- Float_t xpos, ypos, zpos;
-
// Check that FRAME is there otherwise we have no place where to put the TRD
AliModule* FRAME = gAlice->GetModule("FRAME");
if (!FRAME) return;
// Define the chambers
AliTRD::CreateGeometry();
- // Position the the TRD-sectors only in one TRD-volume in the spaceframe
- xpos = 0.;
- ypos = 0.;
- zpos = 0.;
- gMC->Gspos("TRD ",1,"BTR1",xpos,ypos,zpos,0,"ONLY");
-
}
//_____________________________________________________________________________
void AliTRDv1::CreateMaterials()
{
//
- // Create materials for the Transition Radiation Detector version 1
+ // Create materials for the Transition Radiation Detector version 2
//
AliTRD::CreateMaterials();
}
//_____________________________________________________________________________
-void AliTRDv1::Init()
+void AliTRDv1::Diffusion(Float_t driftlength, Float_t *xyz)
{
//
- // Initialise the Transition Radiation Detector after the geometry is built
+ // Applies the diffusion smearing to the position of a single electron
//
- printf("**************************************"
- " TRD "
- "**************************************\n");
- printf("\n Version 1 of TRD initialing, "
- "with openings for PHOS and RICH\n\n");
+ if ((driftlength > 0) &&
+ (driftlength < kDrThick)) {
+ Float_t driftSqrt = TMath::Sqrt(driftlength);
+ Float_t sigmaT = driftSqrt * fDiffusionT;
+ Float_t sigmaL = driftSqrt * fDiffusionL;
+ xyz[0] = gRandom->Gaus(xyz[0], sigmaL);
+ xyz[1] = gRandom->Gaus(xyz[1], sigmaT);
+ xyz[2] = gRandom->Gaus(xyz[2], sigmaT);
+ }
+ else {
+ xyz[0] = 0.0;
+ xyz[1] = 0.0;
+ xyz[2] = 0.0;
+ }
- AliTRD::Init();
+}
+//_____________________________________________________________________________
+void AliTRDv1::Hits2Digits()
+{
+ //
+ // Creates TRD digits from hits. This procedure includes the following:
+ // - Diffusion
+ // - Gas gain including fluctuations
+ // - Pad-response (simple Gaussian approximation)
+ // - Electronics noise
+ // - Electronics gain
+ // - Digitization
+ // - ADC threshold
+ // The corresponding parameter can be adjusted via the various Set-functions.
+ // If these parameters are not explicitly set, default values are used (see
+ // Init-function).
+ // To produce digits from a root-file with TRD-hits use the
+ // slowDigitsCreate.C macro.
//
- // Check that FRAME is there otherwise we have no place where to
- // put TRD
- AliModule* FRAME=gAlice->GetModule("FRAME");
- if(!FRAME) {
- Error("Ctor","TRD needs FRAME to be present\n");
- exit(1);
- } else
- if(FRAME->IsVersion()!=0) {
- Error("Ctor","FRAME version 0 needed with this version of TRD\n");
- exit(1);
+
+ printf("AliTRDv1::Hits2Digits -- Start creating digits\n");
+
+ ///////////////////////////////////////////////////////////////
+ // Parameter
+ ///////////////////////////////////////////////////////////////
+
+ // Converts number of electrons to fC
+ const Float_t el2fC = 1.602E-19 * 1.0E15;
+
+ ///////////////////////////////////////////////////////////////
+
+ Int_t nBytes = 0;
+
+ AliTRDhit *TRDhit;
+
+ // Get the pointer to the hit tree
+ TTree *HitTree = gAlice->TreeH();
+ // Get the pointer to the digits tree
+ TTree *DigitsTree = gAlice->TreeD();
+
+ // Get the number of entries in the hit tree
+ // (Number of primary particles creating a hit somewhere)
+ Int_t nTrack = (Int_t) HitTree->GetEntries();
+
+ Int_t chamBeg = 0;
+ Int_t chamEnd = kNcham;
+ if (fSensChamber) chamEnd = chamBeg = fSensChamber;
+ Int_t planBeg = 0;
+ Int_t planEnd = kNplan;
+ if (fSensPlane) planEnd = planBeg = fSensPlane;
+ Int_t sectBeg = 0;
+ Int_t sectEnd = kNsect;
+ if (fSensSector) sectEnd = sectBeg = fSensSector;
+
+ // Loop through all the chambers
+ for (Int_t icham = chamBeg; icham < chamEnd; icham++) {
+ for (Int_t iplan = planBeg; iplan < planEnd; iplan++) {
+ for (Int_t isect = sectBeg; isect < sectEnd; isect++) {
+
+ Int_t nDigits = 0;
+
+ printf("AliTRDv1::Hits2Digits -- Digitizing chamber %d, plane %d, sector %d\n"
+ ,icham+1,iplan+1,isect+1);
+
+ // Create a detector matrix to keep the signal and track numbers
+ AliTRDmatrix *matrix = new AliTRDmatrix(fRowMax[iplan][icham][isect]
+ ,fColMax[iplan]
+ ,fTimeMax
+ ,isect+1,icham+1,iplan+1);
+
+ // Loop through all entries in the tree
+ for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) {
+
+ gAlice->ResetHits();
+ nBytes += HitTree->GetEvent(iTrack);
+
+ // Get the number of hits in the TRD created by this particle
+ Int_t nHit = fHits->GetEntriesFast();
+
+ // Loop through the TRD hits
+ for (Int_t iHit = 0; iHit < nHit; iHit++) {
+
+ if (!(TRDhit = (AliTRDhit *) fHits->UncheckedAt(iHit)))
+ continue;
+
+ Float_t x = TRDhit->fX;
+ Float_t y = TRDhit->fY;
+ Float_t z = TRDhit->fZ;
+ Float_t q = TRDhit->fQ;
+ Int_t track = TRDhit->fTrack;
+ Int_t plane = TRDhit->fPlane;
+ Int_t sector = TRDhit->fSector;
+ Int_t chamber = TRDhit->fChamber;
+
+ if ((sector != isect+1) ||
+ (plane != iplan+1) ||
+ (chamber != icham+1))
+ continue;
+
+ // Rotate the sectors on top of each other
+ Float_t phi = 2.0 * kPI / (Float_t) kNsect
+ * ((Float_t) sector - 0.5);
+ Float_t xRot = -x * TMath::Cos(phi) + y * TMath::Sin(phi);
+ Float_t yRot = x * TMath::Sin(phi) + y * TMath::Cos(phi);
+ Float_t zRot = z;
+
+ // The hit position in pad coordinates (center pad)
+ // The pad row (z-direction)
+ Int_t rowH = (Int_t) ((zRot - fRow0[iplan][icham][isect]) / fRowPadSize);
+ // The pad column (rphi-direction)
+ Int_t colH = (Int_t) ((yRot - fCol0[iplan] ) / fColPadSize);
+ // The time bucket
+ Int_t timeH = (Int_t) ((xRot - fTime0[iplan] ) / fTimeBinSize);
+
+ // Array to sum up the signal in a box surrounding the
+ // hit postition
+ const Int_t timeBox = 5;
+ const Int_t colBox = 7;
+ const Int_t rowBox = 5;
+ Float_t signalSum[rowBox][colBox][timeBox];
+ for (Int_t iRow = 0; iRow < rowBox; iRow++ ) {
+ for (Int_t iCol = 0; iCol < colBox; iCol++ ) {
+ for (Int_t iTime = 0; iTime < timeBox; iTime++) {
+ signalSum[iRow][iCol][iTime] = 0;
+ }
+ }
+ }
+
+ // Loop over all electrons of this hit
+ Int_t nEl = (Int_t) q;
+ for (Int_t iEl = 0; iEl < nEl; iEl++) {
+
+ // Apply the diffusion smearing
+ Float_t driftlength = xRot - fTime0[iplan];
+ Float_t xyz[3];
+ xyz[0] = xRot;
+ xyz[1] = yRot;
+ xyz[2] = zRot;
+ Diffusion(driftlength,xyz);
+
+ // At this point absorption effects that depend on the
+ // driftlength could be taken into account.
+
+ // The electron position and the distance to the hit position
+ // in pad units
+ // The pad row (z-direction)
+ Int_t rowE = (Int_t) ((xyz[2] - fRow0[iplan][icham][isect]) / fRowPadSize);
+ Int_t rowD = rowH - rowE;
+ // The pad column (rphi-direction)
+ Int_t colE = (Int_t) ((xyz[1] - fCol0[iplan] ) / fColPadSize);
+ Int_t colD = colH - colE;
+ // The time bucket
+ Int_t timeE = (Int_t) ((xyz[0] - fTime0[iplan] ) / fTimeBinSize);
+ Int_t timeD = timeH - timeE;
+
+ // Apply the gas gain including fluctuations
+ Int_t signal = (Int_t) (-fGasGain * TMath::Log(gRandom->Rndm()));
+
+ // The distance of the electron to the center of the pad
+ // in units of pad width
+ Float_t dist = (xyz[1] - fCol0[iplan] - (colE + 0.5) * fColPadSize)
+ / fColPadSize;
+
+ // Sum up the signal in the different pixels
+ // and apply the pad response
+ Int_t rowIdx = rowD + (Int_t) ( rowBox / 2);
+ Int_t colIdx = colD + (Int_t) ( colBox / 2);
+ Int_t timeIdx = timeD + (Int_t) (timeBox / 2);
+ signalSum[rowIdx][colIdx-1][timeIdx] += PadResponse(dist-1.) * signal;
+ signalSum[rowIdx][colIdx ][timeIdx] += PadResponse(dist ) * signal;
+ signalSum[rowIdx][colIdx+1][timeIdx] += PadResponse(dist+1.) * signal;
+
+ }
+
+ // Add the padcluster to the detector matrix
+ for (Int_t iRow = 0; iRow < rowBox; iRow++ ) {
+ for (Int_t iCol = 0; iCol < colBox; iCol++ ) {
+ for (Int_t iTime = 0; iTime < timeBox; iTime++) {
+
+ Int_t rowB = rowH + iRow - (Int_t) ( rowBox / 2);
+ Int_t colB = colH + iCol - (Int_t) ( colBox / 2);
+ Int_t timeB = timeH + iTime - (Int_t) (timeBox / 2);
+
+ Float_t signalB = signalSum[iRow][iCol][iTime];
+ if (signalB > 0.0) {
+ matrix->AddSignal(rowB,colB,timeB,signalB);
+ if (!(matrix->AddTrack(rowB,colB,timeB,track)))
+ printf(" More than three tracks in a pixel!\n");
+ }
+
+ }
+ }
+ }
+
+ }
+
+ }
+
+ // Create the hits for this chamber
+ for (Int_t iRow = 0; iRow < fRowMax[iplan][icham][isect]; iRow++ ) {
+ for (Int_t iCol = 0; iCol < fColMax[iplan] ; iCol++ ) {
+ for (Int_t iTime = 0; iTime < fTimeMax ; iTime++) {
+
+ Float_t signalAmp = matrix->GetSignal(iRow,iCol,iTime);
+
+ // Add the noise
+ signalAmp = TMath::Max(gRandom->Gaus(signalAmp,fNoise),(Float_t) 0.0);
+ // Convert to fC
+ signalAmp *= el2fC;
+ // Convert to mV
+ signalAmp *= fChipGain;
+ // Convert to ADC counts
+ Int_t adc = (Int_t) (signalAmp * (fADCoutRange / fADCinRange));
+
+ // Apply threshold on ADC value
+ if (adc > fADCthreshold) {
+
+ Int_t trackSave[3];
+ for (Int_t ii = 0; ii < 3; ii++) {
+ trackSave[ii] = matrix->GetTrack(iRow,iCol,iTime,ii);
+ }
+
+ Int_t digits[7];
+ digits[0] = matrix->GetSector();
+ digits[1] = matrix->GetChamber();
+ digits[2] = matrix->GetPlane();
+ digits[3] = iRow;
+ digits[4] = iCol;
+ digits[5] = iTime;
+ digits[6] = adc;
+
+ // Add this digit to the TClonesArray
+ AddDigit(trackSave,digits);
+ nDigits++;
+
+ }
+
+ }
+ }
+ }
+
+ printf("AliTRDv1::Hits2Digits -- Number of digits found: %d\n",nDigits);
+
+ // Clean up
+ delete matrix;
+
+ }
}
+ }
- for (Int_t i = 0; i < 80; i++) printf("*");
- printf("\n");
-
- // Identifier of the sensitive volume (amplification region)
- fIdSens = gMC->VolId("UL06");
+ // Fill the digits-tree
+ printf("AliTRDv1::Hits2Digits -- Fill the digits tree\n");
+ DigitsTree->Fill();
+
+}
+
+//_____________________________________________________________________________
+void AliTRDv1::Digits2Clusters()
+{
+
+ //
+ // Method to convert AliTRDdigits created by AliTRDv1::Hits2Digits()
+ // into AliTRDclusters
+ // To produce cluster from a root-file with TRD-digits use the
+ // slowClusterCreate.C macro.
+ //
+
+ printf("AliTRDv1::Digits2Clusters -- Start creating clusters\n");
+
+ AliTRDdigit *TRDdigit;
+ TClonesArray *TRDDigits;
+
+ // Parameters
+ Float_t maxThresh = fClusMaxThresh; // threshold value for maximum
+ Float_t signalThresh = fClusSigThresh; // threshold value for digit signal
+ Int_t clusteringMethod = fClusMethod; // clustering method option (for testing)
+
+ const Float_t epsilon = 0.01; // iteration limit for unfolding procedure
+
+ // Get the pointer to the digits tree
+ TTree *DigitTree = gAlice->TreeD();
+ // Get the pointer to the cluster tree
+ TTree *ClusterTree = gAlice->TreeD();
+
+ // Get the pointer to the digits container
+ TRDDigits = Digits();
+
+ Int_t chamBeg = 0;
+ Int_t chamEnd = kNcham;
+ if (fSensChamber) chamEnd = chamBeg = fSensChamber;
+ Int_t planBeg = 0;
+ Int_t planEnd = kNplan;
+ if (fSensPlane) planEnd = planBeg = fSensPlane;
+ Int_t sectBeg = 0;
+ Int_t sectEnd = kNsect;
+ if (fSensSector) sectEnd = sectBeg = fSensSector;
+
+ // Import the digit tree
+ gAlice->ResetDigits();
+ Int_t nbytes;
+ nbytes += DigitTree->GetEvent(1);
+
+ // Get the number of digits in the detector
+ Int_t nTRDDigits = TRDDigits->GetEntriesFast();
+
+ // *** Start clustering *** in every chamber
+ for (Int_t icham = chamBeg; icham < chamEnd; icham++) {
+ for (Int_t iplan = planBeg; iplan < planEnd; iplan++) {
+ for (Int_t isect = sectBeg; isect < sectEnd; isect++) {
+
+ Int_t nClusters = 0;
+ printf("AliTRDv1::Digits2Clusters -- Finding clusters in chamber %d, plane %d, sector %d\n"
+ ,icham+1,iplan+1,isect+1);
+
+ // Create a detector matrix to keep maxima
+ AliTRDmatrix *digitMatrix = new AliTRDmatrix(fRowMax[iplan][icham][isect]
+ ,fColMax[iplan]
+ ,fTimeMax,isect+1
+ ,icham+1,iplan+1);
+ // Create a matrix to contain maximum flags
+ AliTRDmatrix *maximaMatrix = new AliTRDmatrix(fRowMax[iplan][icham][isect]
+ ,fColMax[iplan]
+ ,fTimeMax
+ ,isect+1,icham+1,iplan+1);
+
+ // Loop through all TRD digits
+ for (Int_t iTRDDigits = 0; iTRDDigits < nTRDDigits; iTRDDigits++) {
+
+ // Get the information for this digit
+ TRDdigit = (AliTRDdigit*) TRDDigits->UncheckedAt(iTRDDigits);
+ Int_t signal = TRDdigit->fSignal;
+ Int_t sector = TRDdigit->fSector;
+ Int_t chamber = TRDdigit->fChamber;
+ Int_t plane = TRDdigit->fPlane;
+ Int_t row = TRDdigit->fRow;
+ Int_t col = TRDdigit->fCol;
+ Int_t time = TRDdigit->fTime;
+
+ Int_t track[3];
+ for (Int_t iTrack = 0; iTrack < 3; iTrack++) {
+ track[iTrack] = TRDdigit->AliDigit::fTracks[iTrack];
+ }
+
+ if ((sector != isect+1) ||
+ (plane != iplan+1) ||
+ (chamber != icham+1))
+ continue;
+
+ // Fill the detector matrix
+ if (signal > signalThresh) {
+ digitMatrix->SetSignal(row,col,time,signal);
+ for (Int_t iTrack = 0; iTrack < 3; iTrack++) {
+ if (track[iTrack] > 0) {
+ digitMatrix->AddTrack(row,col,time,track[iTrack]);
+ }
+ }
+ }
+
+ }
+
+ // Loop chamber and find maxima in digitMatrix
+ for (Int_t row = 0; row < fRowMax[iplan][icham][isect]; row++) {
+ for (Int_t col = 1; col < fColMax[iplan] ; col++) {
+ for (Int_t time = 0; time < fTimeMax ; time++) {
+
+ if (digitMatrix->GetSignal(row,col,time)
+ < digitMatrix->GetSignal(row,col - 1,time)) {
+ // really maximum?
+ if (col > 1) {
+ if (digitMatrix->GetSignal(row,col - 2,time)
+ < digitMatrix->GetSignal(row,col - 1,time)) {
+ // yes, so set maximum flag
+ maximaMatrix->SetSignal(row,col - 1,time,1);
+ }
+ else maximaMatrix->SetSignal(row,col - 1,time,0);
+ }
+ }
+
+ } // time
+ } // col
+ } // row
+
+ // now check maxima and calculate cluster position
+ for (Int_t row = 0; row < fRowMax[iplan][icham][isect]; row++) {
+ for (Int_t col = 1; col < fColMax[iplan] ; col++) {
+ for (Int_t time = 0; time < fTimeMax ; time++) {
+
+ if ((maximaMatrix->GetSignal(row,col,time) > 0)
+ && (digitMatrix->GetSignal(row,col,time) > maxThresh)) {
+
+ Int_t clusters[5] = {0}; // cluster-object data
+
+ Float_t ratio = 0; // ratio resulting from unfolding
+ Float_t padSignal[5] = {0}; // signals on max and neighbouring pads
+ Float_t clusterSignal[3] = {0}; // signals from cluster
+ Float_t clusterPos[3] = {0}; // cluster in ALICE refFrame coords
+ Float_t clusterPads[6] = {0}; // cluster pad info
+
+ // setting values
+ clusters[0] = isect+1; // = isect ????
+ clusters[1] = icham+1; // = ichamber ????
+ clusters[2] = iplan+1; // = iplane ????
+ clusters[3] = time;
+
+ clusterPads[0] = icham+1;
+ clusterPads[1] = isect+1;
+ clusterPads[2] = iplan+1;
+
+ for (Int_t iPad = 0; iPad < 3; iPad++) {
+ clusterSignal[iPad] = digitMatrix->GetSignal(row,col-1+iPad,time);
+ }
+
+ // neighbouring maximum on right side?
+ if (col < fColMax[iplan] - 2) {
+ if (maximaMatrix->GetSignal(row,col + 2,time) > 0) {
+ for (Int_t iPad = 0; iPad < 5; iPad++) {
+ padSignal[iPad] = digitMatrix->GetSignal(row,col-1+iPad,time);
+ }
+
+ // unfold:
+ ratio = Unfold(epsilon, padSignal);
+
+ // set signal on overlapping pad to ratio
+ clusterSignal[2] *= ratio;
+ }
+ }
+
+ switch (clusteringMethod) {
+ case 1:
+ // method 1: simply center of mass
+ clusterPads[3] = row + 0.5;
+ clusterPads[4] = col - 0.5 + (clusterSignal[2] - clusterSignal[0]) /
+ (clusterSignal[1] + clusterSignal[2] + clusterSignal[3]);
+ clusterPads[5] = time + 0.5;
+
+ nClusters++;
+ break;
+ case 2:
+ // method 2: integral gauss fit on 3 pads
+ TH1F *hPadCharges = new TH1F("hPadCharges", "Charges on center 3 pads"
+ , 5, -1.5, 3.5);
+ for (Int_t iCol = -1; iCol <= 3; iCol++) {
+ if (clusterSignal[iCol] < 1) clusterSignal[iCol] = 1;
+ hPadCharges->Fill(iCol, clusterSignal[iCol]);
+ }
+ hPadCharges->Fit("gaus", "IQ", "SAME", -0.5, 2.5);
+ TF1 *fPadChargeFit = hPadCharges->GetFunction("gaus");
+ Double_t colMean = fPadChargeFit->GetParameter(1);
+
+ clusterPads[3] = row + 0.5;
+ clusterPads[4] = col - 1.5 + colMean;
+ clusterPads[5] = time + 0.5;
+
+ delete hPadCharges;
+
+ nClusters++;
+ break;
+ }
+
+ Float_t clusterCharge = clusterSignal[0]
+ + clusterSignal[1]
+ + clusterSignal[2];
+ clusters[4] = (Int_t)clusterCharge;
+
+ Int_t trackSave[3];
+ for (Int_t iTrack = 0; iTrack < 3; iTrack++) {
+ trackSave[iTrack] = digitMatrix->GetTrack(row,col,time,iTrack);
+ }
+
+ // Calculate cluster position in ALICE refFrame coords
+ // and set array clusterPos to calculated values
+ Pads2XYZ(clusterPads, clusterPos);
+
+ // Add cluster to reconstruction tree
+ AddCluster(trackSave,clusters,clusterPos);
+
+ }
- // Identifier of the TRD-spaceframe volumina
- fIdSpace1 = gMC->VolId("B028");
- fIdSpace2 = gMC->VolId("B029");
- fIdSpace3 = gMC->VolId("B030");
+ } // time
+ } // col
+ } // row
+
+ printf("AliTRDv1::Digits2Clusters -- Number of clusters found: %d\n",nClusters);
+
+ delete digitMatrix;
+ delete maximaMatrix;
+
+ } // isect
+ } // iplan
+ } // icham
+
+ // Fill the cluster-tree
+ printf("AliTRDv1::Digits2Clusters -- Total number of clusters found: %d\n"
+ ,fClusters->GetEntries());
+ printf("AliTRDv1::Digits2Clusters -- Fill the cluster tree\n");
+ ClusterTree->Fill();
+
+}
+
+//_____________________________________________________________________________
+void AliTRDv1::Init()
+{
+ //
+ // Initialise Transition Radiation Detector after geometry has been built.
+ // Includes the default settings of all parameter for the digitization.
+ //
+
+ AliTRD::Init();
+
+ printf(" Slow simulator\n");
+ if (fSensPlane)
+ printf(" Only plane %d is sensitive\n",fSensPlane);
+ if (fSensChamber)
+ printf(" Only chamber %d is sensitive\n",fSensChamber);
+ if (fSensSector)
+ printf(" Only sector %d is sensitive\n",fSensSector);
+
+ // First ionization potential (eV) for the gas mixture (90% Xe + 10% CO2)
+ const Float_t kPoti = 12.1;
+ // Maximum energy (50 keV);
+ const Float_t kEend = 50000.0;
+ // Ermilova distribution for the delta-ray spectrum
+ Float_t Poti = TMath::Log(kPoti);
+ Float_t Eend = TMath::Log(kEend);
+ fDeltaE = new TF1("deltae",Ermilova,Poti,Eend,0);
+
+ // Identifier of the sensitive volume (drift region)
+ fIdSens = gMC->VolId("UL05");
// Identifier of the TRD-driftchambers
fIdChamber1 = gMC->VolId("UCIO");
fIdChamber2 = gMC->VolId("UCIM");
fIdChamber3 = gMC->VolId("UCII");
- printf("**************************************"
- " TRD "
- "**************************************\n");
+ // The default parameter for the digitization
+ if (!(fGasGain)) fGasGain = 2.0E3;
+ if (!(fNoise)) fNoise = 3000.;
+ if (!(fChipGain)) fChipGain = 10.;
+ if (!(fADCoutRange)) fADCoutRange = 255.;
+ if (!(fADCinRange)) fADCinRange = 2000.;
+ if (!(fADCthreshold)) fADCthreshold = 1;
+
+ // Transverse and longitudinal diffusion coefficients (Xe/Isobutane)
+ if (!(fDiffusionT)) fDiffusionT = 0.060;
+ if (!(fDiffusionL)) fDiffusionL = 0.017;
+
+ // The default parameter for the clustering
+ if (!(fClusMaxThresh)) fClusMaxThresh = 5.0;
+ if (!(fClusSigThresh)) fClusSigThresh = 2.0;
+ if (!(fClusMethod)) fClusMethod = 1;
+
+ for (Int_t i = 0; i < 80; i++) printf("*");
+ printf("\n");
+
}
//_____________________________________________________________________________
-void AliTRDv1::StepManager()
+Float_t AliTRDv1::PadResponse(Float_t x)
{
//
- // Procedure called at every step in the TRD
- // Fast simulator. If switched on, a hit is produced when a track
- // crosses the border between amplification region and pad plane.
+ // The pad response for the chevron pads.
+ // We use a simple Gaussian approximation which should be good
+ // enough for our purpose.
//
- Int_t vol[3];
- Int_t iIdSens, icSens;
- Int_t iIdSpace, icSpace;
- Int_t iIdChamber, icChamber;
+ // The parameters for the response function
+ const Float_t aa = 0.8872;
+ const Float_t bb = -0.00573;
+ const Float_t cc = 0.454;
+ const Float_t cc2 = cc*cc;
+
+ Float_t pr = aa * (bb + TMath::Exp(-x*x / (2. * cc2)));
+
+ return (pr);
+
+}
+
+//_____________________________________________________________________________
+void AliTRDv1::SetSensPlane(Int_t iplane)
+{
+ //
+ // Defines the hit-sensitive plane (1-6)
+ //
- Int_t secMap1[10] = { 3, 7, 8, 9, 10, 11, 2, 1, 18, 17 };
- Int_t secMap2[ 5] = { 16, 15, 14, 13, 12 };
- Int_t secMap3[ 3] = { 5, 6, 4 };
+ if ((iplane < 0) || (iplane > 6)) {
+ printf("Wrong input value: %d\n",iplane);
+ printf("Use standard setting\n");
+ fSensPlane = 0;
+ fSensSelect = 0;
+ return;
+ }
- Float_t hits[4];
+ fSensSelect = 1;
+ fSensPlane = iplane;
- TLorentzVector p;
+}
+
+//_____________________________________________________________________________
+void AliTRDv1::SetSensChamber(Int_t ichamber)
+{
+ //
+ // Defines the hit-sensitive chamber (1-5)
+ //
+
+ if ((ichamber < 0) || (ichamber > 5)) {
+ printf("Wrong input value: %d\n",ichamber);
+ printf("Use standard setting\n");
+ fSensChamber = 0;
+ fSensSelect = 0;
+ return;
+ }
+
+ fSensSelect = 1;
+ fSensChamber = ichamber;
+
+}
+
+//_____________________________________________________________________________
+void AliTRDv1::SetSensSector(Int_t isector)
+{
+ //
+ // Defines the hit-sensitive sector (1-18)
+ //
+
+ if ((isector < 0) || (isector > 18)) {
+ printf("Wrong input value: %d\n",isector);
+ printf("Use standard setting\n");
+ fSensSector = 0;
+ fSensSelect = 0;
+ return;
+ }
+
+ fSensSelect = 1;
+ fSensSector = isector;
+
+}
+
+//_____________________________________________________________________________
+void AliTRDv1::StepManager()
+{
+ //
+ // Called at every step in the Transition Radiation Detector version 2.
+ // Slow simulator. Every charged track produces electron cluster as hits
+ // along its path across the drift volume. The step size is set acording
+ // to Bethe-Bloch. The energy distribution of the delta electrons follows
+ // a spectrum taken from Ermilova et al.
+ //
+
+ Int_t iIdSens, icSens;
+ Int_t iIdSpace, icSpace;
+ Int_t iIdChamber, icChamber;
+ Int_t vol[3];
+ Int_t iPid;
+
+ Float_t hits[4];
+ Float_t random[1];
+ Float_t charge;
+ Float_t aMass;
+
+ Double_t pTot;
+ Double_t qTot;
+ Double_t eDelta;
+ Double_t betaGamma, pp;
+
+ TLorentzVector pos, mom;
TClonesArray &lhits = *fHits;
- // Writing out hits enabled?
- if (!(fHitsOn)) return;
+ const Double_t kBig = 1.0E+12;
+
+ // Ionization energy
+ const Float_t kWion = 22.04;
+ // Maximum energy for e+ e- g for the step-size calculation
+ const Float_t kPTotMax = 0.002;
+ // Plateau value of the energy-loss for electron in xenon
+ // taken from: Allison + Comb, Ann. Rev. Nucl. Sci. (1980), 30, 253
+ //const Double_t kPlateau = 1.70;
+ // the averaged value (26/3/99)
+ const Float_t kPlateau = 1.55;
+ // dN1/dx|min for the gas mixture (90% Xe + 10% CO2)
+ const Float_t kPrim = 48.0;
+ // First ionization potential (eV) for the gas mixture (90% Xe + 10% CO2)
+ const Float_t kPoti = 12.1;
+
+ // Set the maximum step size to a very large number for all
+ // neutral particles and those outside the driftvolume
+ gMC->SetMaxStep(kBig);
+
+ // Use only charged tracks
+ if (( gMC->TrackCharge() ) &&
+ (!gMC->IsTrackStop() ) &&
+ (!gMC->IsTrackDisappeared())) {
- // Use only charged tracks and count them only once per volume
- if (gMC->TrackCharge() &&
- gMC->IsTrackExiting()) {
-
- // Check on sensitive volume
+ // Inside a sensitive volume?
iIdSens = gMC->CurrentVolID(icSens);
if (iIdSens == fIdSens) {
- gMC->TrackPosition(p);
- for (Int_t i = 0; i < 3; i++) hits[i] = p[i];
- // No charge created
- hits[3] = 0;
-
iIdSpace = gMC->CurrentVolOffID(4,icSpace );
iIdChamber = gMC->CurrentVolOffID(1,icChamber);
+ // Calculate the energy of the delta-electrons
+ eDelta = TMath::Exp(fDeltaE->GetRandom()) - kPoti;
+ eDelta = TMath::Max(eDelta,0.0);
+
+ // The number of secondary electrons created
+ qTot = (Double_t) ((Int_t) (eDelta / kWion) + 1);
+
+ // The hit coordinates and charge
+ gMC->TrackPosition(pos);
+ hits[0] = pos[0];
+ hits[1] = pos[1];
+ hits[2] = pos[2];
+ hits[3] = qTot;
+
// The sector number
- if (iIdSpace == fIdSpace1)
- vol[0] = secMap1[icSpace-1];
- else if (iIdSpace == fIdSpace2)
- vol[0] = secMap2[icSpace-1];
- else if (iIdSpace == fIdSpace3)
- vol[0] = secMap3[icSpace-1];
+ Float_t phi = pos[1] != 0 ? TMath::Atan2(pos[0],pos[1]) : (pos[0] > 0 ? 180. : 0.);
+ vol[0] = ((Int_t) (phi / 20)) + 1;
// The chamber number
// 1: outer left
// The plane number
vol[2] = icChamber - TMath::Nint((Float_t) (icChamber / 7)) * 6;
- new(lhits[fNhits++]) AliTRDhit(fIshunt,gAlice->CurrentTrack(),vol,hits);
+ // Check on selected volumes
+ Int_t addthishit = 1;
+ if (fSensSelect) {
+ if ((fSensPlane) && (vol[2] != fSensPlane )) addthishit = 0;
+ if ((fSensChamber) && (vol[1] != fSensChamber)) addthishit = 0;
+ if ((fSensSector) && (vol[0] != fSensSector )) addthishit = 0;
+ }
+
+ // Add this hit
+ if (addthishit) {
+
+ new(lhits[fNhits++]) AliTRDhit(fIshunt,gAlice->CurrentTrack(),vol,hits);
+
+ // The energy loss according to Bethe Bloch
+ gMC->TrackMomentum(mom);
+ pTot = mom.Rho();
+ iPid = gMC->TrackPid();
+ if ( (iPid > 3) ||
+ ((iPid <= 3) && (pTot < kPTotMax))) {
+ aMass = gMC->TrackMass();
+ betaGamma = pTot / aMass;
+ pp = kPrim * BetheBloch(betaGamma);
+ // Take charge > 1 into account
+ charge = gMC->TrackCharge();
+ if (TMath::Abs(charge) > 1) pp = pp * charge*charge;
+ }
+ // Electrons above 20 Mev/c are at the plateau
+ else {
+ pp = kPrim * kPlateau;
+ }
+
+ // Calculate the maximum step size for the next tracking step
+ if (pp > 0) {
+ do
+ gMC->Rndm(random,1);
+ while ((random[0] == 1.) || (random[0] == 0.));
+ gMC->SetMaxStep( - TMath::Log(random[0]) / pp);
+ }
+
+ }
+ else {
+ // set step size to maximal value
+ gMC->SetMaxStep(kBig);
+ }
}
- }
+ }
+
+}
+
+//_____________________________________________________________________________
+Double_t AliTRDv1::BetheBloch(Double_t bg)
+{
+ //
+ // Parametrization of the Bethe-Bloch-curve
+ // The parametrization is the same as for the TPC and is taken from Lehrhaus.
+ //
+
+ // This parameters have been adjusted to averaged values from GEANT
+ const Double_t kP1 = 7.17960e-02;
+ const Double_t kP2 = 8.54196;
+ const Double_t kP3 = 1.38065e-06;
+ const Double_t kP4 = 5.30972;
+ const Double_t kP5 = 2.83798;
+
+ // This parameters have been adjusted to Xe-data found in:
+ // Allison & Cobb, Ann. Rev. Nucl. Sci. (1980), 30, 253
+ //const Double_t kP1 = 0.76176E-1;
+ //const Double_t kP2 = 10.632;
+ //const Double_t kP3 = 3.17983E-6;
+ //const Double_t kP4 = 1.8631;
+ //const Double_t kP5 = 1.9479;
+
+ if (bg > 0) {
+ Double_t yy = bg / TMath::Sqrt(1. + bg*bg);
+ Double_t aa = TMath::Power(yy,kP4);
+ Double_t bb = TMath::Power((1./bg),kP5);
+ bb = TMath::Log(kP3 + bb);
+ return ((kP2 - aa - bb)*kP1 / aa);
+ }
+ else
+ return 0;
}
+
+//_____________________________________________________________________________
+Double_t Ermilova(Double_t *x, Double_t *)
+{
+ //
+ // Calculates the delta-ray energy distribution according to Ermilova.
+ // Logarithmic scale !
+ //
+
+ Double_t energy;
+ Double_t dpos;
+ Double_t dnde;
+
+ Int_t pos1, pos2;
+
+ const Int_t nV = 31;
+
+ Float_t vxe[nV] = { 2.3026, 2.9957, 3.4012, 3.6889, 3.9120
+ , 4.0943, 4.2485, 4.3820, 4.4998, 4.6052
+ , 4.7005, 5.0752, 5.2983, 5.7038, 5.9915
+ , 6.2146, 6.5221, 6.9078, 7.3132, 7.6009
+ , 8.0064, 8.5172, 8.6995, 8.9872, 9.2103
+ , 9.4727, 9.9035,10.3735,10.5966,10.8198
+ ,11.5129 };
+
+ Float_t vye[nV] = { 80.0 , 31.0 , 23.3 , 21.1 , 21.0
+ , 20.9 , 20.8 , 20.0 , 16.0 , 11.0
+ , 8.0 , 6.0 , 5.2 , 4.6 , 4.0
+ , 3.5 , 3.0 , 1.4 , 0.67 , 0.44
+ , 0.3 , 0.18 , 0.12 , 0.08 , 0.056
+ , 0.04 , 0.023, 0.015, 0.011, 0.01
+ , 0.004 };
+
+ energy = x[0];
+
+ // Find the position
+ pos1 = pos2 = 0;
+ dpos = 0;
+ do {
+ dpos = energy - vxe[pos2++];
+ }
+ while (dpos > 0);
+ pos2--;
+ if (pos2 > nV) pos2 = nV;
+ pos1 = pos2 - 1;
+
+ // Differentiate between the sampling points
+ dnde = (vye[pos1] - vye[pos2]) / (vxe[pos2] - vxe[pos1]);
+
+ return dnde;
+
+}
+
+//_____________________________________________________________________________
+void AliTRDv1::Pads2XYZ(Float_t *pads, Float_t *pos)
+{
+ // Method to convert pad coordinates (row,col,time)
+ // into ALICE reference frame coordinates (x,y,z)
+
+ Int_t chamber = (Int_t) pads[0]; // chamber info (1-5)
+ Int_t sector = (Int_t) pads[1]; // sector info (1-18)
+ Int_t plane = (Int_t) pads[2]; // plane info (1-6)
+
+ Int_t icham = chamber - 1; // chamber info (0-4)
+ Int_t isect = sector - 1; // sector info (0-17)
+ Int_t iplan = plane - 1; // plane info (0-5)
+
+ Float_t padRow = pads[3]; // Pad Row position
+ Float_t padCol = pads[4]; // Pad Column position
+ Float_t timeSlice = pads[5]; // Time "position"
+
+ // calculate (x,y) position in rotated chamber
+ Float_t yRot = fCol0[iplan] + padCol * fColPadSize;
+ Float_t xRot = fTime0[iplan] + timeSlice * fTimeBinSize;
+ // calculate z-position:
+ Float_t z = fRow0[iplan][icham][isect] + padRow * fRowPadSize;
+
+ /**
+ rotate chamber back to original position
+ 1. mirror at y-axis, 2. rotate back to position (-phi)
+ / cos(phi) -sin(phi) \ / -1 0 \ / -cos(phi) -sin(phi) \
+ \ sin(phi) cos(phi) / * \ 0 1 / = \ -sin(phi) cos(phi) /
+ **/
+ //Float_t phi = 2*kPI / kNsect * ((Float_t) sector - 0.5);
+ //Float_t x = -xRot * TMath::Cos(phi) - yRot * TMath::Sin(phi);
+ //Float_t y = -xRot * TMath::Sin(phi) + yRot * TMath::Cos(phi);
+ Float_t phi = 2*kPI / kNsect * ((Float_t) sector - 0.5);
+ Float_t x = -xRot * TMath::Cos(phi) + yRot * TMath::Sin(phi);
+ Float_t y = xRot * TMath::Sin(phi) + yRot * TMath::Cos(phi);
+
+ // Setting values
+ pos[0] = x;
+ pos[1] = y;
+ pos[2] = z;
+
+}
+
+//_____________________________________________________________________________
+Float_t AliTRDv1::Unfold(Float_t eps, Float_t* padSignal)
+{
+ // Method to unfold neighbouring maxima.
+ // The charge ratio on the overlapping pad is calculated
+ // until there is no more change within the range given by eps.
+ // The resulting ratio is then returned to the calling method.
+
+ Int_t itStep = 0; // count iteration steps
+
+ Float_t ratio = 0.5; // start value for ratio
+ Float_t prevRatio = 0; // store previous ratio
+
+ Float_t newLeftSignal[3] = {0}; // array to store left cluster signal
+ Float_t newRightSignal[3] = {0}; // array to store right cluster signal
+
+ // start iteration:
+ while ((TMath::Abs(prevRatio - ratio) > eps) && (itStep < 10)) {
+
+ itStep++;
+ prevRatio = ratio;
+
+ // cluster position according to charge ratio
+ Float_t maxLeft = (ratio*padSignal[2] - padSignal[0]) /
+ (padSignal[0] + padSignal[1] + ratio*padSignal[2]);
+ Float_t maxRight = (padSignal[4] - (1-ratio)*padSignal[2]) /
+ ((1-ratio)*padSignal[2] + padSignal[3] + padSignal[4]);
+
+ // set cluster charge ratio
+ Float_t ampLeft = padSignal[1];
+ Float_t ampRight = padSignal[3];
+
+ // apply pad response to parameters
+ newLeftSignal[0] = ampLeft*PadResponse(-1 - maxLeft);
+ newLeftSignal[1] = ampLeft*PadResponse( 0 - maxLeft);
+ newLeftSignal[2] = ampLeft*PadResponse( 1 - maxLeft);
+
+ newRightSignal[0] = ampRight*PadResponse(-1 - maxRight);
+ newRightSignal[1] = ampRight*PadResponse( 0 - maxRight);
+ newRightSignal[2] = ampRight*PadResponse( 1 - maxRight);
+
+ // calculate new overlapping ratio
+ ratio = newLeftSignal[2]/(newLeftSignal[2] + newRightSignal[0]);
+
+ }
+
+ return ratio;
+
+}
+
// Manager and hits classes for set:TRD version 1 //
////////////////////////////////////////////////////////
+#include <TF1.h>
#include "AliTRD.h"
-
+
+// Energy spectrum of the delta-rays
+Double_t Ermilova(Double_t *x, Double_t *par);
+
class AliTRDv1 : public AliTRD {
public:
AliTRDv1() {}
AliTRDv1(const char *name, const char *title);
- virtual ~AliTRDv1() {}
- virtual void CreateGeometry();
- virtual void CreateMaterials();
- virtual Int_t IsVersion() const {return 1;}
- virtual void SetHits(Int_t ihit = 1) { fHitsOn = ihit; };
- virtual void StepManager();
- virtual void Init();
+ virtual ~AliTRDv1();
+ virtual void CreateGeometry();
+ virtual void CreateMaterials();
+ virtual Int_t IsVersion() const { return 1; };
+ virtual void StepManager();
+ virtual void SetSensPlane(Int_t iplane = 0);
+ virtual void SetSensChamber(Int_t ichamber = 0);
+ virtual void SetSensSector(Int_t isector = 0);
+ virtual void Init();
+ virtual void Hits2Digits();
+ virtual void Digits2Clusters();
+
+ virtual void SetGasGain(Float_t gasgain) { fGasGain = gasgain; };
+ virtual void SetNoise(Float_t noise) { fNoise = noise; };
+ virtual void SetChipGain(Float_t chipgain) { fChipGain = chipgain; };
+ virtual void SetADCoutRange(Float_t range) { fADCoutRange = range; };
+ virtual void SetADCinRange(Float_t range) { fADCinRange = range; };
+ virtual void SetADCthreshold(Int_t thresh) { fADCthreshold = thresh; };
+ virtual void SetDiffusionT(Float_t diff) { fDiffusionT = diff; };
+ virtual void SetDiffusionL(Float_t diff) { fDiffusionL = diff; };
+
+ virtual void SetClusMaxThresh(Float_t thresh) { fClusMaxThresh = thresh; };
+ virtual void SetClusSigThresh(Float_t thresh) { fClusSigThresh = thresh; };
+ virtual void SetClusMethod(Int_t meth) { fClusMethod = meth; };
+
+ virtual Float_t GetGasGain() { return fGasGain; };
+ virtual Float_t GetNoise() { return fNoise; };
+ virtual Float_t GetChipGain() { return fChipGain; };
+ virtual Float_t GetADCoutRange() { return fADCoutRange; };
+ virtual Float_t GetADCinRange() { return fADCinRange; };
+ virtual Int_t GetADCthreshold() { return fADCthreshold; };
+ virtual Float_t GetDiffusionT() { return fDiffusionT; };
+ virtual Float_t GetDiffusionL() { return fDiffusionL; };
+
+ virtual Float_t GetClusMaxThresh() { return fClusMaxThresh; };
+ virtual Float_t GetClusSigThresh() { return fClusSigThresh; };
+ virtual Int_t GetClusMethod() { return fClusMethod; };
protected:
- Int_t fIdSens; // Sensitive volume identifier
+ Int_t fIdSens; // Sensitive volume identifier
+
+ Int_t fIdChamber1; // Driftchamber volume identifier
+ Int_t fIdChamber2; //
+ Int_t fIdChamber3; //
+
+ Int_t fSensSelect; // Switch to select only parts of the detector
+ Int_t fSensPlane; // Sensitive detector plane
+ Int_t fSensChamber; // Sensitive detector chamber
+ Int_t fSensSector; // Sensitive detector sector
+
+ Float_t fGasGain; // Gas gain
+ Float_t fNoise; // Electronics noise
+ Float_t fChipGain; // Electronics gain
+ Float_t fADCoutRange; // ADC output range (number of channels)
+ Float_t fADCinRange; // ADC input range (input charge)
+ Int_t fADCthreshold; // ADC threshold in ADC channel
+ Float_t fDiffusionT; // Diffusion in transverse direction
+ Float_t fDiffusionL; // Diffusion in logitudinal direction
- Int_t fIdSpace1; // Spaceframe volume identifier
- Int_t fIdSpace2; //
- Int_t fIdSpace3; //
+ Float_t fClusMaxThresh; // Threshold value for cluster maximum
+ Float_t fClusSigThresh; // Threshold value for cluster signal
+ Int_t fClusMethod; // Clustering method
- Int_t fIdChamber1; // Driftchamber volume identifier
- Int_t fIdChamber2; //
- Int_t fIdChamber3; //
+private:
+ virtual Double_t BetheBloch(Double_t bg);
+ virtual void Diffusion(Float_t driftlength, Float_t *xyz);
+ virtual Float_t PadResponse(Float_t x);
+ virtual void Pads2XYZ(Float_t *pads, Float_t *pos);
+ virtual Float_t Unfold(Float_t eps, Float_t *padSignal);
- Int_t fHitsOn; // Used to switch hits on
-
- ClassDef(AliTRDv1,1) // Transition Radiation Detector version 1
+ TF1 *fDeltaE; // Energy distribution of the delta-electrons
+
+ ClassDef(AliTRDv1,1) // Transition Radiation Detector version 1 (slow simulator)
};
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Author: The ALICE Off-line Project. *
- * Contributors are mentioned in the code where appropriate. *
- * *
- * Permission to use, copy, modify and distribute this software and its *
- * documentation strictly for non-commercial purposes is hereby granted *
- * without fee, provided that the above copyright notice appears in all *
- * copies and that both the copyright notice and this permission notice *
- * appear in the supporting documentation. The authors make no claims *
- * about the suitability of this software for any purpose. It is *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-
-/*
-$Log$
-Revision 1.14 1999/10/04 14:48:07 fca
-Avoid warnings on non-ansi compiler HP-UX CC
-
-Revision 1.13 1999/09/29 09:24:35 fca
-Introduction of the Copyright and cvs Log
-
-*/
-
-///////////////////////////////////////////////////////////////////////////////
-// //
-// Transition Radiation Detector version 2 -- detailed simulation //
-// //
-//Begin_Html
-/*
-<img src="picts/AliTRDv2Class.gif">
-*/
-//End_Html
-// //
-// //
-///////////////////////////////////////////////////////////////////////////////
-
-#include <stdlib.h>
-
-#include <TMath.h>
-#include <TVector.h>
-#include <TRandom.h>
-
-#include "AliTRDv2.h"
-#include "AliTRDmatrix.h"
-#include "AliRun.h"
-#include "AliMC.h"
-#include "AliConst.h"
-
-ClassImp(AliTRDv2)
-
-//_____________________________________________________________________________
-AliTRDv2::AliTRDv2(const char *name, const char *title)
- :AliTRD(name, title)
-{
- //
- // Standard constructor for Transition Radiation Detector version 2
- //
-
- fIdSens = 0;
-
- fIdSpace1 = 0;
- fIdSpace2 = 0;
- fIdSpace3 = 0;
-
- fIdChamber1 = 0;
- fIdChamber2 = 0;
- fIdChamber3 = 0;
-
- fSensSelect = 0;
- fSensPlane = 0;
- fSensChamber = 0;
- fSensSector = 0;
-
- Int_t iplan;
-
- for (iplan = 0; iplan < kNplan; iplan++) {
- for (Int_t icham = 0; icham < kNcham; icham++) {
- fRowMax[iplan][icham] = 0;
- }
- fColMax[iplan] = 0;
- }
- fTimeMax = 0;
-
- fRowPadSize = 0;
- fColPadSize = 0;
- fTimeBinSize = 0;
-
- fGasGain = 0;
- fNoise = 0;
- fChipGain = 0;
- fADCoutRange = 0;
- fADCinRange = 0;
- fADCthreshold = 0;
-
- fDiffusionT = 0;
- fDiffusionL = 0;
-
- fDeltaE = NULL;
-
- SetBufferSize(128000);
-
-}
-
-//_____________________________________________________________________________
-AliTRDv2::~AliTRDv2()
-{
-
- if (fDeltaE) delete fDeltaE;
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::CreateGeometry()
-{
- //
- // Create the GEANT geometry for the Transition Radiation Detector - Version 2
- // This version covers the full azimuth.
- //
- // Author: Christoph Blume (C.Blume@gsi.de) 20/07/99
- //
-
- Float_t xpos, ypos, zpos;
-
- // Check that FRAME is there otherwise we have no place where to put the TRD
- AliModule* FRAME = gAlice->GetModule("FRAME");
- if (!FRAME) return;
-
- // Define the chambers
- AliTRD::CreateGeometry();
-
- // Position the the TRD-sectors in all TRD-volumes in the spaceframe
- xpos = 0.;
- ypos = 0.;
- zpos = 0.;
- gMC->Gspos("TRD ",1,"BTR1",xpos,ypos,zpos,0,"ONLY");
- gMC->Gspos("TRD ",2,"BTR2",xpos,ypos,zpos,0,"ONLY");
- gMC->Gspos("TRD ",3,"BTR3",xpos,ypos,zpos,0,"ONLY");
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::CreateMaterials()
-{
- //
- // Create materials for the Transition Radiation Detector version 2
- //
-
- AliTRD::CreateMaterials();
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::Diffusion(Float_t driftlength, Float_t *xyz)
-{
- //
- // Applies the diffusion smearing to the position of a single electron
- //
-
- if ((driftlength > 0) &&
- (driftlength < kDrThick)) {
- Float_t driftSqrt = TMath::Sqrt(driftlength);
- Float_t sigmaT = driftSqrt * fDiffusionT;
- Float_t sigmaL = driftSqrt * fDiffusionL;
- xyz[0] = gRandom->Gaus(xyz[0], sigmaL);
- xyz[1] = gRandom->Gaus(xyz[1], sigmaT);
- xyz[2] = gRandom->Gaus(xyz[2], sigmaT);
- }
- else {
- xyz[0] = 0.0;
- xyz[1] = 0.0;
- xyz[2] = 0.0;
- }
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::Hits2Digits()
-{
- //
- // Creates TRD digits from hits. This procedure includes the following:
- // - Diffusion
- // - Gas gain including fluctuations
- // - Pad-response (simple Gaussian approximation)
- // - Electronics noise
- // - Electronics gain
- // - Digitization
- // - ADC threshold
- // The corresponding parameter can be adjusted via the various Set-functions.
- // If these parameters are not explicitly set, default values are used (see
- // Init-function).
- // To produce digits from a galice.root file with TRD-hits use the
- // digitsCreate.C macro.
- //
-
- printf(" Start creating digits\n");
-
- ///////////////////////////////////////////////////////////////
- // Parameter
- ///////////////////////////////////////////////////////////////
-
- // Converts number of electrons to fC
- const Float_t el2fC = 1.602E-19 * 1.0E15;
-
- ///////////////////////////////////////////////////////////////
-
- Int_t nBytes = 0;
-
- AliTRDhit *TRDhit;
-
- Int_t iplan;
- Int_t iRow;
-
- // Position of pad 0,0,0
- //
- // chambers seen from the top:
- // +----------------------------+
- // | |
- // | | ^
- // | | rphi|
- // | | |
- // |0 | |
- // +----------------------------+ +------>
- // z
- // chambers seen from the side: ^
- // +----------------------------+ time|
- // | | |
- // |0 | |
- // +----------------------------+ +------>
- // z
- //
- // The pad row (z-direction)
- Float_t row0[kNplan][kNcham];
- for (iplan = 0; iplan < kNplan; iplan++) {
- row0[iplan][0] = -fClengthI[iplan]/2. - fClengthM[iplan] - fClengthO[iplan]
- + kCcthick;
- row0[iplan][1] = -fClengthI[iplan]/2. - fClengthM[iplan]
- + kCcthick;
- row0[iplan][2] = -fClengthI[iplan]/2.
- + kCcthick;
- row0[iplan][3] = fClengthI[iplan]/2.
- + kCcthick;
- row0[iplan][4] = fClengthI[iplan]/2. + fClengthM[iplan]
- + kCcthick;
- }
- // The pad column (rphi-direction)
- Float_t col0[kNplan];
- for (iplan = 0; iplan < kNplan; iplan++) {
- col0[iplan] = -fCwidth[iplan]/2. + kCcthick;
- }
- // The time bucket
- Float_t time0[kNplan];
- for (iplan = 0; iplan < kNplan; iplan++) {
- time0[iplan] = kRmin + kCcframe/2. + kDrZpos - 0.5 * kDrThick
- + iplan * (kCheight + kCspace);
- }
-
- // Get the pointer to the hit tree
- TTree *HitTree = gAlice->TreeH();
- // Get the pointer to the digits tree
- TTree *DigitsTree = gAlice->TreeD();
-
- // Get the number of entries in the hit tree
- // (Number of primary particles creating a hit somewhere)
- Int_t nTrack = (Int_t) HitTree->GetEntries();
-
- Int_t chamBeg = 0;
- Int_t chamEnd = kNcham;
- if (fSensChamber) chamEnd = chamBeg = fSensChamber;
- Int_t planBeg = 0;
- Int_t planEnd = kNplan;
- if (fSensPlane) planEnd = planBeg = fSensPlane;
- Int_t sectBeg = 0;
- Int_t sectEnd = kNsect;
- if (fSensSector) sectEnd = sectBeg = fSensSector;
-
- // Loop through all the chambers
- for (Int_t icham = chamBeg; icham < chamEnd; icham++) {
- for (iplan = planBeg; iplan < planEnd; iplan++) {
- for (Int_t isect = sectBeg; isect < sectEnd; isect++) {
-
- printf(" Digitizing chamber %d, plane %d, sector %d\n"
- ,icham+1,iplan+1,isect+1);
-
- // Create a detector matrix to keep the signal and track numbers
- AliTRDmatrix *matrix = new AliTRDmatrix(fRowMax[iplan][icham]
- ,fColMax[iplan]
- ,fTimeMax
- ,isect+1,icham+1,iplan+1);
-
- // Loop through all entries in the tree
- for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) {
-
- gAlice->ResetHits();
- nBytes += HitTree->GetEvent(iTrack);
-
- // Get the number of hits in the TRD created by this particle
- Int_t nHit = fHits->GetEntriesFast();
-
- // Loop through the TRD hits
- for (Int_t iHit = 0; iHit < nHit; iHit++) {
-
- if (!(TRDhit = (AliTRDhit *) fHits->UncheckedAt(iHit)))
- continue;
-
- Float_t x = TRDhit->fX;
- Float_t y = TRDhit->fY;
- Float_t z = TRDhit->fZ;
- Float_t q = TRDhit->fQ;
- Int_t track = TRDhit->fTrack;
- Int_t plane = TRDhit->fPlane;
- Int_t sector = TRDhit->fSector;
- Int_t chamber = TRDhit->fChamber;
-
- if ((sector != isect+1) ||
- (plane != iplan+1) ||
- (chamber != icham+1))
- continue;
-
- // Rotate the sectors on top of each other
- Float_t phi = 2.0 * kPI / (Float_t) kNsect
- * ((Float_t) sector - 0.5);
- Float_t xRot = -x * TMath::Cos(phi) + y * TMath::Sin(phi);
- Float_t yRot = x * TMath::Sin(phi) + y * TMath::Cos(phi);
- Float_t zRot = z;
-
- // The hit position in pad coordinates (center pad)
- // The pad row (z-direction)
- Int_t rowH = (Int_t) ((zRot - row0[iplan][icham]) / fRowPadSize);
- // The pad column (rphi-direction)
- Int_t colH = (Int_t) ((yRot - col0[iplan] ) / fColPadSize);
- // The time bucket
- Int_t timeH = (Int_t) ((xRot - time0[iplan] ) / fTimeBinSize);
-
- // Array to sum up the signal in a box surrounding the
- // hit postition
- const Int_t timeBox = 5;
- const Int_t colBox = 7;
- const Int_t rowBox = 5;
- Float_t signalSum[rowBox][colBox][timeBox];
- for (iRow = 0; iRow < rowBox; iRow++ ) {
- for (Int_t iCol = 0; iCol < colBox; iCol++ ) {
- for (Int_t iTime = 0; iTime < timeBox; iTime++) {
- signalSum[iRow][iCol][iTime] = 0;
- }
- }
- }
-
- // Loop over all electrons of this hit
- Int_t nEl = (Int_t) q;
- for (Int_t iEl = 0; iEl < nEl; iEl++) {
-
- // Apply the diffusion smearing
- Float_t driftlength = xRot - time0[iplan];
- Float_t xyz[3];
- xyz[0] = xRot;
- xyz[1] = yRot;
- xyz[2] = zRot;
- Diffusion(driftlength,xyz);
-
- // At this point absorption effects that depend on the
- // driftlength could be taken into account.
-
- // The electron position and the distance to the hit position
- // in pad units
- // The pad row (z-direction)
- Int_t rowE = (Int_t) ((xyz[2] - row0[iplan][icham]) / fRowPadSize);
- Int_t rowD = rowH - rowE;
- // The pad column (rphi-direction)
- Int_t colE = (Int_t) ((xyz[1] - col0[iplan] ) / fColPadSize);
- Int_t colD = colH - colE;
- // The time bucket
- Int_t timeE = (Int_t) ((xyz[0] - time0[iplan] ) / fTimeBinSize);
- Int_t timeD = timeH - timeE;
-
- // Apply the gas gain including fluctuations
- Int_t signal = (Int_t) (-fGasGain * TMath::Log(gRandom->Rndm()));
-
- // The distance of the electron to the center of the pad
- // in units of pad width
- Float_t dist = (xyz[1] - col0[iplan] - (colE + 0.5) * fColPadSize)
- / fColPadSize;
-
- // Sum up the signal in the different pixels
- // and apply the pad response
- Int_t rowIdx = rowD + (Int_t) ( rowBox / 2);
- Int_t colIdx = colD + (Int_t) ( colBox / 2);
- Int_t timeIdx = timeD + (Int_t) (timeBox / 2);
- signalSum[rowIdx][colIdx-1][timeIdx] += PadResponse(dist-1.) * signal;
- signalSum[rowIdx][colIdx ][timeIdx] += PadResponse(dist ) * signal;
- signalSum[rowIdx][colIdx+1][timeIdx] += PadResponse(dist+1.) * signal;
-
- }
-
- // Add the padcluster to the detector matrix
- for (iRow = 0; iRow < rowBox; iRow++ ) {
- for (Int_t iCol = 0; iCol < colBox; iCol++ ) {
- for (Int_t iTime = 0; iTime < timeBox; iTime++) {
-
- Int_t rowB = rowH + iRow - (Int_t) ( rowBox / 2);
- Int_t colB = colH + iCol - (Int_t) ( colBox / 2);
- Int_t timeB = timeH + iTime - (Int_t) (timeBox / 2);
-
- Float_t signalB = signalSum[iRow][iCol][iTime];
- if (signalB > 0.0) {
- matrix->AddSignal(rowB,colB,timeB,signalB);
- if (!(matrix->AddTrack(rowB,colB,timeB,track)))
- printf("More than three tracks in a pixel!\n");
- }
-
- }
- }
- }
-
- }
-
- }
-
- // Create the hits for this chamber
- for (Int_t iRow = 0; iRow < fRowMax[iplan][icham]; iRow++ ) {
- for (Int_t iCol = 0; iCol < fColMax[iplan] ; iCol++ ) {
- for (Int_t iTime = 0; iTime < fTimeMax ; iTime++) {
-
- Float_t signalAmp = matrix->GetSignal(iRow,iCol,iTime);
-
- // Add the noise
- signalAmp = TMath::Max(gRandom->Gaus(signalAmp,fNoise),(Float_t) 0.0);
- // Convert to fC
- signalAmp *= el2fC;
- // Convert to mV
- signalAmp *= fChipGain;
- // Convert to ADC counts
- Int_t adc = (Int_t) (signalAmp * (fADCoutRange / fADCinRange));
-
- // Apply threshold on ADC value
- if (adc > fADCthreshold) {
-
- Int_t trackSave[3];
- for (Int_t ii = 0; ii < 3; ii++) {
- trackSave[ii] = matrix->GetTrack(iRow,iCol,iTime,ii);
- }
-
- Int_t digits[7];
- digits[0] = matrix->GetSector();
- digits[1] = matrix->GetChamber();
- digits[2] = matrix->GetPlane();
- digits[3] = iRow;
- digits[4] = iCol;
- digits[5] = iTime;
- digits[6] = adc;
-
- // Add this digit to the TClonesArray
- AddDigit(trackSave,digits);
-
- }
-
- }
- }
- }
-
- // Clean up
- delete matrix;
-
- }
- }
- }
-
- // Fill the digits-tree
- DigitsTree->Fill();
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::Init()
-{
- //
- // Initialise Transition Radiation Detector after geometry has been built.
- // Includes the default settings of all parameter for the digitization.
- //
-
- printf("**************************************"
- " TRD "
- "**************************************\n");
- printf("\n Version 2 of TRD initialing, "
- "symmetric TRD\n\n");
-
- AliTRD::Init();
-
-
- //
- // Check that FRAME is there otherwise we have no place where to
- // put TRD
- AliModule* FRAME=gAlice->GetModule("FRAME");
- if(!FRAME) {
- Error("Ctor","TRD needs FRAME to be present\n");
- exit(1);
- } else
- if(FRAME->IsVersion()!=1) {
- Error("Ctor","FRAME version 1 needed with this version of TRD\n");
- exit(1);
- }
-
- if (fSensPlane)
- printf(" Only plane %d is sensitive\n",fSensPlane);
- if (fSensChamber)
- printf(" Only chamber %d is sensitive\n",fSensChamber);
- if (fSensSector)
- printf(" Only sector %d is sensitive\n",fSensSector);
-
- for (Int_t i = 0; i < 80; i++) printf("*");
- printf("\n");
-
- // First ionization potential (eV) for the gas mixture (90% Xe + 10% CO2)
- const Float_t kPoti = 12.1;
- // Maximum energy (50 keV);
- const Float_t kEend = 50000.0;
- // Ermilova distribution for the delta-ray spectrum
- Float_t Poti = TMath::Log(kPoti);
- Float_t Eend = TMath::Log(kEend);
- fDeltaE = new TF1("deltae",Ermilova,Poti,Eend,0);
-
- // Identifier of the sensitive volume (drift region)
- fIdSens = gMC->VolId("UL05");
-
- // Identifier of the TRD-spaceframe volumina
- fIdSpace1 = gMC->VolId("B028");
- fIdSpace2 = gMC->VolId("B029");
- fIdSpace3 = gMC->VolId("B030");
-
- // Identifier of the TRD-driftchambers
- fIdChamber1 = gMC->VolId("UCIO");
- fIdChamber2 = gMC->VolId("UCIM");
- fIdChamber3 = gMC->VolId("UCII");
-
- // The default pad dimensions
- if (!(fRowPadSize)) fRowPadSize = 4.5;
- if (!(fColPadSize)) fColPadSize = 1.0;
- if (!(fTimeBinSize)) fTimeBinSize = 0.1;
-
- // The maximum number of pads
- for (Int_t iplan = 0; iplan < kNplan; iplan++) {
- // Rows
- fRowMax[iplan][0] = 1 + TMath::Nint((fClengthO[iplan] - 2. * kCcthick)
- / fRowPadSize - 0.5);
- fRowMax[iplan][1] = 1 + TMath::Nint((fClengthM[iplan] - 2. * kCcthick)
- / fRowPadSize - 0.5);
- fRowMax[iplan][2] = 1 + TMath::Nint((fClengthI[iplan] - 2. * kCcthick)
- / fRowPadSize - 0.5);
- fRowMax[iplan][3] = 1 + TMath::Nint((fClengthM[iplan] - 2. * kCcthick)
- / fRowPadSize - 0.5);
- fRowMax[iplan][4] = 1 + TMath::Nint((fClengthO[iplan] - 2. * kCcthick)
- / fRowPadSize - 0.5);
- // Columns
- fColMax[iplan] = 1 + TMath::Nint((fCwidth[iplan] - 2. * kCcthick)
- / fColPadSize - 0.5);
- }
- // Time buckets
- fTimeMax = 1 + TMath::Nint(kDrThick / fTimeBinSize - 0.5);
-
- // The default parameter for the digitization
- if (!(fGasGain)) fGasGain = 2.0E3;
- if (!(fNoise)) fNoise = 3000.;
- if (!(fChipGain)) fChipGain = 10.;
- if (!(fADCoutRange)) fADCoutRange = 255.;
- if (!(fADCinRange)) fADCinRange = 2000.;
- if (!(fADCthreshold)) fADCthreshold = 0;
-
- // Transverse and longitudinal diffusion coefficients (Xe/Isobutane)
- if (!(fDiffusionT)) fDiffusionT = 0.060;
- if (!(fDiffusionL)) fDiffusionL = 0.017;
-
- printf("**************************************"
- " TRD "
- "**************************************\n");
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::MakeBranch(Option_t* option)
-{
- //
- // Create Tree branches for the TRD digits.
- //
-
- Int_t buffersize = 4000;
- Char_t branchname[10];
-
- sprintf(branchname,"%s",GetName());
-
- AliDetector::MakeBranch(option);
-
- Char_t *D = strstr(option,"D");
- if (fDigits && gAlice->TreeD() && D) {
- gAlice->TreeD()->Branch(branchname,&fDigits, buffersize);
- printf("Making Branch %s for digits\n",branchname);
- }
-
-}
-
-//_____________________________________________________________________________
-Float_t AliTRDv2::PadResponse(Float_t x)
-{
- //
- // The pad response for the chevron pads.
- // We use a simple Gaussian approximation which should be good
- // enough for our purpose.
- //
-
- // The parameters for the response function
- const Float_t aa = 0.8872;
- const Float_t bb = -0.00573;
- const Float_t cc = 0.454;
- const Float_t cc2 = cc*cc;
-
- Float_t pr = aa * (bb + TMath::Exp(-x*x / (2. * cc2)));
-
- //TF1 *funPR = new TF1("funPR","[0]*([1]+exp(-x*x /(2.*[2])))",-3,3);
- //funPR->SetParameter(0,aa );
- //funPR->SetParameter(1,bb );
- //funPR->SetParameter(2,cc2);
- //
- //Float_t pr = funPR->Eval(distance);
- //
- //delete funPR;
-
- return (pr);
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::SetSensPlane(Int_t iplane)
-{
- //
- // Defines the hit-sensitive plane (1-6)
- //
-
- if ((iplane < 0) || (iplane > 6)) {
- printf("Wrong input value: %d\n",iplane);
- printf("Use standard setting\n");
- fSensPlane = 0;
- fSensSelect = 0;
- return;
- }
-
- fSensSelect = 1;
- fSensPlane = iplane;
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::SetSensChamber(Int_t ichamber)
-{
- //
- // Defines the hit-sensitive chamber (1-5)
- //
-
- if ((ichamber < 0) || (ichamber > 5)) {
- printf("Wrong input value: %d\n",ichamber);
- printf("Use standard setting\n");
- fSensChamber = 0;
- fSensSelect = 0;
- return;
- }
-
- fSensSelect = 1;
- fSensChamber = ichamber;
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::SetSensSector(Int_t isector)
-{
- //
- // Defines the hit-sensitive sector (1-18)
- //
-
- if ((isector < 0) || (isector > 18)) {
- printf("Wrong input value: %d\n",isector);
- printf("Use standard setting\n");
- fSensSector = 0;
- fSensSelect = 0;
- return;
- }
-
- fSensSelect = 1;
- fSensSector = isector;
-
-}
-
-//_____________________________________________________________________________
-void AliTRDv2::StepManager()
-{
- //
- // Called at every step in the Transition Radiation Detector version 2.
- // Slow simulator. Every charged track produces electron cluster as hits
- // along its path across the drift volume. The step size is set acording
- // to Bethe-Bloch. The energy distribution of the delta electrons follows
- // a spectrum taken from Ermilova et al.
- //
-
- Int_t iIdSens, icSens;
- Int_t iIdSpace, icSpace;
- Int_t iIdChamber, icChamber;
- Int_t vol[3];
- Int_t iPid;
-
- Int_t secMap1[10] = { 3, 7, 8, 9, 10, 11, 2, 1, 18, 17 };
- Int_t secMap2[ 5] = { 16, 15, 14, 13, 12 };
- Int_t secMap3[ 3] = { 5, 6, 4 };
-
- Float_t hits[4];
- Float_t random[1];
- Float_t charge;
- Float_t aMass;
-
- Double_t pTot;
- Double_t qTot;
- Double_t eDelta;
- Double_t betaGamma, pp;
-
- TLorentzVector pos, mom;
- TClonesArray &lhits = *fHits;
-
- const Double_t kBig = 1.0E+12;
-
- // Ionization energy
- const Float_t kWion = 22.04;
- // Maximum energy for e+ e- g for the step-size calculation
- const Float_t kPTotMax = 0.002;
- // Plateau value of the energy-loss for electron in xenon
- // taken from: Allison + Comb, Ann. Rev. Nucl. Sci. (1980), 30, 253
- //const Double_t kPlateau = 1.70;
- // the averaged value (26/3/99)
- const Float_t kPlateau = 1.55;
- // dN1/dx|min for the gas mixture (90% Xe + 10% CO2)
- const Float_t kPrim = 48.0;
- // First ionization potential (eV) for the gas mixture (90% Xe + 10% CO2)
- const Float_t kPoti = 12.1;
-
- // Set the maximum step size to a very large number for all
- // neutral particles and those outside the driftvolume
- gMC->SetMaxStep(kBig);
-
- // Use only charged tracks
- if (( gMC->TrackCharge() ) &&
- (!gMC->IsTrackStop() ) &&
- (!gMC->IsTrackDisappeared())) {
-
- // Inside a sensitive volume?
- iIdSens = gMC->CurrentVolID(icSens);
- if (iIdSens == fIdSens) {
-
- iIdSpace = gMC->CurrentVolOffID(4,icSpace );
- iIdChamber = gMC->CurrentVolOffID(1,icChamber);
-
- // Calculate the energy of the delta-electrons
- eDelta = TMath::Exp(fDeltaE->GetRandom()) - kPoti;
- eDelta = TMath::Max(eDelta,0.0);
-
- // The number of secondary electrons created
- qTot = (Double_t) ((Int_t) (eDelta / kWion) + 1);
-
- // The hit coordinates and charge
- gMC->TrackPosition(pos);
- hits[0] = pos[0];
- hits[1] = pos[1];
- hits[2] = pos[2];
- hits[3] = qTot;
-
- // The sector number
- if (iIdSpace == fIdSpace1)
- vol[0] = secMap1[icSpace-1];
- else if (iIdSpace == fIdSpace2)
- vol[0] = secMap2[icSpace-1];
- else if (iIdSpace == fIdSpace3)
- vol[0] = secMap3[icSpace-1];
-
- // The chamber number
- // 1: outer left
- // 2: middle left
- // 3: inner
- // 4: middle right
- // 5: outer right
- if (iIdChamber == fIdChamber1)
- vol[1] = (hits[2] < 0 ? 1 : 5);
- else if (iIdChamber == fIdChamber2)
- vol[1] = (hits[2] < 0 ? 2 : 4);
- else if (iIdChamber == fIdChamber3)
- vol[1] = 3;
-
- // The plane number
- vol[2] = icChamber - TMath::Nint((Float_t) (icChamber / 7)) * 6;
-
- // Check on selected volumes
- Int_t addthishit = 1;
- if (fSensSelect) {
- if ((fSensPlane) && (vol[2] != fSensPlane )) addthishit = 0;
- if ((fSensChamber) && (vol[1] != fSensChamber)) addthishit = 0;
- if ((fSensSector) && (vol[0] != fSensSector )) addthishit = 0;
- }
-
- // Add this hit
- if (addthishit) {
-
- new(lhits[fNhits++]) AliTRDhit(fIshunt,gAlice->CurrentTrack(),vol,hits);
-
- // The energy loss according to Bethe Bloch
- gMC->TrackMomentum(mom);
- pTot = mom.Rho();
- iPid = gMC->TrackPid();
- if ( (iPid > 3) ||
- ((iPid <= 3) && (pTot < kPTotMax))) {
- aMass = gMC->TrackMass();
- betaGamma = pTot / aMass;
- pp = kPrim * BetheBloch(betaGamma);
- // Take charge > 1 into account
- charge = gMC->TrackCharge();
- if (TMath::Abs(charge) > 1) pp = pp * charge*charge;
- }
- // Electrons above 20 Mev/c are at the plateau
- else {
- pp = kPrim * kPlateau;
- }
-
- // Calculate the maximum step size for the next tracking step
- if (pp > 0) {
- do
- gMC->Rndm(random,1);
- while ((random[0] == 1.) || (random[0] == 0.));
- gMC->SetMaxStep( - TMath::Log(random[0]) / pp);
- }
-
- }
- else {
- // set step size to maximal value
- gMC->SetMaxStep(kBig);
- }
-
- }
-
- }
-
-}
-
-//_____________________________________________________________________________
-Double_t AliTRDv2::BetheBloch(Double_t bg)
-{
- //
- // Parametrization of the Bethe-Bloch-curve
- // The parametrization is the same as for the TPC and is taken from Lehrhaus.
- //
-
- // This parameters have been adjusted to averaged values from GEANT
- const Double_t kP1 = 7.17960e-02;
- const Double_t kP2 = 8.54196;
- const Double_t kP3 = 1.38065e-06;
- const Double_t kP4 = 5.30972;
- const Double_t kP5 = 2.83798;
-
- // This parameters have been adjusted to Xe-data found in:
- // Allison & Cobb, Ann. Rev. Nucl. Sci. (1980), 30, 253
- //const Double_t kP1 = 0.76176E-1;
- //const Double_t kP2 = 10.632;
- //const Double_t kP3 = 3.17983E-6;
- //const Double_t kP4 = 1.8631;
- //const Double_t kP5 = 1.9479;
-
- if (bg > 0) {
- Double_t yy = bg / TMath::Sqrt(1. + bg*bg);
- Double_t aa = TMath::Power(yy,kP4);
- Double_t bb = TMath::Power((1./bg),kP5);
- bb = TMath::Log(kP3 + bb);
- return ((kP2 - aa - bb)*kP1 / aa);
- }
- else
- return 0;
-
-}
-
-//_____________________________________________________________________________
-Double_t Ermilova(Double_t *x, Double_t *)
-{
- //
- // Calculates the delta-ray energy distribution according to Ermilova.
- // Logarithmic scale !
- //
-
- Double_t energy;
- Double_t dpos;
- Double_t dnde;
-
- Int_t pos1, pos2;
-
- const Int_t nV = 31;
-
- Float_t vxe[nV] = { 2.3026, 2.9957, 3.4012, 3.6889, 3.9120
- , 4.0943, 4.2485, 4.3820, 4.4998, 4.6052
- , 4.7005, 5.0752, 5.2983, 5.7038, 5.9915
- , 6.2146, 6.5221, 6.9078, 7.3132, 7.6009
- , 8.0064, 8.5172, 8.6995, 8.9872, 9.2103
- , 9.4727, 9.9035,10.3735,10.5966,10.8198
- ,11.5129 };
-
- Float_t vye[nV] = { 80.0 , 31.0 , 23.3 , 21.1 , 21.0
- , 20.9 , 20.8 , 20.0 , 16.0 , 11.0
- , 8.0 , 6.0 , 5.2 , 4.6 , 4.0
- , 3.5 , 3.0 , 1.4 , 0.67 , 0.44
- , 0.3 , 0.18 , 0.12 , 0.08 , 0.056
- , 0.04 , 0.023, 0.015, 0.011, 0.01
- , 0.004 };
-
- energy = x[0];
-
- // Find the position
- pos1 = pos2 = 0;
- dpos = 0;
- do {
- dpos = energy - vxe[pos2++];
- }
- while (dpos > 0);
- pos2--;
- if (pos2 > nV) pos2 = nV;
- pos1 = pos2 - 1;
-
- // Differentiate between the sampling points
- dnde = (vye[pos1] - vye[pos2]) / (vxe[pos2] - vxe[pos1]);
-
- return dnde;
-
-}
+++ /dev/null
-#ifndef TRDv2_H
-#define TRDv2_H
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice */
-
-/* $Id$ */
-
-////////////////////////////////////////////////////////
-// Manager and hits classes for set:TRD version 2 //
-////////////////////////////////////////////////////////
-
-#include <TF1.h>
-#include "AliTRD.h"
-
-// Energy spectrum of the delta-rays
-Double_t Ermilova(Double_t *x, Double_t *par);
-
-class AliTRDv2 : public AliTRD {
-
-public:
- AliTRDv2() {}
- AliTRDv2(const char *name, const char *title);
- virtual ~AliTRDv2();
- virtual void CreateGeometry();
- virtual void CreateMaterials();
- virtual Int_t IsVersion() const {return 2;}
- virtual void MakeBranch(Option_t* option);
- virtual void StepManager();
- virtual void SetSensPlane(Int_t iplane = 0);
- virtual void SetSensChamber(Int_t ichamber = 0);
- virtual void SetSensSector(Int_t isector = 0);
- virtual void Init();
- virtual void Hits2Digits();
-
- virtual void SetRowPadSize(Float_t size) { fRowPadSize = size; };
- virtual void SetColPadSize(Float_t size) { fColPadSize = size; };
- virtual void SetTimeBinSize(Float_t size) { fTimeBinSize = size; };
-
- virtual void SetGasGain(Float_t gasgain) { fGasGain = gasgain; };
- virtual void SetNoise(Float_t noise) { fNoise = noise; };
- virtual void SetChipGain(Float_t chipgain) { fChipGain = chipgain; };
- virtual void SetADCoutRange(Float_t range) { fADCoutRange = range; };
- virtual void SetADCinRange(Float_t range) { fADCinRange = range; };
- virtual void SetADCthreshold(Int_t thresh) { fADCthreshold = thresh; };
- virtual void SetDiffusionT(Float_t diff) { fDiffusionT = diff; };
- virtual void SetDiffusionL(Float_t diff) { fDiffusionL = diff; };
-
- virtual Float_t GetRowPadSize() { return fRowPadSize; };
- virtual Float_t GetColPadSize() { return fColPadSize; };
- virtual Float_t GetTimeBinSize() { return fTimeBinSize; };
-
- virtual Float_t GetGasGain() { return fGasGain; };
- virtual Float_t GetNoise() { return fNoise; };
- virtual Float_t GetChipGain() { return fChipGain; };
- virtual Float_t GetADCoutRange() { return fADCoutRange; };
- virtual Float_t GetADCinRange() { return fADCinRange; };
- virtual Int_t GetADCthreshold() { return fADCthreshold; };
- virtual Float_t GetDiffusionT() { return fDiffusionT; };
- virtual Float_t GetDiffusionL() { return fDiffusionL; };
-
- virtual Int_t GetRowMax(Int_t iplan, Int_t icham) { return fRowMax[iplan-1][icham-1]; };
- virtual Int_t GetColMax(Int_t iplan) { return fColMax[iplan-1]; };
- virtual Int_t GetTimeMax() { return fTimeMax; };
-
-protected:
- Int_t fIdSens; // Sensitive volume identifier
-
- Int_t fIdSpace1; // Spaceframe volume identifier
- Int_t fIdSpace2; //
- Int_t fIdSpace3; //
-
- Int_t fIdChamber1; // Driftchamber volume identifier
- Int_t fIdChamber2; //
- Int_t fIdChamber3; //
-
- Int_t fSensSelect; // Switch to select only parts of the detector
- Int_t fSensPlane; // Sensitive detector plane
- Int_t fSensChamber; // Sensitive detector chamber
- Int_t fSensSector; // Sensitive detector sector
-
- Int_t fRowMax[kNplan][kNcham]; // Number of pad-rows
- Int_t fColMax[kNplan]; // Number of pad-columns
- Int_t fTimeMax; // Number of time buckets
-
- Float_t fRowPadSize; // Pad size in z-direction
- Float_t fColPadSize; // Pad size in rphi-direction
- Float_t fTimeBinSize; // Size of the time buckets
-
- Float_t fGasGain; // Gas gain
- Float_t fNoise; // Electronics noise
- Float_t fChipGain; // Electronics gain
- Float_t fADCoutRange; // ADC output range (number of channels)
- Float_t fADCinRange; // ADC input range (input charge)
- Int_t fADCthreshold; // ADC threshold in ADC channel
- Float_t fDiffusionT; // Diffusion in transverse direction
- Float_t fDiffusionL; // Diffusion in logitudinal direction
-
-private:
- virtual Double_t BetheBloch(Double_t bg);
- virtual void Diffusion(Float_t driftlength, Float_t *xyz);
- virtual Float_t PadResponse(Float_t x);
-
- TF1 *fDeltaE; // Energy distribution of the delta-electrons
-
- ClassDef(AliTRDv2,1) // Transition Radiation Detector version 2 (slow simulator)
-
-};
-
-#endif
{
gMC->Gsatt("*", "seen", -1);
gMC->Gsatt("alic", "seen", 0);
- gROOT->Macro("ViewTRD.C");
+ AliTRD *TRD = gAlice->GetModule("TRD");
+ if (TRD->Hole())
+ gROOT->Macro("ViewTRDhole.C");
+ else
+ gROOT->Macro("ViewTRDfull.C");
gMC->Gdopt("hide", "on");
gMC->Gdopt("shad", "on");
gMC->Gsatt("*", "fill", 7);
# C++ sources
-SRCS = AliTRD.cxx AliTRDv0.cxx AliTRDv1.cxx AliTRDv2.cxx \
- AliTRDmatrix.cxx AliTRDpixel.cxx
+SRCS = AliTRD.cxx AliTRDv0.cxx AliTRDv1.cxx \
+ AliTRDpixel.cxx AliTRDmatrix.cxx
# C++ Headers
#pragma link off all classes;
#pragma link off all functions;
-#pragma link C++ class AliTRD;
+#pragma link C++ class AliTRD-;
#pragma link C++ class AliTRDv0;
#pragma link C++ class AliTRDv1;
-#pragma link C++ class AliTRDv2;
#pragma link C++ class AliTRDhit;
#pragma link C++ class AliTRDdigit;
+#pragma link C++ class AliTRDcluster;
#pragma link C++ class AliTRDpixel;
#pragma link C++ class AliTRDmatrix;
--- /dev/null
+{
+
+ gMC->Gsatt("B071","SEEN", 0);
+ gMC->Gsatt("B074","SEEN", 0);
+ gMC->Gsatt("B075","SEEN", 0);
+ gMC->Gsatt("B077","SEEN", 0);
+ gMC->Gsatt("BTR1","SEEN", 0);
+ gMC->Gsatt("BTR2","SEEN", 0);
+ gMC->Gsatt("BTR3","SEEN", 0);
+ gMC->Gsatt("TRD1","SEEN", 0);
+ gMC->Gsatt("UCII","SEEN", 0);
+ gMC->Gsatt("UCIM","SEEN", 0);
+ gMC->Gsatt("UCIO","SEEN", 0);
+ gMC->Gsatt("UL02","SEEN", 1);
+ gMC->Gsatt("UL05","SEEN", 1);
+ gMC->Gsatt("UL06","SEEN", 1);
+
+}
--- /dev/null
+{
+
+ gMC->Gsatt("B071","SEEN", 0);
+ gMC->Gsatt("B074","SEEN", 0);
+ gMC->Gsatt("B075","SEEN", 0);
+ gMC->Gsatt("B077","SEEN", 0);
+ gMC->Gsatt("B078","SEEN", 0);
+ gMC->Gsatt("B079","SEEN", 0);
+ gMC->Gsatt("BTR1","SEEN", 0);
+ gMC->Gsatt("BTR2","SEEN", 0);
+ gMC->Gsatt("BTR3","SEEN", 0);
+ gMC->Gsatt("TRD1","SEEN", 0);
+ gMC->Gsatt("TRD2","SEEN", 0);
+ gMC->Gsatt("TRD3","SEEN", 0);
+ gMC->Gsatt("UCII","SEEN", 0);
+ gMC->Gsatt("UCIM","SEEN", 0);
+ gMC->Gsatt("UCIO","SEEN", 0);
+ gMC->Gsatt("UL02","SEEN", 1);
+ gMC->Gsatt("UL05","SEEN", 1);
+ gMC->Gsatt("UL06","SEEN", 1);
+
+}
--- /dev/null
+{
+
+/////////////////////////////////////////////////////////////////////////
+//
+// Example macro for the analysis of the TRD cluster
+//
+/////////////////////////////////////////////////////////////////////////
+
+ // Dynamically link some shared libs
+ if (gClassTable->GetID("AliRun") < 0) {
+ gROOT->LoadMacro("loadlibs.C");
+ loadlibs();
+ }
+
+ // Input file name
+ Char_t *alifile = "galice_c_v0.root";
+
+ // Event number
+ Int_t nEvent = 0;
+
+ // Define the objects
+ AliTRDv1 *TRD;
+ TClonesArray *TRDCluster;
+ AliTRDcluster *OneTRDcluster;
+
+ TH1F *hZ = new TH1F("hZ","Cluster z-position",700,-350.0,350.0);
+
+ // Connect the AliRoot file containing Geometry, Kine, Hits, and Digits
+ TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(alifile);
+ if (!gafl) {
+ cout << "Open the ALIROOT-file " << alifile << endl;
+ gafl = new TFile(alifile);
+ }
+ else {
+ cout << alifile << " is already open" << endl;
+ }
+
+ // Get AliRun object from file or create it if not on file
+ if (!gAlice) {
+ gAlice = (AliRun*) gafl->Get("gAlice");
+ if (gAlice)
+ cout << "AliRun object found on file" << endl;
+ else
+ gAlice = new AliRun("gAlice","Alice test program");
+ }
+
+ // Import the Trees for the event nEvent in the file
+ Int_t nparticles = gAlice->GetEvent(nEvent);
+ cout << "nparticles = " << nparticles << endl;
+ if (nparticles <= 0) break;
+
+ // Get the pointer to the tree
+ TTree *ClusterTree = gAlice->TreeD();
+
+ // Get the pointer to the detector classes
+ TRD = (AliTRDv1 *) gAlice->GetDetector("TRD");
+ // Get the pointer to the hit container
+ if (TRD) TRDCluster = TRD->Cluster();
+
+ // Reconstruct the address
+ ClusterTree->GetBranch("TRDcluster")->SetAddress(&TRDCluster);
+
+ Int_t nEntries = ClusterTree->GetEntries();
+ cout << "Number of entries in cluster tree = " << nEntries << endl;
+
+ // Loop through all entries in the tree
+ Int_t nbytes;
+ for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) {
+
+ cout << "iEntry = " << iEntry << endl;
+
+ // Import the tree
+ gAlice->ResetDigits();
+ nbytes += ClusterTree->GetEvent(iEntry);
+
+ // Get the number of digits in the detector
+ Int_t nTRDCluster = TRDCluster->GetEntriesFast();
+ cout << " nTRDCluster = " << nTRDCluster << endl;
+
+ // Loop through all TRD digits
+ for (Int_t iTRDCluster = 0; iTRDCluster < nTRDCluster; iTRDCluster++) {
+
+ // Get the information for this digit
+ OneTRDcluster = (AliTRDcluster*) TRDCluster->UncheckedAt(iTRDCluster);
+ hZ->Fill(OneTRDcluster->fZ);
+
+ }
+
+ }
+
+ hZ->Draw();
+
+}
--- /dev/null
+{
+
+/////////////////////////////////////////////////////////////////////////
+//
+// Creates cluster from the hit information (fast simulator).
+// An additional hit-tree is added to the input file.
+//
+/////////////////////////////////////////////////////////////////////////
+
+ // Dynamically link some shared libs
+ if (gClassTable->GetID("AliRun") < 0) {
+ gROOT->LoadMacro("loadlibs.C");
+ loadlibs();
+ }
+
+ // Input (and output) file name
+ Char_t *alifile = "galice_c_v0.root";
+
+ // Event number
+ Int_t nEvent = 0;
+
+ // Connect the AliRoot file containing Geometry, Kine, Hits, and Digits
+ TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(alifile);
+ if (!gafl) {
+ cout << "Open the ALIROOT-file " << alifile << endl;
+ gafl = new TFile(alifile,"UPDATE");
+ }
+ else {
+ cout << alifile << " is already open" << endl;
+ }
+
+ // Get AliRun object from file or create it if not on file
+ if (!gAlice) {
+ gAlice = (AliRun*) gafl->Get("gAlice");
+ if (gAlice)
+ cout << "AliRun object found on file" << endl;
+ else
+ gAlice = new AliRun("gAlice","Alice test program");
+ }
+
+ // Import the Trees for the event nEvent in the file
+ Int_t nparticles = gAlice->GetEvent(nEvent);
+ if (nparticles <= 0) break;
+
+ // Get the pointer to the detector classes
+ AliTRDv0 *TRD = (AliTRDv0*) gAlice->GetDetector("TRD");
+
+ // Create the clusters
+ TRD->Hits2Clusters();
+
+ // Write the new tree into the input file
+ cout << "Entries in digits tree = " << gAlice->TreeD()->GetEntries() << endl;
+ Char_t treeName[7];
+ sprintf(treeName,"TreeD%d",nEvent);
+ gAlice->TreeD()->Write(treeName);
+
+}
--- /dev/null
+{
+
+/////////////////////////////////////////////////////////////////////////
+//
+// Creates cluster from the digit information. An additional hit-tree
+// is added to the input file.
+//
+/////////////////////////////////////////////////////////////////////////
+
+ // Dynamically link some shared libs
+ if (gClassTable->GetID("AliRun") < 0) {
+ gROOT->LoadMacro("loadlibs.C");
+ loadlibs();
+ }
+
+ // Input (and output) file name
+ Char_t *alifile = "galice_c_v1.root";
+
+ // Event number
+ Int_t nEvent = 0;
+
+ // Connect the AliRoot file containing Geometry, Kine, Hits, and Digits
+ TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(alifile);
+ if (!gafl) {
+ cout << "Open the ALIROOT-file " << alifile << endl;
+ gafl = new TFile(alifile,"UPDATE");
+ }
+ else {
+ cout << alifile << " is already open" << endl;
+ }
+
+ // Get AliRun object from file or create it if not on file
+ if (!gAlice) {
+ gAlice = (AliRun*) gafl->Get("gAlice");
+ if (gAlice)
+ cout << "AliRun object found on file" << endl;
+ else
+ gAlice = new AliRun("gAlice","Alice test program");
+ }
+
+ // Import the Trees for the event nEvent in the file
+ Int_t nparticles = gAlice->GetEvent(nEvent);
+ if (nparticles <= 0) break;
+
+ // Get the pointer to the detector classes
+ AliTRDv1 *TRD = (AliTRDv1*) gAlice->GetDetector("TRD");
+
+ // Create the clusters
+ TRD->Digits2Clusters();
+
+ // Write the new tree into the input file
+ cout << "Entries in digits tree = " << gAlice->TreeD()->GetEntries() << endl;
+ Char_t treeName[7];
+ sprintf(treeName,"TreeD%d",nEvent);
+ gAlice->TreeD()->Write(treeName);
+
+}
--- /dev/null
+{
+
+/////////////////////////////////////////////////////////////////////////
+//
+// Example macro for the analysis of the TRD digits and the use
+// of the AliTRDmatrix class.
+//
+/////////////////////////////////////////////////////////////////////////
+
+ // Dynamically link some shared libs
+ if (gClassTable->GetID("AliRun") < 0) {
+ gROOT->LoadMacro("loadlibs.C");
+ loadlibs();
+ }
+
+ // Input file name
+//Char_t *alifile = "galice_d_v1.root";
+ Char_t *alifile = "galice_c_v1.root";
+
+ // Event number
+ Int_t nEvent = 0;
+
+ // Define the objects
+ AliTRDv1 *TRD;
+ TClonesArray *TRDDigits;
+ AliTRDdigit *OneTRDDigit;
+
+ // Connect the AliRoot file containing Geometry, Kine, Hits, and Digits
+ TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(alifile);
+ if (!gafl) {
+ cout << "Open the ALIROOT-file " << alifile << endl;
+ gafl = new TFile(alifile);
+ }
+ else {
+ cout << alifile << " is already open" << endl;
+ }
+
+ // Get AliRun object from file or create it if not on file
+ if (!gAlice) {
+ gAlice = (AliRun*) gafl->Get("gAlice");
+ if (gAlice)
+ cout << "AliRun object found on file" << endl;
+ else
+ gAlice = new AliRun("gAlice","Alice test program");
+ }
+
+ // Import the Trees for the event nEvent in the file
+ Int_t nparticles = gAlice->GetEvent(nEvent);
+ if (nparticles <= 0) break;
+
+ // Get the pointer to the hit-tree
+ TTree *DigitsTree = gAlice->TreeD();
+
+ // Get the pointer to the detector classes
+ TRD = (AliTRDv1*) gAlice->GetDetector("TRD");
+ // Get the pointer to the hit container
+ if (TRD) TRDDigits = TRD->Digits();
+
+ DigitsTree->GetBranch("TRD")->SetAddress(&TRDDigits);
+
+ // Define the detector matrix for one chamber (Sector 6, Chamber 3, Plane 1)
+ const Int_t iSec = 6;
+ const Int_t iCha = 3;
+ const Int_t iPla = 1;
+ Int_t rowMax = TRD->GetRowMax(iPla,iCha,iSec);
+ Int_t colMax = TRD->GetColMax(iPla);
+ Int_t timeMax = TRD->GetTimeMax();
+ cout << " rowMax = " << rowMax
+ << " colMax = " << colMax
+ << " timeMax = " << timeMax << endl;
+ AliTRDmatrix *TRDMatrix = new AliTRDmatrix(rowMax,colMax,timeMax,iSec,iCha,iPla);
+
+ Int_t nEntries = DigitsTree->GetEntries();
+ cout << "Number of entries in digits tree = " << nEntries << endl;
+
+ // Loop through all entries in the tree
+ Int_t nbytes;
+ for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) {
+
+ cout << "iEntry = " << iEntry << endl;
+
+ // Import the tree
+ gAlice->ResetDigits();
+ nbytes += DigitsTree->GetEvent(iEntry);
+
+ // Get the number of digits in the detector
+ Int_t nTRDDigits = TRDDigits->GetEntriesFast();
+ cout << " nTRDDigits = " << nTRDDigits << endl;
+
+ // Loop through all TRD digits
+ for (Int_t iTRDDigits = 0; iTRDDigits < nTRDDigits; iTRDDigits++) {
+
+ // Get the information for this digit
+ OneTRDDigit = (AliTRDdigit*) TRDDigits->UncheckedAt(iTRDDigits);
+ Int_t signal = OneTRDDigit->fSignal;
+ Int_t sector = OneTRDDigit->fSector;
+ Int_t chamber = OneTRDDigit->fChamber;
+ Int_t plane = OneTRDDigit->fPlane;
+ Int_t row = OneTRDDigit->fRow;
+ Int_t col = OneTRDDigit->fCol;
+ Int_t time = OneTRDDigit->fTime;
+
+ // Fill the detector matrix
+ if (signal > 1) {
+ TRDMatrix->SetSignal(row,col,time,signal);
+ }
+
+ }
+
+ }
+
+ // Display the detector matrix
+ TRDMatrix->Draw();
+ TRDMatrix->DrawRow(18);
+ TRDMatrix->DrawCol(58);
+ TRDMatrix->DrawTime(20);
+
+}
--- /dev/null
+{
+
+/////////////////////////////////////////////////////////////////////////
+//
+// Creates the digits from the hit information. An additional hit-tree
+// is added to the input file.
+//
+/////////////////////////////////////////////////////////////////////////
+
+ // Dynamically link some shared libs
+ if (gClassTable->GetID("AliRun") < 0) {
+ gROOT->LoadMacro("loadlibs.C");
+ loadlibs();
+ }
+
+ // Input (and output) file name
+ Char_t *alifile = "galice_d_v1.root";
+
+ // Event number
+ Int_t nEvent = 0;
+
+ // Connect the AliRoot file containing Geometry, Kine, and Hits
+ TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(alifile);
+ if (!gafl) {
+ cout << "Open the ALIROOT-file " << alifile << endl;
+ gafl = new TFile(alifile,"UPDATE");
+ }
+ else {
+ cout << alifile << " is already open" << endl;
+ }
+
+ // Get AliRun object from file or create it if not on file
+ if (!gAlice) {
+ gAlice = (AliRun*) gafl->Get("gAlice");
+ if (gAlice)
+ cout << "AliRun object found on file" << endl;
+ else
+ gAlice = new AliRun("gAlice","Alice test program");
+ }
+
+ // Import the Trees for the event nEvent in the file
+ Int_t nparticles = gAlice->GetEvent(nEvent);
+ if (nparticles <= 0) break;
+
+ // Get the pointer to the detector class
+ AliTRDv1 *TRD = (AliTRDv1*) gAlice->GetDetector("TRD");
+
+ // Create the digits and fill the digits-tree
+ TRD->Hits2Digits();
+
+ // Write the new tree into the input file
+ cout << "Entries in digits tree = " << gAlice->TreeD()->GetEntries() << endl;
+ Char_t treeName[7];
+ sprintf(treeName,"TreeD%d",nEvent);
+ gAlice->TreeD()->Write(treeName);
+
+}
AliTRD *TRD = new AliTRDv1("TRD","TRD version 0");
// Select the gas mixture (0: 97% Xe + 3% isobutane, 1: 90% Xe + 10% CO2)
TRD->SetGasMix(0);
-TRD->SetHits(1);
}
if(iFMD) {
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff