/*
$Log$
+Revision 1.37 2001/06/12 07:17:18 kowal2
+Hits2SDigits method implemented (summable digits)
+
+Revision 1.36 2001/05/16 14:57:25 alibrary
+New files for folders and Stack
+
+Revision 1.35 2001/05/08 16:02:22 kowal2
+Updated material specifications
+
+Revision 1.34 2001/05/08 15:00:15 hristov
+Corrections for tracking in arbitrary magnenetic field. Changes towards a concept of global Alice track. Back propagation of reconstructed tracks (Yu.Belikov)
+
+Revision 1.33 2001/04/03 12:40:43 kowal2
+Removed printouts
+
+Revision 1.32 2001/03/12 17:47:36 hristov
+Changes needed on Sun with CC 5.0
+
+Revision 1.31 2001/03/12 08:21:50 kowal2
+Corrected C++ bug in the material definitions
+
+Revision 1.30 2001/03/01 17:34:47 kowal2
+Correction due to the accuracy problem
+
+Revision 1.29 2001/02/28 16:34:40 kowal2
+Protection against nonphysical values of the avalanche size,
+10**6 is the maximum
+
+Revision 1.28 2001/01/26 19:57:19 hristov
+Major upgrade of AliRoot code
+
+Revision 1.27 2001/01/13 17:29:33 kowal2
+Sun compiler correction
+
+Revision 1.26 2001/01/10 07:59:43 kowal2
+Corrections to load points from the noncompressed hits.
+
Revision 1.25 2000/11/02 07:25:31 kowal2
Changes due to the new hit structure.
Memory leak removed.
#include <TObjectTable.h>
#include "TParticle.h"
#include "AliTPC.h"
-#include <TFile.h>
+#include <TFile.h>
+#include <TROOT.h>
+#include <TSystem.h>
#include "AliRun.h"
#include <iostream.h>
+#include <stdlib.h>
#include <fstream.h>
#include "AliMC.h"
#include "AliMagF.h"
//MI changes
fDigitsArray = 0;
fClustersArray = 0;
- fTPCParam=0;
+ fDefaults = 0;
fTrackHits = 0;
fHitType = 2;
fTPCParam = 0;
//MI change
fDigitsArray = 0;
fClustersArray= 0;
+ fDefaults = 0;
//
fTrackHits = new AliTPCTrackHits; //MI - 13.09.2000
fTrackHits->SetHitPrecision(0.002);
// Set TPC parameters
//
- if (!strcmp(title,"Default")) {
- fTPCParam = new AliTPCParamSR;
+
+ if (!strcmp(title,"Default")) {
+ fTPCParam = new AliTPCParamSR;
} else {
cerr<<"AliTPC warning: in Config.C you must set non-default parameters\n";
fTPCParam=0;
//-----------------------------------------------------------------
// This is a track finder.
//-----------------------------------------------------------------
- AliTPCtracker::Clusters2Tracks(fTPCParam,of);
+ AliTPCtracker tracker(fTPCParam);
+ tracker.Clusters2Tracks(gFile,of);
}
//_____________________________________________________________________________
AliMixture(37, "Mylar",amat,zmat,density,-3,wmat);
- // G10 60% SiO2 + 40% epoxy, I use A and Z for SiO2
+ // SiO2 - used later for the glass fiber
amat[0]=28.086;
amat[1]=15.9994;
wmat[0]=1.;
wmat[1]=2.;
- density = 1.7;
- AliMixture(38,"SiO2",amat,zmat,2.2,-2,wmat); //SiO2 - quartz
-
- gMC->Gfmate((*fIdmate)[38],namate,amat[0],zmat[0],rho,X0,absl,buf,nbuf);
-
- AliMaterial(39,"G10",amat[0],zmat[0],density,999.,999.);
+ AliMixture(38,"SiO2",amat,zmat,2.2,-2,wmat); //SiO2 - quartz (rho=2.2)
// Al
// Si
amat[0] = 28.086;
- zmat[0] = 14.,
+ zmat[0] = 14.;
density = 2.33;
AliMixture(44,"Plexiglas",amat,zmat,density,-3,wmat);
+ // Epoxy - C14 H20 O3
+
+
+ amat[0]=12.011;
+ amat[1]=1.;
+ amat[2]=15.9994;
+
+ zmat[0]=6.;
+ zmat[1]=1.;
+ zmat[2]=8.;
+
+ wmat[0]=14.;
+ wmat[1]=20.;
+ wmat[2]=3.;
+
+ density=1.25;
+
+ AliMixture(45,"Epoxy",amat,zmat,density,-3,wmat);
+
+ // Carbon
+
+ amat[0]=12.011;
+ zmat[0]=6.;
+ density= 2.265;
+
+ AliMaterial(46,"C",amat[0],zmat[0],density,999.,999.);
+
+ // get epoxy
+
+ gMC->Gfmate((*fIdmate)[45],namate,amat[1],zmat[1],rho,X0,absl,buf,nbuf);
+
+ // Carbon fiber
+
+ wmat[0]=0.644; // by weight!
+ wmat[1]=0.356;
+
+ density=0.5*(1.25+2.265);
+
+ AliMixture(47,"Cfiber",amat,zmat,density,2,wmat);
+
+ // get SiO2
+
+ gMC->Gfmate((*fIdmate)[38],namate,amat[0],zmat[0],rho,X0,absl,buf,nbuf);
+
+ wmat[0]=0.725; // by weight!
+ wmat[1]=0.275;
+
+ density=1.7;
+
+ AliMixture(39,"G10",amat,zmat,density,2,wmat);
+
+
//----------------------------------------------------------
AliMedium(11,"Si",41,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
AliMedium(12,"G10",39,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
AliMedium(13,"Plexiglas",44,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
+ AliMedium(14,"Epoxy",45,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);
+ AliMedium(15,"Cfiber",47,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
}
-void AliTPC::Digits2Clusters(TFile *of)
+void AliTPC::Digits2Clusters(TFile *of, Int_t eventnumber)
{
//-----------------------------------------------------------------
// This is a simple cluster finder.
//-----------------------------------------------------------------
- AliTPCclusterer::Digits2Clusters(fTPCParam,of);
+ AliTPCclusterer::Digits2Clusters(fTPCParam,of,eventnumber);
}
extern Double_t SigmaY2(Double_t, Double_t, Double_t);
extern Double_t SigmaZ2(Double_t, Double_t);
//_____________________________________________________________________________
-void AliTPC::Hits2Clusters(TFile *of)
+void AliTPC::Hits2Clusters(TFile *of, Int_t eventn)
{
//--------------------------------------------------------
// TPC simple cluster generator from hits
return;
}
- if(fTPCParam == 0){
- printf("AliTPCParam MUST be created firstly\n");
- return;
- }
+ //if(fDefaults == 0) SetDefaults();
Float_t sigmaRphi,sigmaZ,clRphi,clZ;
//
TParticle *particle; // pointer to a given particle
AliTPChit *tpcHit; // pointer to a sigle TPC hit
- TClonesArray *particles; //pointer to the particle list
Int_t sector;
Int_t ipart;
Float_t xyz[5];
TTree *tH = gAlice->TreeH();
Stat_t ntracks = tH->GetEntries();
- particles=gAlice->Particles();
//Switch to the output file
of->cd();
+ char cname[100];
+
+ sprintf(cname,"TreeC_TPC_%d",eventn);
+
fTPCParam->Write(fTPCParam->GetTitle());
AliTPCClustersArray carray;
carray.Setup(fTPCParam);
continue;
}
ipart=tpcHit->Track();
- particle=(TParticle*)particles->UncheckedAt(ipart);
+ particle=gAlice->Particle(ipart);
pl=particle->Pz();
pt=particle->Pt();
if(pt < 1.e-9) pt=1.e-9;
cerr<<"Number of made clusters : "<<nclusters<<" \n";
- carray.GetTree()->Write();
+ carray.GetTree()->Write(cname);
savedir->cd(); //switch back to the input file
//
TParticle *particle; // pointer to a given particle
AliTPChit *tpcHit; // pointer to a sigle TPC hit
- TClonesArray *particles; //pointer to the particle list
Int_t sector,nhits;
Int_t ipart;
const Int_t kcmaxhits=30000;
TTree *tH = gAlice->TreeH();
Stat_t ntracks = tH->GetEntries();
- particles=gAlice->Particles();
- Int_t npart = particles->GetEntriesFast();
+ Int_t npart = gAlice->GetNtrack();
//------------------------------------------------------------
// Loop over tracks
sector=tpcHit->fSector; // sector number
if(sector != isec) continue;
ipart=tpcHit->Track();
- if (ipart<npart) particle=(TParticle*)particles->UncheckedAt(ipart);
+ if (ipart<npart) particle=gAlice->Particle(ipart);
//find row number
xxxx->Delete();
}
+//___________________________________________
+void AliTPC::SDigits2Digits(Int_t eventnumber)
+{
+
+
+ cerr<<"Digitizing TPC...\n";
+
+ Hits2Digits(eventnumber);
+
+
+ //write results
+
+ // char treeName[100];
+
+ // sprintf(treeName,"TreeD_%s_%d",fTPCParam->GetTitle(),eventnumber);
+
+ // GetDigitsArray()->GetTree()->Write(treeName);
+}
+//__________________________________________________________________
+void AliTPC::SetDefaults(){
+
+
+ cerr<<"Setting default parameters...\n";
+
+ // Set response functions
+
+ AliTPCParamSR *param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60");
+ AliTPCPRF2D * prfinner = new AliTPCPRF2D;
+ AliTPCPRF2D * prfouter = new AliTPCPRF2D;
+ AliTPCRF1D * rf = new AliTPCRF1D(kTRUE);
+ rf->SetGauss(param->GetZSigma(),param->GetZWidth(),1.);
+ rf->SetOffset(3*param->GetZSigma());
+ rf->Update();
+
+ TDirectory *savedir=gDirectory;
+ TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2d.root");
+ if (!f->IsOpen()) {
+ cerr<<"Can't open $ALICE_ROOT/TPC/AliTPCprf2d.root !\n" ;
+ exit(3);
+ }
+ prfinner->Read("prf_07504_Gati_056068_d02");
+ prfouter->Read("prf_10006_Gati_047051_d03");
+ f->Close();
+ savedir->cd();
+
+ param->SetInnerPRF(prfinner);
+ param->SetOuterPRF(prfouter);
+ param->SetTimeRF(rf);
+ // set fTPCParam
+
+ SetParam(param);
+
+
+ fDefaults = 1;
+
+}
//__________________________________________________________________
-void AliTPC::Hits2Digits()
+void AliTPC::Hits2Digits(Int_t eventnumber)
{
//----------------------------------------------------
- // Loop over all sectors
+ // Loop over all sectors for a single event
//----------------------------------------------------
- if(fTPCParam == 0){
- printf("AliTPCParam MUST be created firstly\n");
- return;
- }
+
+ if(fDefaults == 0) SetDefaults(); // check if the parameters are set
+
+ //setup TPCDigitsArray
+
+ if(GetDigitsArray()) delete GetDigitsArray();
+
+ AliTPCDigitsArray *arr = new AliTPCDigitsArray;
+ arr->SetClass("AliSimDigits");
+ arr->Setup(fTPCParam);
+ arr->MakeTree(fDigitsFile);
+ SetDigitsArray(arr);
+
+ fDigitsSwitch=0; // standard digits
+
+ cerr<<"Digitizing TPC...\n";
+
+ for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) Hits2DigitsSector(isec);
+
+ // write results
+
+ char treeName[100];
+
+ sprintf(treeName,"TreeD_%s_%d",fTPCParam->GetTitle(),eventnumber);
+
+ GetDigitsArray()->GetTree()->Write(treeName);
+
+
+}
+
+
+
+//__________________________________________________________________
+void AliTPC::Hits2SDigits(Int_t eventnumber)
+{
+
+ //-----------------------------------------------------------
+ // summable digits - 16 bit "ADC", no noise, no saturation
+ //-----------------------------------------------------------
+
+ //----------------------------------------------------
+ // Loop over all sectors for a single event
+ //----------------------------------------------------
+
+
+ if(fDefaults == 0) SetDefaults();
+
+ //setup TPCDigitsArray
+
+ if(GetDigitsArray()) delete GetDigitsArray();
+
+ AliTPCDigitsArray *arr = new AliTPCDigitsArray;
+ arr->SetClass("AliSimDigits");
+ arr->Setup(fTPCParam);
+ arr->MakeTree(fDigitsFile);
+ SetDigitsArray(arr);
+
+ cerr<<"Digitizing TPC...\n";
+
+ fDigitsSwitch=1; // summable digits
for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) Hits2DigitsSector(isec);
+
+ // write results
+
+ char treeName[100];
+
+ sprintf(treeName,"TreeS_%s_%d",fTPCParam->GetTitle(),eventnumber);
+
+ GetDigitsArray()->GetTree()->Write(treeName);
+
}
+
//_____________________________________________________________________________
void AliTPC::Hits2DigitsSector(Int_t isec)
{
// Get the access to the track hits
//-------------------------------------------------------
+ // check if the parameters are set - important if one calls this method
+ // directly, not from the Hits2Digits
+
+ if(fDefaults == 0) SetDefaults();
TTree *tH = gAlice->TreeH(); // pointer to the hits tree
Stat_t ntracks = tH->GetEntries();
TObjArray **row;
- printf("*** Processing sector number %d ***\n",isec);
+ printf("*** Processing sector number %d ***\n",isec);
Int_t nrows =fTPCParam->GetNRow(isec);
Int_t i;
- if (fDigitsArray->GetTree()==0) fDigitsArray->MakeTree();
+ if (fDigitsArray->GetTree()==0) fDigitsArray->MakeTree(fDigitsFile);
for (i=0;i<nrows;i++){
fDigitsArray->StoreRow(isec,i);
- Int_t ndig = dig->GetDigitSize();
+ Int_t ndig = dig->GetDigitSize();
+
- printf("*** Sector, row, compressed digits %d %d %d ***\n",isec,i,ndig);
+ printf("*** Sector, row, compressed digits %d %d %d ***\n",isec,i,ndig);
fDigitsArray->ClearRow(isec,i);
Int_t gi =it*nofPads+ip; // global index
- q = gRandom->Gaus(q,fTPCParam->GetNoise()*fTPCParam->GetNoiseNormFac());
+ if(fDigitsSwitch == 0){
- q = (Int_t)q;
+ q = gRandom->Gaus(q,fTPCParam->GetNoise()*fTPCParam->GetNoiseNormFac());
- if(q <=zerosup) continue; // do not fill zeros
- if(q > fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat(); // saturation
+ q = (Int_t)q;
+
+ if(q <=zerosup) continue; // do not fill zeros
+ if(q > fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat(); // saturation
+
+ }
+
+ else {
+ q *= 16.;
+ q = (Int_t)q;
+ }
//
// "real" signal or electronic noise (list = -1)?
if((gRandom->Rndm(0)) < attProb) continue; // electron lost!
xyz[0]=tpcHit->X();
xyz[1]=tpcHit->Y();
- xyz[2]=tpcHit->Z();
- xyz[3]= (Float_t) (-gasgain*TMath::Log(gRandom->Rndm()));
+ xyz[2]=tpcHit->Z();
+ //
+ // protection for the nonphysical avalanche size (10**6 maximum)
+ //
+ Double_t rn=TMath::Max(gRandom->Rndm(0),1.93e-22);
+ xyz[3]= (Float_t) (-gasgain*TMath::Log(rn));
index[0]=1;
TransportElectron(xyz,index); //MI change -august
//
Int_t i;
//
- printf("\n");
- for(i=0;i<35;i++) printf("*");
- printf(" TPC_INIT ");
- for(i=0;i<35;i++) printf("*");
- printf("\n");
- //
- for(i=0;i<80;i++) printf("*");
- printf("\n");
+ if(fDebug) {
+ printf("\n%s: ",ClassName());
+ for(i=0;i<35;i++) printf("*");
+ printf(" TPC_INIT ");
+ for(i=0;i<35;i++) printf("*");
+ printf("\n%s: ",ClassName());
+ //
+ for(i=0;i<80;i++) printf("*");
+ printf("\n");
+ }
}
//_____________________________________________________________________________
-void AliTPC::MakeBranch(Option_t* option)
+void AliTPC::MakeBranch(Option_t* option, const char *file)
{
//
// Create Tree branches for the TPC.
char branchname[10];
sprintf(branchname,"%s",GetName());
- AliDetector::MakeBranch(option);
+ AliDetector::MakeBranch(option,file);
- char *d = strstr(option,"D");
+ const char *d = strstr(option,"D");
if (fDigits && gAlice->TreeD() && d) {
- gAlice->TreeD()->Branch(branchname,&fDigits, buffersize);
- printf("Making Branch %s for digits\n",branchname);
+ MakeBranchInTree(gAlice->TreeD(),
+ branchname, &fDigits, buffersize, file);
}
- if (fHitType&2) MakeBranch2(option); // MI change 14.09.2000
+ if (fHitType&2) MakeBranch2(option,file); // MI change 14.09.2000
}
//_____________________________________________________________________________
//add nonisochronity (not implemented yet)
}
-//_____________________________________________________________________________
-void AliTPC::Streamer(TBuffer &R__b)
-{
- //
- // Stream an object of class AliTPC.
- //
- if (R__b.IsReading()) {
- Version_t R__v = R__b.ReadVersion(); if (R__v) { }
- AliDetector::Streamer(R__b);
- R__b >> fTPCParam;
- if (R__v < 2) return;
- R__b >> fNsectors;
- R__b >> fHitType;
-
- } else {
- R__b.WriteVersion(AliTPC::IsA());
- AliDetector::Streamer(R__b);
- R__b << fTPCParam;
- R__b << fNsectors;
- R__b << fHitType;
- }
-}
ClassImp(AliTPCdigit)
//________________________________________________________________________
// Additional code because of the AliTPCTrackHits
-void AliTPC::MakeBranch2(Option_t *option)
+void AliTPC::MakeBranch2(Option_t *option,const char *file)
{
//
// Create a new branch in the current Root Tree
sprintf(branchname,"%s2",GetName());
//
// Get the pointer to the header
- char *cH = strstr(option,"H");
+ const char *cH = strstr(option,"H");
//
if (fTrackHits && gAlice->TreeH() && cH) {
- // gAlice->TreeH()->Branch(branchname,&fTrackHits, fBufferSize);
AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHits",&fTrackHits,
gAlice->TreeH(),fBufferSize,1);
gAlice->TreeH()->GetListOfBranches()->Add(branch);
- printf("* AliDetector::MakeBranch * Making Branch %s for trackhits\n",branchname);
+ if (GetDebug()>1)
+ printf("* AliDetector::MakeBranch * Making Branch %s for trackhits\n",branchname);
+ const char folder [] = "RunMC/Event/Data";
+ if (GetDebug())
+ printf("%15s: Publishing %s to %s\n",ClassName(),branchname,folder);
+ Publish(folder,&fTrackHits,branchname);
+ if (file) {
+ TBranch *b = gAlice->TreeH()->GetBranch(branchname);
+ TDirectory *wd = gDirectory;
+ b->SetFile(file);
+ TIter next( b->GetListOfBranches());
+ while ((b=(TBranch*)next())) {
+ b->SetFile(file);
+ }
+ wd->cd();
+ if (GetDebug()>1)
+ cout << "Diverting branch " << branchname << " to file " << file << endl;
+ }
}
}
{
//
// add hit to the list
- TClonesArray &particles = *(gAlice->Particles());
Int_t rtrack;
if (fIshunt) {
int primary = gAlice->GetPrimary(track);
- ((TParticle *)particles[primary])->SetBit(kKeepBit);
+ gAlice->Particle(primary)->SetBit(kKeepBit);
rtrack=primary;
} else {
rtrack=track;
{
//
Int_t a = 0;
- if(fHitType==1) return AliDetector::LoadPoints(a);
+ /* if(fHitType==1) return AliDetector::LoadPoints(a);
LoadPoints2(a);
+ */
+ if(fHitType==1) AliDetector::LoadPoints(a);
+ else LoadPoints2(a);
+
// LoadPoints3(a);
}
} // end of loop over hits
xxxx->Delete();
-
-
}
+//_______________________________________________________________________________
+void AliTPC::Digits2Reco(Int_t eventnumber)
+{
+ // empty for a time being
+}