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