* provided "as is" without express or implied warranty. *
**************************************************************************/
-/*
-$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
-
-*/
+/* $Id$ */
///////////////////////////////////////////////////////////////////////////////
// //
// //
///////////////////////////////////////////////////////////////////////////////
+#include <stdlib.h>
+
+#include <TLorentzVector.h>
#include <TMath.h>
#include <TRandom.h>
-#include <TVector.h>
+#include <TVector.h>
+#include <TVirtualMC.h>
-#include "AliTRDv0.h"
+#include "AliConst.h"
#include "AliRun.h"
+#include "AliTRDgeometry.h"
+#include "AliTRDhit.h"
+#include "AliTRDv0.h"
#include "AliMC.h"
-#include "AliConst.h"
-
+
ClassImp(AliTRDv0)
+
+//_____________________________________________________________________________
+AliTRDv0::AliTRDv0():AliTRD()
+{
+ //
+ // AliTRDv0 default constructor
+ //
+
+ fHitsOn = 0;
+
+}
//_____________________________________________________________________________
AliTRDv0::AliTRDv0(const char *name, const char *title)
// Standard constructor for Transition Radiation Detector version 0
//
- fIdSens = 0;
fHitsOn = 0;
- fIdChamber1 = 0;
- fIdChamber2 = 0;
- fIdChamber3 = 0;
+}
- fRphiSigma = 0;
- fRphiDist = 0;
+//_____________________________________________________________________________
+AliTRDv0::~AliTRDv0()
+{
+ //
+ // AliTRDv0 destructor
+ //
}
-
+
//_____________________________________________________________________________
void AliTRDv0::CreateGeometry()
{
//
// Check that FRAME is there otherwise we have no place where to put the TRD
- AliModule* FRAME = gAlice->GetModule("FRAME");
- if (!FRAME) return;
+ AliModule* frame = gAlice->GetModule("FRAME");
+ if (!frame) return;
// Define the chambers
AliTRD::CreateGeometry();
}
-//_____________________________________________________________________________
-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()
{
//
- // Initialise Transition Radiation Detector after geometry is built
+ // Initialize Transition Radiation Detector after geometry is built
//
AliTRD::Init();
- // Identifier of the sensitive volume (amplification region)
- fIdSens = gMC->VolId("UL06");
-
- // Identifier of the TRD-driftchambers
- fIdChamber1 = gMC->VolId("UCIO");
- fIdChamber2 = gMC->VolId("UCIM");
- fIdChamber3 = gMC->VolId("UCII");
-
- // Parameter for Hits2Cluster
-
- // Position resolution in rphi-direction
- fRphiSigma = 0.02;
- // Minimum distance of non-overlapping cluster
- fRphiDist = 1.0;
-
- printf(" Fast simulator\n");
+ printf(" Fast simulator\n\n");
for (Int_t i = 0; i < 80; i++) printf("*");
printf("\n");
// crosses the border between amplification region and pad plane.
//
- Int_t vol[3];
- Int_t iIdSens, icSens;
- Int_t iIdChamber, icChamber;
+ Int_t pla = 0;
+ Int_t cha = 0;
+ Int_t sec = 0;
- Float_t hits[4];
+ Float_t hits[3];
+ Int_t det;
TLorentzVector p;
- TClonesArray &lhits = *fHits;
+
+ // Use pad plane as sensitive volume
+ TString cIdSens = "L";
+ TString cIdCurrent;
+ Char_t cIdChamber[3];
+ cIdChamber[2] = 0;
+
+ const Int_t kNplan = AliTRDgeometry::Nplan();
// Writing out hits enabled?
if (!(fHitsOn)) return;
// Use only charged tracks and count them only once per volume
if (gMC->TrackCharge() &&
- gMC->IsTrackExiting()) {
+ gMC->IsTrackEntering()) {
// Check on sensitive volume
- iIdSens = gMC->CurrentVolID(icSens);
- if (iIdSens == fIdSens) {
+ cIdCurrent = gMC->CurrentVolName();
+ if (cIdSens == cIdCurrent[1]) {
gMC->TrackPosition(p);
for (Int_t i = 0; i < 3; i++) hits[i] = p[i];
- // No charge created
- hits[3] = 0;
-
- // The sector number
- 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
- // 2: middle 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)
- 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;
-
- new(lhits[fNhits++]) AliTRDhit(fIshunt,gAlice->CurrentTrack(),vol,hits);
+
+ // The sector number (0 - 17)
+ // The numbering goes clockwise and starts at y = 0
+ Float_t phi = kRaddeg*TMath::ATan2(hits[0],hits[1]);
+ if (phi < 90.)
+ phi = phi + 270.;
+ else
+ phi = phi - 90.;
+ sec = ((Int_t) (phi / 20));
+
+ // The plane and chamber number
+ cIdChamber[0] = cIdCurrent[2];
+ cIdChamber[1] = cIdCurrent[3];
+ Int_t idChamber = atoi(cIdChamber);
+ cha = ((Int_t) idChamber / kNplan);
+ pla = ((Int_t) idChamber % kNplan);
+ det = fGeometry->GetDetector(pla,cha,sec);
+
+ AddHit(gAlice->GetMCApp()->GetCurrentTrackNumber(),det,hits,0,kTRUE);
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff