/**************************************************************************
* 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.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 -- fast simulator //
// //
//Begin_Html
/*
*/
//End_Html
// //
// //
///////////////////////////////////////////////////////////////////////////////
#include
#include
#include
#include "AliTRDv0.h"
#include "AliRun.h"
#include "AliMC.h"
#include "AliConst.h"
ClassImp(AliTRDv0)
//_____________________________________________________________________________
AliTRDv0::AliTRDv0(const char *name, const char *title)
:AliTRD(name, title)
{
//
// Standard constructor for Transition Radiation Detector version 0
//
fIdSens = 0;
fHitsOn = 0;
fIdChamber1 = 0;
fIdChamber2 = 0;
fIdChamber3 = 0;
fRphiSigma = 0;
fRphiDist = 0;
}
//_____________________________________________________________________________
void AliTRDv0::CreateGeometry()
{
//
// Create the GEANT geometry for the Transition Radiation Detector - Version 0
// This version covers the full azimuth.
//
// 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();
}
//_____________________________________________________________________________
void AliTRDv0::CreateMaterials()
{
//
// Create materials for the Transition Radiation Detector
//
AliTRD::CreateMaterials();
}
//_____________________________________________________________________________
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
//
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");
for (Int_t i = 0; i < 80; i++) printf("*");
printf("\n");
}
//_____________________________________________________________________________
void AliTRDv0::StepManager()
{
//
// 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.
//
Int_t vol[3];
Int_t iIdSens, icSens;
Int_t iIdChamber, icChamber;
Float_t hits[4];
TLorentzVector p;
TClonesArray &lhits = *fHits;
// Writing out hits enabled?
if (!(fHitsOn)) return;
// Use only charged tracks and count them only once per volume
if (gMC->TrackCharge() &&
gMC->IsTrackExiting()) {
// Check on 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;
// 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);
}
}
}