1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
17 Revision 1.5 2002/10/23 07:24:57 alibrary
18 Introducing Riostream.h
20 Revision 1.4 2001/01/26 20:00:53 hristov
21 Major upgrade of AliRoot code
23 Revision 1.3 2000/10/11 09:19:12 egangler
24 Corrected some bugs - it should compile now
26 Revision 1.2 2000/06/15 07:58:49 morsch
27 Code from MUON-dev joined
29 Revision 1.1.2.1 2000/04/18 09:11:15 morsch
30 Implementation of MUON Chamber Prototype Class
31 Both read digits from raw data or use the Monte-Carlo.
32 Rachid GUERNANE, IPN Lyon guernane@ipnl.in2p3.fr
37 Implementation of MUON Chamber Prototype Class
38 Both read digits from raw data or use the Monte-Carlo.
39 1-Feb-2000 Rachid GUERNANE, IPN Lyon guernane@ipnl.in2p3.fr
43 ////////////////////////////////////////////////
44 // Manager and hits classes for set:PROTO //
45 ////////////////////////////////////////////////
49 #include <TRotMatrix.h>
55 #include <TObjArray.h>
57 #include <TParticle.h>
63 #include <TDirectory.h>
64 #include <TObjectTable.h>
67 #include <AliDetector.h>
69 #include "AliMUONChamber.h"
70 #include "AliMUONproto.h"
71 #include "AliMUONHit.h"
73 #include "AliMUONClusterFinder.h"
75 #include "Riostream.h"
76 #include "AliCallf77.h"
78 //#include "chainalice2.h"
79 #include "AliMUONSegmentationV0.h"
80 //#include "AliMUONSegResV11.h"
82 ClassImp(AliMUONproto)
84 //___________________________________________
85 AliMUONproto::AliMUONproto()
88 cout << "\n Calling AliMUONproto constructor..." << endl;
95 //___________________________________________
96 AliMUONproto::AliMUONproto(const char *name, const char *title)
98 // AliMUON defines everything, but the chamber for NCH=1
101 // z-Positions of Chambers
102 const Float_t zch[1] = {975.};
105 const Float_t dmi[1] = {0.};
108 const Float_t dma[1] = {29.2};
111 // Default Parameters for ALICE2 prototype
114 (*fChambers)[0] = new AliMUONChamber();
115 AliMUONChamber* chamber = (AliMUONChamber*) (*fChambers)[0];
117 chamber->SetZ(zch[0]);
119 chamber->InitGeo(zch[0]);
120 chamber->SetRInner(dmi[0]/2);
121 chamber->SetROuter(dma[0]/2);
123 for (int i = 0; i <= 99; i++) {
129 #ifdef WE_FORGET_FOR_THE_MOMENT
130 //___________________________________________
131 void AliMUONproto::GetRawDigits(Int_t evnb, Int_t *lptr, Int_t ilen) {
145 AliMUON* MUON = (AliMUON*)gAlice->GetModule("MUON");
146 AliMUONSegmentationV0* seg = (AliMUONSegmentationV0*) Chamber(0).SegmentationModel(1);
151 for (Int_t i = 0; i < 10; i++) {
162 equiplen = lptr[ip++];
164 if (equiplen < 0 ) return;
166 if (itype != (int)(0XCACCA008)) {
168 if (ip < ilen-2) goto nwtype;
173 if ((serial == 190) && (equip == 1)) {
174 for (loop = 0; loop < equiplen; loop++) {
175 nochnl = (lptr[ip] & 0x7ff000 ) >> 12;
176 val = lptr[ip] & 0x3ff;
177 // fill digits from raw data according to cathode connexions
178 if (group[nochnl][2][1]!=0)
179 digits[0] = group[nochnl][2][1] - 12;
180 else if (group[nochnl][1][1]!=0)
181 digits[0] = group[nochnl][1][1] - 12;
183 digits[0] = group[nochnl][0][1] - 12;
184 digits[1] = group[nochnl][0][0];
185 if (digits[0] != seg->Ix(digits[0], digits[1]))
186 printf("Pb pour ix=%d,iy=%d\n", digits[0], digits[1]);
190 if (digits[2] >= fThreshold[nochnl])
191 MUON->AddDigits(ich, tracks, charges, digits);
199 if (ip < ilen-2) goto nwtype;
202 gAlice->TreeD()->Fill();
205 gAlice->TreeD()->Fill();
209 sprintf(hname, "TreeD%d", evnb);
210 gAlice->TreeD()->Write(hname,TObject::kOverwrite);
212 gAlice->TreeD()->Reset();
217 //___________________________________________
218 void AliMUONproto::SetPadSize(Int_t id, Int_t isec, Float_t p1, Float_t p2)
221 ((AliMUONChamber*)(*fChambers)[i])->SetPadSize(isec,p1,p2);
224 //___________________________________________
225 void AliMUONproto::SetChargeSlope(Int_t id, Float_t p1)
228 ((AliMUONChamber*) (*fChambers)[i])->SetChargeSlope(p1);
231 //___________________________________________
232 void AliMUONproto::SetChargeSpread(Int_t id, Float_t p1, Float_t p2)
235 ((AliMUONChamber*) (*fChambers)[i])->SetChargeSpread(p1,p2);
238 //___________________________________________
239 void AliMUONproto::SetMaxAdc(Int_t id, Float_t p1)
242 ((AliMUONChamber*) (*fChambers)[i])->SetMaxAdc(p1);
245 //___________________________________________
246 void AliMUONproto::CreateGeometry()
248 Int_t *idtmed = fIdtmed->GetArray()-1099;
252 Int_t idGas=idtmed[1100];
254 AliMUONChamber *iChamber;
255 //********************************************************************
256 // Prototype ALICE2 **
257 //********************************************************************
258 iChamber=(AliMUONChamber*) (*fChambers)[0];
262 const Float_t X_POS = 11*.975/2; //half x
263 const Float_t Y_POS = 18*.55/2; //half y
264 const Float_t Z_POS = 0.325;
270 gMC->Gsvolu("C01G", "BOX", idGas, bpar, 3);
271 gMC->Gspos("C01G", 1, "ALIC", -bpar[0], bpar[1], zpos, 0, "ONLY");
274 //___________________________________________
275 void AliMUONproto::CreateMaterials()
277 // *** DEFINITION OF AVAILABLE MUON MATERIALS ***
280 Float_t ag1[3] = { 39.95,12.01,16. };
281 Float_t zg1[3] = { 18.,6.,8. };
282 Float_t wg1[3] = { .8,.0667,.13333 };
283 Float_t dg1 = .001821;
285 Float_t epsil = .001; // Tracking precision,
286 Float_t tmaxfd = -20.; // Maximum angle due to field deflection
288 Int_t ISXFLD = gAlice->Field()->Integ();
289 Float_t SXMGMX = gAlice->Field()->Max();
292 // --- Define the various materials for GEANT ---
293 AliMixture(22, "ArCO2 80%$", ag1, zg1, dg1, 3, wg1);
297 AliMedium(1, "ARG_CO2 ", 22, 1, ISXFLD, SXMGMX, tmaxfd, fMaxStepGas,
298 fMaxDestepAlu, epsil, stmin);
300 //AliMedium(1, "AIR_CH_US ", 15, 1, ISXFLD, SXMGMX, tmaxfd, -1., -.3, epsil, stmin);
303 //___________________________________________
304 void AliMUONproto::Init()
306 printf("\n\n\n Start Init for Prototype ALICE2 - CPC chamber type\n\n\n");
309 // Initialize Tracking Chambers
312 ( (AliMUONChamber*) (*fChambers)[0])->Init();
315 // Set the chamber (sensitive region) GEANT identifier
316 ((AliMUONChamber*)(*fChambers)[0])->SetGid(gMC->VolId("C01G"));
318 printf("\n\n\n Finished Init for Prototype ALICE2 - CPC chamber type\n\n\n");
322 //___________________________________________
323 void AliMUONproto::StepManager()
332 Float_t destep, step;
334 static Float_t eloss, eloss2, xhit, yhit, tlength;
335 const Float_t big=1.e10;
338 static Float_t hits[14];
340 TClonesArray &lhits = *fHits;
343 // Set maximum step size for gas
344 // numed=gMC->GetMedium();
346 // Only charged tracks
347 if( !(gMC->TrackCharge()) ) return;
349 // Only gas gap inside chamber
350 // Tag chambers and record hits when track enters
352 id=gMC->CurrentVolID(copy);
354 for (Int_t i=1; i<=kNCH; i++) {
355 if(id==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()){
360 if (idvol == -1) return;
362 // Get current particle id (ipart), track position (pos) and momentum (mom)
363 gMC->TrackPosition(pos);
364 gMC->TrackMomentum(mom);
366 ipart = gMC->TrackPid();
367 //Int_t ipart1 = gMC->IdFromPDG(ipart);
368 //printf("ich, ipart %d %d \n",vol[0],ipart1);
371 // momentum loss and steplength in last step
372 destep = gMC->Edep();
373 step = gMC->TrackStep();
376 // record hits when track enters ...
377 if( gMC->IsTrackEntering()) {
378 gMC->SetMaxStep(fMaxStepGas);
379 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
380 Double_t rt = TMath::Sqrt(tc);
381 Double_t pmom = TMath::Sqrt(tc+mom[2]*mom[2]);
382 Double_t tx=mom[0]/pmom;
383 Double_t ty=mom[1]/pmom;
384 Double_t tz=mom[2]/pmom;
385 Double_t s=((AliMUONChamber*)(*fChambers)[idvol])
388 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
389 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
390 hits[0] = Float_t(ipart); // Geant3 particle type
391 hits[1] = pos[0]+s*tx; // X-position for hit
392 hits[2] = pos[1]+s*ty; // Y-position for hit
393 hits[3] = pos[2]+s*tz; // Z-position for hit
394 hits[4] = theta; // theta angle of incidence
395 hits[5] = phi; // phi angle of incidence
396 hits[8] = (Float_t) fNPadHits; // first padhit
397 hits[9] = -1; // last pad hit
400 hits[10] = mom[3]; // hit momentum P
401 hits[11] = mom[0]; // Px/P
402 hits[12] = mom[1]; // Py/P
403 hits[13] = mom[2]; // Pz/P
406 // phi angle of incidence
412 // Only if not trigger chamber
415 // Initialize hit position (cursor) in the segmentation model
416 ((AliMUONChamber*) (*fChambers)[idvol])
417 ->SigGenInit(pos[0], pos[1], pos[2]);
420 //printf("In the Trigger Chamber #%d\n",idvol-9);
426 // Calculate the charge induced on a pad (disintegration) in case
428 // Mip left chamber ...
429 if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
430 gMC->SetMaxStep(big);
434 // Only if not trigger chamber
437 MakePadHits(0.5*(xhit+pos[0]),0.5*(yhit+pos[1]),eloss,0.0,idvol);
442 if (fNPadHits > (Int_t)hits[8]) {
444 hits[9]= (Float_t) fNPadHits;
448 AliMUONHit(fIshunt,gAlice->CurrentTrack(),vol,hits);
451 // Check additional signal generation conditions
452 // defined by the segmentation
453 // model (boundary crossing conditions)
455 (((AliMUONChamber*) (*fChambers)[idvol])
456 ->SigGenCond(pos[0], pos[1], pos[2]))
458 ((AliMUONChamber*) (*fChambers)[idvol])
459 ->SigGenInit(pos[0], pos[1], pos[2]);
460 // printf("\n-> MakePadHits, reason special %d",ipart);
462 MakePadHits(0.5*(xhit+pos[0]),0.5*(yhit+pos[1]),eloss,0.0,idvol);
468 // nothing special happened, add up energy loss
475 //___________________________________________
476 void AliMUONproto::BuildGeometry()
481 const int kColorMUON = kBlue;
483 Top=gAlice->GetGeometry()->GetNode("alice");
487 float dx, dy, dz, zpos;
489 const Float_t cz[1]={975.};
497 new TBRIK("C_MUON101","Mother volume for Proto.","void",dx*2,dy*2,dz*2);
498 TBRIK* PROTO = new TBRIK("PROT", "Proto. Sens. Region","void",dx,dy,dz);
500 Node = new TNode("MUON101","ChamberNode","C_MUON101",0,0,zpos,"");
501 Node->SetLineColor(kColorMUON);
502 Node->SetVisibility(0);
505 Node = new TNode("MUON201", "Proto. Sens. Region Node", PROTO, -dx, dy, dz);
506 Node->SetLineColor(kColorMUON);
510 //_____________________________________________________________________________
511 void AliMUONproto::FindClusters(Int_t nev, Int_t last_entry)
515 // Loop on chambers and on cathode planes
517 for (Int_t icat = 0; icat < 2; icat++) {
518 gAlice->ResetDigits();
519 gAlice->TreeD()->GetEvent(last_entry+icat);
521 for (Int_t ich = 0; ich < kNTrackingCh; ich++) {
522 AliMUONChamber* iChamber=(AliMUONChamber*) (*fChambers)[ich];
523 TClonesArray *MUONdigits = this->DigitsAddress(ich);
524 if (MUONdigits == 0) continue;
526 // Get ready the current chamber stuff
529 AliMUONResponse* response = iChamber->ResponseModel();
530 AliMUONSegmentation* seg = iChamber->SegmentationModel(icat+1);
531 AliMUONClusterFinder* rec = iChamber->ReconstructionModel();
534 rec->SetSegmentation(seg);
535 rec->SetResponse(response);
536 rec->SetDigits(MUONdigits);
537 rec->SetChamber(ich);
538 rec->FindRawClusters();
542 fRch=RawClustAddress(ich);
548 gAlice->TreeR()->Fill();
554 sprintf(hname,"TreeR%d",nev);
555 gAlice->TreeR()->Write(hname);
556 gAlice->TreeR()->Reset();
559 #ifdef WE_FORGRET_THIS_SHIT
560 //_____________________________________________________________________________
561 void AliMUONproto::SetThreshold()
564 ifstream inputFile("/home/alice/guernane/aliroot/pro/MUON/crped190.dat", ios::in);
566 if (inputFile.fail()) {
567 cout << "Error opening file" << endl;
581 inputFile >> Ntrigger;
582 inputFile >> Nchannel;
586 Float_t ped0[Nchannel];
587 Float_t sig0[Nchannel];
588 Float_t ped1[Nchannel];
589 Float_t sig1[Nchannel];
592 for (Int_t i = 0; i < Nchannel-1; i++) {
599 for (Int_t i = 0; i < Nchannel-1; i++) {
600 inputFile >> ichannel;
601 inputFile >> ped0[i];
602 inputFile >> sig0[i];
603 inputFile >> ped1[i];
604 inputFile >> sig1[i];
605 fThreshold[i] = fNsigma*sig1[i];