Moving to the new VMC naming convention
[u/mrichter/AliRoot.git] / MUON / AliMUONproto.cxx
CommitLineData
a9e2aefa 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
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 **************************************************************************/
2ab0c725 15
88cb7938 16/* $Id$ */
a9e2aefa 17
18/*
19Implementation of MUON Chamber Prototype Class
20Both read digits from raw data or use the Monte-Carlo.
211-Feb-2000 Rachid GUERNANE, IPN Lyon guernane@ipnl.in2p3.fr
22*/
23
24
25////////////////////////////////////////////////
26// Manager and hits classes for set:PROTO //
27////////////////////////////////////////////////
28
29#include <TTUBE.h>
30#include <TBRIK.h>
31#include <TRotMatrix.h>
32#include <TNode.h>
33#include <TTree.h>
34#include <TRandom.h>
35#include <TObject.h>
36#include <TVector.h>
37#include <TObjArray.h>
38#include <TMinuit.h>
39#include <TParticle.h>
40#include <TROOT.h>
41#include <TFile.h>
42#include <TNtuple.h>
43#include <TCanvas.h>
44#include <TPad.h>
45#include <TDirectory.h>
46#include <TObjectTable.h>
47#include <AliPDG.h>
48#include <TArrayI.h>
49#include <AliDetector.h>
50
4bc2de9c 51#include "AliMUONChamber.h"
a9e2aefa 52#include "AliMUONproto.h"
4bc2de9c 53#include "AliMUONHit.h"
a9e2aefa 54#include "TTUBE.h"
55#include "AliMUONClusterFinder.h"
56#include "AliRun.h"
70479d0e 57#include "Riostream.h"
a9e2aefa 58#include "AliCallf77.h"
59#include "AliConst.h"
4bc2de9c 60//#include "chainalice2.h"
61#include "AliMUONSegmentationV0.h"
62//#include "AliMUONSegResV11.h"
a9e2aefa 63
64ClassImp(AliMUONproto)
65
66//___________________________________________
67AliMUONproto::AliMUONproto()
68 : AliMUON()
69{
70 cout << "\n Calling AliMUONproto constructor..." << endl;
71
72 //
73 //
74 //
75}
76
77//___________________________________________
78AliMUONproto::AliMUONproto(const char *name, const char *title)
79 : AliMUON(name,title)
80// AliMUON defines everything, but the chamber for NCH=1
81{
82
83// z-Positions of Chambers
84 const Float_t zch[1] = {975.};
85//
86// inner diameter
87 const Float_t dmi[1] = {0.};
88//
89// outer diameter
90 const Float_t dma[1] = {29.2};
91//
92//
93// Default Parameters for ALICE2 prototype
94
95//
96 (*fChambers)[0] = new AliMUONChamber();
97 AliMUONChamber* chamber = (AliMUONChamber*) (*fChambers)[0];
98 chamber->SetGid(0);
4bc2de9c 99 chamber->SetZ(zch[0]);
a9e2aefa 100//
101 chamber->InitGeo(zch[0]);
102 chamber->SetRInner(dmi[0]/2);
103 chamber->SetROuter(dma[0]/2);
104
105 for (int i = 0; i <= 99; i++) {
106 fThreshold[i] = 0.;
107 }
108
109}
110
4bc2de9c 111#ifdef WE_FORGET_FOR_THE_MOMENT
a9e2aefa 112//___________________________________________
113void AliMUONproto::GetRawDigits(Int_t evnb, Int_t *lptr, Int_t ilen) {
114
115 Int_t ip = 0;
116 Int_t equip = 0;
117 Int_t nochnl;
118 Int_t loop;
119 Int_t val;
120 Int_t itype;
121 Int_t id;
122 Int_t serial;
123 Int_t equiplen;
124
125 Int_t digits[5];
126
127 AliMUON* MUON = (AliMUON*)gAlice->GetModule("MUON");
4bc2de9c 128 AliMUONSegmentationV0* seg = (AliMUONSegmentationV0*) Chamber(0).SegmentationModel(1);
a9e2aefa 129
130 Int_t tracks[10];
131 Int_t charges[10];
132
133 for (Int_t i = 0; i < 10; i++) {
134 tracks[i] = 0;
135 charges[i] = 0;
136 }
137
138 Int_t ich = 0;
139
140 nwtype:
141
142 itype = lptr[ip++];
143 id = lptr[ip++];
144 equiplen = lptr[ip++];
145
146 if (equiplen < 0 ) return;
147
148 if (itype != (int)(0XCACCA008)) {
149 ip += equiplen;
150 if (ip < ilen-2) goto nwtype;
151 }
152 else {
153 serial = id >> 16;
154 equip = id & 0x1;
155 if ((serial == 190) && (equip == 1)) {
156 for (loop = 0; loop < equiplen; loop++) {
157 nochnl = (lptr[ip] & 0x7ff000 ) >> 12;
158 val = lptr[ip] & 0x3ff;
159 // fill digits from raw data according to cathode connexions
160 if (group[nochnl][2][1]!=0)
161 digits[0] = group[nochnl][2][1] - 12;
162 else if (group[nochnl][1][1]!=0)
163 digits[0] = group[nochnl][1][1] - 12;
164 else
165 digits[0] = group[nochnl][0][1] - 12;
166 digits[1] = group[nochnl][0][0];
167 if (digits[0] != seg->Ix(digits[0], digits[1]))
168 printf("Pb pour ix=%d,iy=%d\n", digits[0], digits[1]);
169 digits[2] = val;
170 digits[3] = 0;
171 digits[4] = 0;
172 if (digits[2] >= fThreshold[nochnl])
173 MUON->AddDigits(ich, tracks, charges, digits);
174
175 ip++;
176 }
177 }
178 else
179 ip += equiplen;
180
181 if (ip < ilen-2) goto nwtype;
182 }
183
184 gAlice->TreeD()->Fill();
185 MUON->ResetDigits();
186
187 gAlice->TreeD()->Fill();
188 MUON->ResetDigits();
189
190 char hname[30];
191 sprintf(hname, "TreeD%d", evnb);
2ab0c725 192 gAlice->TreeD()->Write(hname,TObject::kOverwrite);
a9e2aefa 193 // reset tree
194 gAlice->TreeD()->Reset();
195
196}
4bc2de9c 197#endif
a9e2aefa 198
199//___________________________________________
200void AliMUONproto::SetPadSize(Int_t id, Int_t isec, Float_t p1, Float_t p2)
201{
202 Int_t i=2*(id-1);
203 ((AliMUONChamber*)(*fChambers)[i])->SetPadSize(isec,p1,p2);
204}
205
206//___________________________________________
207void AliMUONproto::SetChargeSlope(Int_t id, Float_t p1)
208{
209 Int_t i=2*(id-1);
210 ((AliMUONChamber*) (*fChambers)[i])->SetChargeSlope(p1);
211}
212
213//___________________________________________
214void AliMUONproto::SetChargeSpread(Int_t id, Float_t p1, Float_t p2)
215{
216 Int_t i=2*(id-1);
217 ((AliMUONChamber*) (*fChambers)[i])->SetChargeSpread(p1,p2);
218}
219
220//___________________________________________
221void AliMUONproto::SetMaxAdc(Int_t id, Float_t p1)
222{
223 Int_t i=2*(id-1);
224 ((AliMUONChamber*) (*fChambers)[i])->SetMaxAdc(p1);
225}
226
227//___________________________________________
228void AliMUONproto::CreateGeometry()
229{
230 Int_t *idtmed = fIdtmed->GetArray()-1099;
231//
232 Float_t zpos;
233 Float_t bpar[3];
234 Int_t idGas=idtmed[1100];
235
236 AliMUONChamber *iChamber;
237//********************************************************************
238// Prototype ALICE2 **
239//********************************************************************
240 iChamber=(AliMUONChamber*) (*fChambers)[0];
241
4bc2de9c 242 zpos=iChamber->Z();
a9e2aefa 243//
244 const Float_t X_POS = 11*.975/2; //half x
245 const Float_t Y_POS = 18*.55/2; //half y
246 const Float_t Z_POS = 0.325;
247
248 bpar[0] = X_POS;
249 bpar[1] = Y_POS;
250 bpar[2] = Z_POS;
251
252 gMC->Gsvolu("C01G", "BOX", idGas, bpar, 3);
253 gMC->Gspos("C01G", 1, "ALIC", -bpar[0], bpar[1], zpos, 0, "ONLY");
254}
255
256//___________________________________________
257void AliMUONproto::CreateMaterials()
258{
259 // *** DEFINITION OF AVAILABLE MUON MATERIALS ***
260 //
261 // Ar-CO2 gas
262 Float_t ag1[3] = { 39.95,12.01,16. };
263 Float_t zg1[3] = { 18.,6.,8. };
264 Float_t wg1[3] = { .8,.0667,.13333 };
265 Float_t dg1 = .001821;
266
267 Float_t epsil = .001; // Tracking precision,
268 Float_t tmaxfd = -20.; // Maximum angle due to field deflection
269 Float_t stmin = -.8;
270 Int_t ISXFLD = gAlice->Field()->Integ();
271 Float_t SXMGMX = gAlice->Field()->Max();
272
273 //
274 // --- Define the various materials for GEANT ---
275 AliMixture(22, "ArCO2 80%$", ag1, zg1, dg1, 3, wg1);
276
277 //
278 // Ar/CO2
279 AliMedium(1, "ARG_CO2 ", 22, 1, ISXFLD, SXMGMX, tmaxfd, fMaxStepGas,
280 fMaxDestepAlu, epsil, stmin);
281 // Air
282 //AliMedium(1, "AIR_CH_US ", 15, 1, ISXFLD, SXMGMX, tmaxfd, -1., -.3, epsil, stmin);
283}
284
285//___________________________________________
286void AliMUONproto::Init()
287{
288 printf("\n\n\n Start Init for Prototype ALICE2 - CPC chamber type\n\n\n");
289
290 //
291 // Initialize Tracking Chambers
292 //
293
294 ( (AliMUONChamber*) (*fChambers)[0])->Init();
295
296 //
297 // Set the chamber (sensitive region) GEANT identifier
a9e2aefa 298 ((AliMUONChamber*)(*fChambers)[0])->SetGid(gMC->VolId("C01G"));
299
300 printf("\n\n\n Finished Init for Prototype ALICE2 - CPC chamber type\n\n\n");
301
302}
303
304//___________________________________________
305void AliMUONproto::StepManager()
306{
307 Int_t copy, id;
308 static Int_t idvol;
309 static Int_t vol[2];
310 Int_t ipart;
311 TLorentzVector pos;
312 TLorentzVector mom;
313 Float_t theta,phi;
314 Float_t destep, step;
315
316 static Float_t eloss, eloss2, xhit, yhit, tlength;
317 const Float_t big=1.e10;
318
319 // modifs perso
320 static Float_t hits[14];
321
322 TClonesArray &lhits = *fHits;
323
324 //
325 // Set maximum step size for gas
326 // numed=gMC->GetMedium();
327 //
328 // Only charged tracks
329 if( !(gMC->TrackCharge()) ) return;
330 //
331 // Only gas gap inside chamber
332 // Tag chambers and record hits when track enters
333 idvol=-1;
334 id=gMC->CurrentVolID(copy);
335
4bc2de9c 336 for (Int_t i=1; i<=kNCH; i++) {
a9e2aefa 337 if(id==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()){
338 vol[0]=i;
339 idvol=i-1;
340 }
341 }
342 if (idvol == -1) return;
343 //
344 // Get current particle id (ipart), track position (pos) and momentum (mom)
345 gMC->TrackPosition(pos);
346 gMC->TrackMomentum(mom);
347
348 ipart = gMC->TrackPid();
349 //Int_t ipart1 = gMC->IdFromPDG(ipart);
350 //printf("ich, ipart %d %d \n",vol[0],ipart1);
351
352 //
353 // momentum loss and steplength in last step
354 destep = gMC->Edep();
355 step = gMC->TrackStep();
356
357 //
358 // record hits when track enters ...
359 if( gMC->IsTrackEntering()) {
360 gMC->SetMaxStep(fMaxStepGas);
361 Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
362 Double_t rt = TMath::Sqrt(tc);
363 Double_t pmom = TMath::Sqrt(tc+mom[2]*mom[2]);
364 Double_t tx=mom[0]/pmom;
365 Double_t ty=mom[1]/pmom;
366 Double_t tz=mom[2]/pmom;
367 Double_t s=((AliMUONChamber*)(*fChambers)[idvol])
4bc2de9c 368 ->ResponseModel()
a9e2aefa 369 ->Pitch()/tz;
370 theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
371 phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
372 hits[0] = Float_t(ipart); // Geant3 particle type
373 hits[1] = pos[0]+s*tx; // X-position for hit
374 hits[2] = pos[1]+s*ty; // Y-position for hit
375 hits[3] = pos[2]+s*tz; // Z-position for hit
376 hits[4] = theta; // theta angle of incidence
377 hits[5] = phi; // phi angle of incidence
4bc2de9c 378 hits[8] = (Float_t) fNPadHits; // first padhit
a9e2aefa 379 hits[9] = -1; // last pad hit
380
381 // modifs perso
382 hits[10] = mom[3]; // hit momentum P
383 hits[11] = mom[0]; // Px/P
384 hits[12] = mom[1]; // Py/P
385 hits[13] = mom[2]; // Pz/P
386 // fin modifs perso
387
388 // phi angle of incidence
389 tlength = 0;
390 eloss = 0;
391 eloss2 = 0;
392 xhit = pos[0];
393 yhit = pos[1];
394 // Only if not trigger chamber
395 if(idvol<10) {
396 //
397 // Initialize hit position (cursor) in the segmentation model
398 ((AliMUONChamber*) (*fChambers)[idvol])
399 ->SigGenInit(pos[0], pos[1], pos[2]);
400 } else {
401 //geant3->Gpcxyz();
402 //printf("In the Trigger Chamber #%d\n",idvol-9);
403 }
404 }
405 eloss2+=destep;
406
407 //
408 // Calculate the charge induced on a pad (disintegration) in case
409 //
410 // Mip left chamber ...
411 if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
412 gMC->SetMaxStep(big);
413 eloss += destep;
414 tlength += step;
415
416 // Only if not trigger chamber
417 if(idvol<10) {
418 if (eloss > 0)
4bc2de9c 419 MakePadHits(0.5*(xhit+pos[0]),0.5*(yhit+pos[1]),eloss,0.0,idvol);
a9e2aefa 420 }
421
422 hits[6]=tlength;
423 hits[7]=eloss2;
4bc2de9c 424 if (fNPadHits > (Int_t)hits[8]) {
a9e2aefa 425 hits[8]= hits[8]+1;
4bc2de9c 426 hits[9]= (Float_t) fNPadHits;
a9e2aefa 427 }
428
429 new(lhits[fNhits++])
642f15cf 430 AliMUONHit(fIshunt,gAlice->GetCurrentTrackNumber(),vol,hits);
a9e2aefa 431 eloss = 0;
432 //
433 // Check additional signal generation conditions
434 // defined by the segmentation
435 // model (boundary crossing conditions)
436 } else if
437 (((AliMUONChamber*) (*fChambers)[idvol])
438 ->SigGenCond(pos[0], pos[1], pos[2]))
439 {
440 ((AliMUONChamber*) (*fChambers)[idvol])
441 ->SigGenInit(pos[0], pos[1], pos[2]);
442// printf("\n-> MakePadHits, reason special %d",ipart);
443 if (eloss > 0)
4bc2de9c 444 MakePadHits(0.5*(xhit+pos[0]),0.5*(yhit+pos[1]),eloss,0.0,idvol);
a9e2aefa 445 xhit = pos[0];
446 yhit = pos[1];
447 eloss = destep;
448 tlength += step ;
449 //
450 // nothing special happened, add up energy loss
451 } else {
452 eloss += destep;
453 tlength += step ;
454 }
455}
456
457//___________________________________________
458void AliMUONproto::BuildGeometry()
459{
460 TNode *Node;
461 TNode *Top;
462
463 const int kColorMUON = kBlue;
464 //
465 Top=gAlice->GetGeometry()->GetNode("alice");
466 //
467 //
468 //
469 float dx, dy, dz, zpos;
470
471 const Float_t cz[1]={975.};
472
473 zpos=cz[0];
474
475 dx=11*.975/2;
476 dy=18*.55/2;
477 dz=0.325;
478
479 new TBRIK("C_MUON101","Mother volume for Proto.","void",dx*2,dy*2,dz*2);
480 TBRIK* PROTO = new TBRIK("PROT", "Proto. Sens. Region","void",dx,dy,dz);
481 Top->cd();
482 Node = new TNode("MUON101","ChamberNode","C_MUON101",0,0,zpos,"");
483 Node->SetLineColor(kColorMUON);
484 Node->SetVisibility(0);
485 fNodes->Add(Node);
486 Node->cd();
487 Node = new TNode("MUON201", "Proto. Sens. Region Node", PROTO, -dx, dy, dz);
488 Node->SetLineColor(kColorMUON);
489
490}
491
492//_____________________________________________________________________________
493void AliMUONproto::FindClusters(Int_t nev, Int_t last_entry)
494{
495
496//
497// Loop on chambers and on cathode planes
498//
499 for (Int_t icat = 0; icat < 2; icat++) {
500 gAlice->ResetDigits();
501 gAlice->TreeD()->GetEvent(last_entry+icat);
502
4bc2de9c 503 for (Int_t ich = 0; ich < kNTrackingCh; ich++) {
a9e2aefa 504 AliMUONChamber* iChamber=(AliMUONChamber*) (*fChambers)[ich];
505 TClonesArray *MUONdigits = this->DigitsAddress(ich);
506 if (MUONdigits == 0) continue;
507 //
508 // Get ready the current chamber stuff
509 //
510
4bc2de9c 511 AliMUONResponse* response = iChamber->ResponseModel();
512 AliMUONSegmentation* seg = iChamber->SegmentationModel(icat+1);
513 AliMUONClusterFinder* rec = iChamber->ReconstructionModel();
a9e2aefa 514
515 if (seg) {
516 rec->SetSegmentation(seg);
517 rec->SetResponse(response);
518 rec->SetDigits(MUONdigits);
519 rec->SetChamber(ich);
520 rec->FindRawClusters();
521 }
522
523 TClonesArray *fRch;
524 fRch=RawClustAddress(ich);
525 fRch->Sort();
526 // it seems to work
527
528 } // for ich
529 // fill the tree
530 gAlice->TreeR()->Fill();
531
532 ResetRawClusters();
533 } // for icat
534
535 char hname[30];
536 sprintf(hname,"TreeR%d",nev);
537 gAlice->TreeR()->Write(hname);
538 gAlice->TreeR()->Reset();
539}
540
4bc2de9c 541#ifdef WE_FORGRET_THIS_SHIT
a9e2aefa 542//_____________________________________________________________________________
543void AliMUONproto::SetThreshold()
544{
545
546 ifstream inputFile("/home/alice/guernane/aliroot/pro/MUON/crped190.dat", ios::in);
547
548 if (inputFile.fail()) {
549 cout << "Error opening file" << endl;
550 exit(2);
551 }
552
553 char buff[32];
554 Int_t Serial;
555 Int_t Ntrigger;
556 Int_t Nchannel;
557 Int_t i1;
558 Int_t i2;
559
560 inputFile >> buff;
561
562 inputFile >> Serial;
563 inputFile >> Ntrigger;
564 inputFile >> Nchannel;
565 inputFile >> i1;
566 inputFile >> i2;
567
568 Float_t ped0[Nchannel];
569 Float_t sig0[Nchannel];
570 Float_t ped1[Nchannel];
571 Float_t sig1[Nchannel];
572 Int_t ichannel;
573
574 for (Int_t i = 0; i < Nchannel-1; i++) {
575 ped0[i] = 0;
576 sig0[i] = 0;
577 ped1[i] = 0;
578 sig1[i] = 0;
579 }
580
581 for (Int_t i = 0; i < Nchannel-1; i++) {
582 inputFile >> ichannel;
583 inputFile >> ped0[i];
584 inputFile >> sig0[i];
585 inputFile >> ped1[i];
586 inputFile >> sig1[i];
587 fThreshold[i] = fNsigma*sig1[i];
588 }
589
590 inputFile.close();
591}
6b677e96 592#endif