BmnRoot 26.10.0-18-g4e3414f9
Software framework of the BM@N experiment
Loading...
Searching...
No Matches
ResidOk_r7.C
Go to the documentation of this file.
1#if !defined(__CLING__) || defined(__MAKECLING__)
2
3// CBM includes
4#include "BmnFieldConst.h"
5#include "BmnFieldMap.h"
6#include "BmnFileSource.h"
7#include "BmnNewFieldMap.h"
8#include "CbmKF.h"
9#include "CbmKFParticle.h"
11#include "CbmKFTrack.h"
12#include "CbmMCTrack.h"
13#include "CbmStsCluster.h"
14#include "CbmStsDigi.h"
15#include "CbmStsDigiPar.h"
16#include "CbmStsHit.h"
17#include "CbmStsKFTrackFitter.h"
18#include "CbmStsPoint.h"
19#include "CbmStsSector.h"
20#include "CbmStsSensor.h"
21#include "CbmStsStation.h"
22#include "CbmStsTrack.h"
23#include "CbmTrackMatch.h"
24#include "CbmVertex.h"
25#include "FairParAsciiFileIo.h"
26#include "FairParRootFileIo.h"
27#include "FairRunAna.h"
28#include "FairRuntimeDb.h"
29#include "FairTrackParam.h"
30#include "TKey.h"
31#include "UniRun.h"
32
33// ROOT includes
34#include <TBranch.h>
35#include <TClonesArray.h>
36#include <TDatabasePDG.h>
37// #include <TDecompLU.h>
38// #include <TDecompSVD.h>
39#include "BmnEventHeader.h"
40#include "TCanvas.h"
41
42#include <Riostream.h>
43#include <TFile.h>
44#include <TGeoCone.h>
45#include <TGeoManager.h>
46#include <TGeoVolume.h>
47#include <TH1.h>
48#include <TH2.h>
49#include <TMatrixD.h>
50#include <TROOT.h>
51#include <TRandom.h>
52#include <TSystem.h>
53#include <TTree.h>
54#include <TVector3.h>
55#include <TVectorD.h>
56#include <map>
57#include <set>
58#include <vector>
59
60#endif
61
62// #include "Chain.C"
63Bool_t GetBmnGeom(const char* fileName);
64void ApplyAlignment();
66
67const Int_t ip_r7[6][2] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11};
68
69// ----------------------------------------------------------------------------
70
71void ResidOk_r7(TString inFile = "", TString histFile = "", Double_t Mscale = 1.366666667)
72{
73
74 Mscale = 1250 / 900.0;
75 // Get hit residuals w.r.t. tracks
76 // (exclude each hit from the track and refit it to avoid bias)
77
78 // TObjArray *Hlist = new TObjArray();
79 Bool_t MC = false;
80
81 vector<Int_t> effstaUp = {0, 0, 0, 0, 0, 0, 0};
82 vector<Int_t> effstaDown = {0, 0, 0, 0, 0, 0, 0};
83
84 // TString inFile="1589_testSI.root";//"1882_lorentz.root";
85 TString outFile = "aaa.root";
86 // TString histFile = "resid.histo_1589_testSI.root";
87 TFile out(histFile, "RECREATE");
88
89 Float_t xext, yext, zhit, dXext, dYext, prec, chi2, xsmo, ysmo, dXsmo, dYsmo, dXZ, dYZ, precPT, precPX, precNOf,
90 trkAng, Qf, Qb;
91 Int_t flag, ista, nhits, nsecGEM, evNo, nsecGEMhot;
92 Float_t sectors[8];
93 Int_t nsec;
94 TTree* txyz = new TTree("hits", "Hits");
95 txyz->Branch("evNo", &evNo, "evNo/I");
96 txyz->Branch("xext", &xext, "xext/F");
97 txyz->Branch("yext", &yext, "yext/F");
98 txyz->Branch("zhit", &zhit, "zhit/F");
99 txyz->Branch("dXext", &dXext, "dXext/F");
100 txyz->Branch("dYext", &dYext, "dYext/F");
101 txyz->Branch("xsmo", &xsmo, "xsmo/F");
102 txyz->Branch("ysmo", &ysmo, "ysmo/F");
103 txyz->Branch("dXsmo", &dXsmo, "dXsmo/F");
104 txyz->Branch("dYsmo", &dYsmo, "dYsmo/F");
105 // txyz->Branch("dXZ",&dXZ,"dXZ/F");
106 // txyz->Branch("dYZ",&dYZ,"dYZ/F");
107 txyz->Branch("prec", &prec, "prec/F");
108 // txyz->Branch("precPT",&precPT,"precPT/F");
109 // txyz->Branch("precPX",&precPX,"precPX/F");
110 // txyz->Branch("precNOf",&precNOf,"precNOf/F");
111 txyz->Branch("chi2", &chi2, "chi2/F");
112 txyz->Branch("ista", &ista, "ista/I");
113 txyz->Branch("nhits", &nhits, "nhits/I");
114 txyz->Branch("nsec", &nsec, "nsec/I");
115 txyz->Branch("nsecGEM", &nsecGEM, "nsecGEM/I");
116 // txyz->Branch("nsecGEMhot",&nsecGEMhot,"nsecGEMhot/I");
117 txyz->Branch("trkAng", &trkAng, "trkAng/F");
118 txyz->Branch("Qf", &Qf, "Qf/F");
119 txyz->Branch("Qb", &Qb, "Qb/F");
120
121 TString parFile = inFile;
122 // fRun->Init();
123
124 // FairRootManager* ioman = FairRootManager::Instance();
125 // if ( ! ioman ) Fatal("Init", "No FairRootManager");
126 // TClonesArray* fHit = (TClonesArray*) ioman->GetObject("StsHit");//StsPoint
127
128 // fRun->SetOutputFile(outFile);
129
130 TList* parFileList = new TList();
131
132 TObjString stsDigiFile = "sts_v1_BMN_SI_GEM.digi.par"; // sts_v16BMN_gem.digi.par
133 // TObjString stsDigiFile2 = "/home/pokat/bmnM/bmnR/AZ_data/sts/GEMS_v16_AZ_3mm.root";
134 parFileList->Add(&stsDigiFile);
135
136 GetBmnGeom("SIGEMS_r7.root"); //"GEMS_RunWinter2016_old.root");
137
139
140 FairRunAna* fRun = new FairRunAna();
141 fRun->SetSource(new BmnFileSource(inFile));
142 fRun->SetOutputFile(outFile);
143
144 FairRuntimeDb* rtdb = fRun->GetRuntimeDb();
145
146 FairParRootFileIo* parIo1 = new FairParRootFileIo();
147 FairParAsciiFileIo* parIo2 = new FairParAsciiFileIo();
148 parIo1->open(parFile.Data());
149 parIo2->open(parFileList, "in");
150 gSystem->Unlink(parIo2->getFilename()); // drop concatenated all_<pid>.par
151 rtdb->setFirstInput(parIo1);
152 rtdb->setSecondInput(parIo2);
153 rtdb->setOutput(parIo1);
154 rtdb->saveOutput();
155
156 // ----- Field map -------------------------------------------
157 /*
158 Int_t run_period=1, run_number=065;
159 UniRun* pCurrentRun = UniRun::GetRun(run_period, run_number);
160 if (pCurrentRun == 0) {
161 exit(-2);
162 }
163 */
164 /* Double_t fieldScale = 0.01;
165 double map_current = 900.0;
166 //int* current_current = pCurrentRun->GetFieldCurrent();
167 //int* current_current = NULL;
168 int current_current[1] = {500};
169 if (current_current == NULL) {
170 fieldScale = 0.01;
171
172 } else
173 fieldScale = (*current_current) / map_current;
174 // cout<<" SCALE: "<<fieldScale<<endl;
175 fieldScale=TMath::Max(fieldScale,0.01);
176 */
177 // if(!MC){
178 BmnNewFieldMap* magField = new BmnNewFieldMap("field_sp41v4_ascii_Extrap.dat");
179 magField->SetScale(
180 Mscale); // 1.366666667);//1.777777778);//0.666666667);//2.0);//1.777777778);//2.0);//fieldScale);1.333333333
181 magField->Init();
182 fRun->SetField(magField);
183 cout << "FIELD SCALE: " << magField->GetScale() << endl;
184
185 // BmnFieldConst* magField = new BmnFieldConst();
186 // magField->SetFieldRegion(-3000., 3000., -3000., 3000., -3000., 3000);
187 // magField->SetField(0, 0.01 , 0.);
188 // magField->Init();
189 // ------------------------------------------------------------------------
190
191 cout << "FIELD SCALE: " << magField->GetScale() << endl;
192 /*} else {
193 BmnFieldConst* magField = new BmnFieldConst();
194 magField->SetFieldRegion(-3000., 3000., -3000., 3000., -3000., 3000);
195 magField->SetField(0, 5.9 , 0.);
196
197 magField->Init();
198 fRun->SetField(magField);
199 }
200 */
201
202 FairTask* pKF = new CbmKF();
203 fRun->AddTask(pKF);
204
205 fRun->Init();
206
207 // TClonesArray *pSTS = (TClonesArray *) FairRootManager::Instance()->GetObject("StsPoint");
208 // TClonesArray *trMC = (TClonesArray *) FairRootManager::Instance()->GetObject("MCTrack");
209 // AZ TChain *rec = (TChain*) gROOT->FindObject("cbmout");
210 TChain* rec = FairRootManager::Instance()->GetInChain();
211 TClonesArray* stsTrack = (TClonesArray*)FairRootManager::Instance()->GetObject("StsTrack");
212 // TClonesArray *trMatch = (TClonesArray *) FairRootManager::Instance()->GetObject("StsTrackMatch");
213 TClonesArray* stsHit = (TClonesArray*)FairRootManager::Instance()->GetObject("StsHit");
214 TClonesArray* stsDigi = (TClonesArray*)FairRootManager::Instance()->GetObject("StsDigi");
215 TClonesArray* stsCluster = (TClonesArray*)FairRootManager::Instance()->GetObject("StsCluster");
216 CbmVertex* primVtx = (CbmVertex*)FairRootManager::Instance()->GetObject("PrimaryVertex.");
217 TClonesArray* BmnH = (TClonesArray*)FairRootManager::Instance()->GetObject("EventHeaderBmn");
218 CbmStsKFTrackFitter fitter;
219 fitter.Init();
220
221 Int_t events = rec->GetEntries();
222 Int_t events0 = events;
223 // if (nRead != 0) events = TMath::Min (events, nRead);
224 // cout << " Number of events = " << events << " N of StsHits: "<<stsHit->GetEntries()<< endl;
225
226 // Get STS geometry parameter container
227 // CbmGeoStsPar *geoPar = (CbmGeoStsPar*) rtdb->getContainer("CbmGeoStsPar");
228 FairRuntimeDb* db = fRun->GetRuntimeDb();
229 if (!db)
230 Fatal("SetParContainers", "No runtime database");
231 // Get STS digitisation parameter container
232 CbmStsDigiPar* digiPar = (CbmStsDigiPar*)db->getContainer("CbmStsDigiPar");
233 digiPar->init(parIo2);
234
235 // cout<<"digiPAR: "<<digiPar->GetNStations()<<endl;
236
237 CbmStsDigiScheme* digiScheme = new CbmStsDigiScheme();
238 digiScheme->Init(NULL, digiPar);
239
240 // Int_t nStat = digiScheme->GetNStations();
241 // Double_t *zStat = new Double_t [nStat];
242 // for (Int_t i = 0; i < nStat; ++i) zStat[i] = digiScheme->GetStationByNr(i+1)->GetZ();
243
244 FairTrackParam param, pIter1;
245
246 for (Int_t iev = 0; iev < events; ++iev) {
247 if (iev % 10000 == 0)
248 cout << iev << endl;
249 // if (iev < nBeg) continue;
250
251 FairRootManager::Instance()->ReadEvent(iev);
252 // continue; //!!!!
253 evNo = iev;
254 Int_t nPoints = stsHit->GetEntriesFast(); // pSTS->GetEntriesFast();
255 Int_t nTracks = stsTrack->GetEntriesFast(); // trSTS->GetEntriesFast();
256 // CbmStsDigi *digi = (CbmStsDigi*) stsDigi->UncheckedAt(iev);
257 // cout << ", Mult: " << stsTrack->GetEntriesFast() << " hits: "<<nPoints<<"\n";
258 // if( nTracks<1 ) continue;
259 //-----------------control tracks from PV---------------------------
260 /* Int_t T0 = ((BmnEventHeader *)BmnH->At(0))->GetTrigT0();
261 Int_t BC2 = ((BmnEventHeader *)BmnH->At(0))->GetTrigBC2();
262 Int_t VETO = ((BmnEventHeader *)BmnH->At(0))->GetTrigVETO();
263 Int_t BD = ((BmnEventHeader *)BmnH->At(0))->GetTrigBD();
264
265 if(T0!=1 || BC2!=1 || VETO!=0) continue;
266
267 //cout<<"to: "<<T0<<" bc2: "<<BC2<<" veto: "<<VETO<<" bd: "<<BD<<endl;
268
269 //if( nTracks<1 ) continue;
270 TVector3 PVpos(primVtx->GetX(), primVtx->GetY(), primVtx->GetZ());
271 TVector3 PVexp(0, 0, -24.0);
272 Int_t nhitsLarger3=0, fromPV=0;
273 if(nTracks<2) continue;
274 if ((PVpos - PVexp).Mag() > 5) continue; //1.8
275
276
277 for (Int_t i = 0; i < nTracks; ++i) {
278 CbmStsTrack *tr = (CbmStsTrack *)stsTrack->UncheckedAt(i);
279 Int_t nhits1 = tr->GetNStsHits();
280 fromPV++;
281 if(nhits1>3) nhitsLarger3++;
282 }
283
284 if(fromPV<2 && nhitsLarger3<2) continue;
285 */
286
287 for (Int_t i = 0; i < nTracks; ++i) {
288 CbmStsTrack* tr = (CbmStsTrack*)stsTrack->UncheckedAt(i);
289 // cout<<"\n"<<endl;
290 nhits = tr->GetNStsHits();
291 /* if(nhits<4) continue;
292
293
294 CbmKFTrack kftr = CbmKFTrack(*tr);
295 FairTrackParam param;
296 kftr.Extrapolate(primVtx->GetZ());
297 kftr.GetTrackParam(param);
298
299 // distance betwee track fit to PV && PV postion
300 Double_t dx = param.GetX() - primVtx->GetX();
301 Double_t dy = param.GetY() - primVtx->GetY();
302 Double_t dz = param.GetZ() - primVtx->GetZ();
303
304 Double_t dist = TMath::Sqrt (dx * dx + dy * dy + dz * dz);
305
306 if(dist>0.5) continue; */
307 //-----------------control tracks from PV---------------------------
308
309 // if (nhits !=6) continue;
310 // if (nhits !=6 || stsDigi->GetEntries()>80) continue; // not enough hits
311 // if (nhits !=6 || stsCluster->GetEntries()!=12) continue; // not enough hits
312 TVector3 trackPos[2];
313 tr->GetParamFirst()->Position(trackPos[0]);
314 tr->GetParamLast()->Position(trackPos[1]);
315
316 CbmStsTrack* track = new CbmStsTrack(*tr);
317 // Int_t pdgDeuteron=1000010020;
318 // fitter.DoFit(track, pdgDeuteron );
319
320 prec = 1. / track->GetParamFirst()->GetQp();
321 chi2 = track->GetChi2() / track->GetNDF();
322 TVector3 mom, mom2;
323 // cout<<" prec: "<<1/track->GetParamFirst()->GetQp() <<"x,y,z,tx,ty:
324 // "<<track->GetParamFirst()->GetX()<<","<<track->GetParamFirst()->GetY()<<","<<track->GetParamFirst()->GetZ()<<","
325 // <<track->GetParamFirst()->GetTx()<<","<<track->GetParamFirst()->GetTy()<<endl;
326
327 // cout<<" prec: "<<1/track->GetParamLast()->GetQp() <<"x,y,z,tx,ty:
328 // "<<track->GetParamLast()->GetX()<<","<<track->GetParamLast()->GetY()<<","<<track->GetParamLast()->GetZ()<<","
329 // <<track->GetParamLast()->GetTx()<<","<<track->GetParamLast()->GetTy()<<endl;
330
331 // track->GetParamFirst()->Momentum(mom);
332 // precPT=mom.Pt();
333
334 // CbmKFTrack track2 = CbmKFTrack(*track);
335 // track2.SetPID(pdgDeuteron);
336 // CbmKFParticle *deuteron= new CbmKFParticle(&track2);
337
338 // FairTrackParam param_ex;
339 // track2.Extrapolate(-100);
340 // track2.Smooth(primVtx->GetZ());
341 // track2.GetTrackParam(param_ex);
342
343 // param_ex.Momentum(mom2);
344 // precNOf=mom2.Pt();
345 // precPX=mom2.Px();
346
347 // txyz->Fill();
348 /*
349 CbmStsTrack *tr = (CbmStsTrack*) trSTS->UncheckedAt(i);
350 CbmTrackMatch *trM = (CbmTrackMatch*) trMatch->UncheckedAt(i);
351 Int_t id = trM->GetMCTrackId();
352 if (id < 0) continue;
353 //if (trM->GetNofTrueHits() <= trM->GetNofFakeHits() + trM->GetNofWrongHits()) continue;
354
355 CbmMCTrack *mcTr = (CbmMCTrack*) trMC->UncheckedAt(id);
356 TVector3 trackPos[2];
357 tr->GetParamFirst()->Position(trackPos[0]);
358 tr->GetParamLast()->Position(trackPos[1]);
359 prec = 1./tr->GetParamFirst()->GetQp();
360 chi2 = tr->GetChi2() / tr->GetNDF();
361 nhits = tr->GetNStsHits();
362 */
363
364 //-----------------------GP extrapolate----------------------
365
366 //-----------------------end---------------------------------
367 // cout<<""<<endl;
368 nsec = 0;
369 nsecGEM = 0;
370 nsecGEMhot = 0;
371 std::vector<int> v1 = {1, 1, 1, 1, 1, 1};
372 std::vector<int> v2 = {1, 1, 1, 3, 3, 3};
373 std::vector<int> v3 = {0, 0, 0, 0, 0, 0};
374 // for (Int_t Niter = 0; Niter < 2; ++Niter)
375 /*for (Int_t iHit = 0; iHit < nhits; ++iHit) {
376 CbmStsHit *hit = (CbmStsHit*) stsHit->UncheckedAt(track->GetStsHitIndex(iHit));
377 Int_t secn=hit->GetSectorNr();
378 Int_t stan=hit->GetStationNr();
379 v3[iHit]=secn;
380 } */
381
382 // if(v3!=v1 && v3!=v2) continue;
383 // cout<<"\n";
384 for (Int_t iHit = 0; iHit < nhits; ++iHit) {
385 CbmStsHit* hit = (CbmStsHit*)stsHit->UncheckedAt(track->GetStsHitIndex(iHit));
386
387 if (!hit) {
388 cout << "-E- CbmStsMatchTracks::Exec: " << "No StsHit " << iHit << " for track " << i << endl;
389 // warn = kTRUE;
390 continue;
391 }
392 TVector3 pos, residXY;
393 hit->Position(pos);
394 // if(hit->GetStationNr()==1) pos[2]+=0.022;
395 // if(hit->GetStationNr()==2) pos[2]-=0.022;
396 // if(hit->GetStationNr()==3) pos[2]+=0.022;
397 // cout<<" nhits: "<<nhits<<" iHit: "<<iHit<<" posZ: "<<pos[2]<<" index:
398 // "<<track->GetStsHitIndex(iHit)<<endl;
399
400 // if(Niter==1){
401 // pos[0]= pos[0]-param.GetTx()/2;
402 // }
403
404 // residXY=GetStrightResid( trackPos[0], trackPos[1], pos );
405
406 // dXZ=residXY[0];
407 // dYZ=residXY[1];
408
409 // Remove this hit and refit
410 CbmStsTrack trackC(*track);
411 TArrayI& hits = *trackC.GetStsHits();
412 hits.Set(nhits - 1);
413 Int_t indx = 0;
414 for (Int_t j = 0; j < nhits; ++j) {
415 if (j == iHit)
416 continue; // UNincludes hit
417 // if (j == -1) continue; //includes hit
418 hits[indx++] = track->GetStsHitIndex(j);
419 }
420 CbmKFTrack kftr = CbmKFTrack(trackC);
421 fitter.SetKFHits(kftr, &trackC);
422 kftr.Fit(kTRUE); // downstream
423 kftr.GetTrackParam(param);
424 // cout << track->GetParamFirst()->GetX() << " " << track->GetParamFirst()->GetZ() << " "
425 // << param.GetX() << " " << param.GetZ() << endl;
426 // cout<<"last sta: "<<ista<< " prec: "<<1/param.GetQp() <<"x,y,z,tx,ty:
427 // "<<param.GetX()<<","<<param.GetY()<<","<<param.GetZ()<<","<<param.GetTx()<<","<<param.GetTy()<<endl;
428 trackC.SetParamLast(param);
429 kftr.Fit(kFALSE); // upstream
430 kftr.GetTrackParam(param);
431 trackC.SetParamFirst(param);
432 // cout<<"first sta: "<<ista<< " prec: "<<1/param.GetQp() <<"x,y,z,tx,ty:
433 // "<<param.GetX()<<","<<param.GetY()<<","<<param.GetZ()<<","<<param.GetTx()<<","<<param.GetTy()<<endl;
434
435 // Get residuals
436 CbmKFTrack tr = CbmKFTrack(trackC);
437 tr.Extrapolate(pos[2]);
438 ista = hit->GetStationNr();
439 tr.GetTrackParam(param);
440 // cout<<"Extrap sta: "<<ista<< " prec: "<<1/param.GetQp() <<"x,y,z,tx,ty:
441 // "<<param.GetX()<<","<<param.GetY()<<","<<param.GetZ()<<","<<param.GetTx()<<","<<param.GetTy()<<endl;
442 xext = param.GetX();
443 yext = param.GetY();
444 zhit = pos.Z();
445 // cout<<"pos z: "<<pos[2]<<" zhit: "<<zhit<<endl;
446 dXext = pos[0] - xext;
447 dYext = pos[1] - yext;
448 CbmKFTrack tr1 = CbmKFTrack(trackC);
449 fitter.SetKFHits(tr1, &trackC);
450 tr1.Smooth(pos[2]);
451 tr1.GetTrackParam(param);
452 // cout<<"SMO sta: "<<ista<< " prec: "<<1/param.GetQp() <<"x,y,z,tx,ty:
453 // "<<param.GetX()<<","<<param.GetY()<<","<<param.GetZ()<<","<<param.GetTx()<<","<<param.GetTy()<<endl;
454 // if(Niter==0) tr1.GetTrackParam(pIter1);
455 xsmo = param.GetX();
456 ysmo = param.GetY();
457
458 // cout<<"xsmo: "<< xsmo<<" ysmo: "<<ysmo<<endl;
459 dXsmo = pos[0] - xsmo;
460 dYsmo = pos[1] - ysmo;
461 nsec = hit->GetSectorNr();
462 /* if (ista == 1) {
463 nsec = hit->GetSectorNr();
464 } else
465 nsec = -100;
466
467 if (ista > 1) {
468 // if(hit->GetSectorNr()<=2) nsecGEM=1;
469 // else nsecGEM=2;
470
471 nsecGEM = hit->GetSectorNr();
472
473 // if(nsecGEM%2==0) continue;
474 } */
475 // cout<<"dxsmo: "<< dXsmo<<" dysmo: "<<dYsmo*10<<" posX: "<<pos[0]<<" posY: "<<pos[1]<< "xsmo:
476 // "<<xsmo<< " ysmo: "<<ysmo<<" Z: "<<pos[2]<<endl; cout<<"Z: "<<pos[2]<<endl;
477
478 trkAng = TMath::ATan(param.GetTx());
479 // stsCluster->GetQtot();
480
481 // Double_t adcHit=((CbmStsCluster)(stsCluster->UncheckedAt(hit->fDigiF)))->GetQtot();
482 CbmStsCluster* stsclF = (CbmStsCluster*)stsCluster->UncheckedAt(hit->fDigiF);
483 CbmStsCluster* stsclB = (CbmStsCluster*)stsCluster->UncheckedAt(hit->fDigiB);
484 Qf = stsclF->GetQtot();
485 Qb = stsclB->GetQtot();
486 // cout<<" adcHIT: "<<stscl->GetQtot()<<endl;
487
488 // if(Niter==1) txyz->Fill();
489 txyz->Fill();
490
491 // hits->Draw("dXsmo:xsmo","ista==6 && abs(dXsmo)<1 && nhits>4 && xsmo>-100 && prec<100","colz")
492 /*if(nhits>4 && abs(xsmo)>-100 && prec<100){
493 if(ista==1 && TMath::Abs(dXsmo)<1 ) h1->Fill(xsmo,dXsmo);
494 else if(ista==2 && TMath::Abs(dXsmo)<1 ) h2->Fill(xsmo,dXsmo);
495 else if(ista==3 && TMath::Abs(dXsmo)<1 ) h3->Fill(xsmo,dXsmo);
496 else if(ista==4 && TMath::Abs(dXsmo)<1 ) h4->Fill(xsmo,dXsmo);
497 else if(ista==5 && TMath::Abs(dXsmo)<1 )h5->Fill(xsmo,dXsmo);
498 else if(ista==6 && TMath::Abs(dXsmo)<1 ) h6->Fill(xsmo,dXsmo);
499
500 }
501 */
502
503 } // for (Int_t iHit = 0; iHit < nhits;
504
505 } // for (Int_t i = 0; i < nTracks;
506
507 } // for (Int_t iev = 0; iev < events;
508 // cout<<" size: "<<effstaUp.size()<<endl;
509 // for (size_t nstaeff = 0; nstaeff < effstaUp.size(); nstaeff++)
510 // cout << " sta: " << nstaeff << " eff: " << (Double_t)effstaUp[nstaeff] / (Double_t)effstaDown[nstaeff] <<
511 // endl;
512 // TCanvas* cHist2 = new TCanvas("cHist2","Canvas 2",600,600,800,800);
513 // cHist2->Divide(2);
514
515 // cHist2->cd(1); /*gPad->SetLogz();*/ gPad->SetGrid(); h1->Draw("colz");
516 // cHist2->cd(2); /*gPad->SetLogz();*/gPad->SetGrid(); h2->Draw("colz");
517 // cHist2->cd(1); /*gPad->SetLogz();*/gPad->SetGrid(); h3->Draw("colz");
518 // cHist2->cd(2); /*gPad->SetLogz();*/gPad->SetGrid(); h4->Draw("colz");
519 // cHist2->cd(3); /*gPad->SetLogz();*/ gPad->SetGrid();h5->Draw("colz");
520 // cHist2->cd(4); /*gPad->SetLogz();*/gPad->SetGrid(); h6->Draw("colz");
521
522 // cHist2->Update();
523 // out.Add(cHist2);
524
525 // out.Add(cHist3);
526 // Hlist->Write();
527 // Hlist->Add(txyz);
528 // Hlist->Add(hist2);
529 // cHist2->cd(1); hist1->Draw();
530 // cHist2->cd(2); hist2->Draw();
531 // cHist2->Update();
532 // Hlist->Write();
533
534 out.Write();
535 out.Close();
536}
537// ----------------------------------------------------------------------------
538
539Bool_t GetBmnGeom(const char* fileName)
540{
541 // ---> Get TGeoManager and top node (TOP)
542
543 TFile* geoFile = new TFile(fileName, "READ");
544 if (!geoFile->IsOpen()) {
545 cout << "-E- Could not open file!" << endl;
546 return kFALSE;
547 }
548 geoFile->FindObjectAny("SIGEMS"); //"GEMS_geom SIGEMS"); //AZ
549 return kTRUE; // AZ
550
551 // ---> Get TOP node from file
552 TList* keyList = geoFile->GetListOfKeys();
553 TIter next(keyList);
554 TKey* key = NULL;
555 TGeoVolume* top = NULL;
556 while ((key = (TKey*)next())) {
557 TString className(key->GetClassName());
558 if (className.BeginsWith("TGeoVolume")) {
559 top = dynamic_cast<TGeoVolume*>(key->ReadObj());
560 std::cout << "Found volume " << top->GetName() << endl;
561 break;
562 }
563 }
564 if (!top) {
565 cout << "-E- Could not find volume object in file" << endl;
566 return kFALSE;
567 }
568 // cout << gGeoManager << endl;
569 gGeoManager->GetListOfVolumes()->ls();
570 gGeoManager->SetTopVolume(top);
571
572 return kTRUE;
573}
574
576{
577
578 // Apply alignment
579
580 const int ip_r7[6][2] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11};
581
582 TGeoNode* st_node = NULL;
583 // Get STS node
584 TGeoNode* sts = NULL;
585 gGeoManager->CdTop();
586 TGeoNode* cave = gGeoManager->GetCurrentNode();
587 for (Int_t iNode = 0; iNode < cave->GetNdaughters(); iNode++) {
588 TGeoNode* node = cave->GetDaughter(iNode);
589 TString name = node->GetName();
590 if (name.Contains("GEMS_0")) {
591 sts = node;
592 gGeoManager->CdDown(iNode);
593 break;
594 }
595 }
596 if (!sts) {
597 cout << "-E- CbmBmnStsDigitize::ApplyAlignment: Cannot find top GEM node" << endl;
598 // return kFALSE;
599 }
600
601 // Loop over stations in STS
602 Int_t statNr = 0;
603 Int_t nStat0 = sts->GetNdaughters();
604 Int_t nMod0sts = 0;
605
606 cout << "STATIONS : " << nStat0 << endl;
607 /*
608 Double_t xAlign[26] =
609 {3.32829,-2.67635,2.83055,-3.16945, //I Si
610 2.84975, -3.11505, //II Si
611 8.91718,2.91888,-3.09552,-9.09292, 8.91618,2.91618,-3.08382,-9.08382, //III Si
612 -40.716, 40.914, 41.085, -40.598, 42.278, -39.429, -38.45, 43.25, 44.076, -37.596, -36.8092, 44.8207};//GEM
613 //{-40.8258, 40.8258, 40.845, -40.845, 40.8608, -40.8608, -40.8578, 40.8578, 40.845, -40.845, -40.8298, 40.8298};
614 //for 4 GeV
615
616 Double_t yAlign[26] =
617 {4.46297,5.07597, -1.55503,-1.05503, //I Si
618 1.8115,1.42, //II Si
619 7.84856,7.80306,7.92256,7.85506,-4.25594,-4.25594, -4.25594,-4.25594, // III Si
620 22.5211, 22.2759, 22.574, 22.6124, 22.8827, 22.516, 22.5385, 22.9377, 22.5828, 22.3481, 22.7875, 23.007};//GEM
621 //{22.2753, 22.4706, 22.4851, 22.6689, 22.7851, 22.5916, 22.3308, 23.1519, 22.4786, 22.4484, 22.5971, 23.2184};
622 //for 4 GeV
623
624 Double_t zAlign3[9] =
625 {14.365-0.025, 17.69-0.05, 26.0775, //3 planes of Si
626 42.2,67.77,115.569,137.8325,163.54,186.103}; //6 planes of GEM
627
628 Double_t zAlign23[26] =
629 {-0.315,0.315,0.365,-0.365, // I Si
630 -0.33,0.33, //II Si
631 0.8325,-0.7675, 0.7525, -0.8175, -0.4275, 1.1825, -0.2875, 1.3225, //III Si
632 0,0, 0,0, -0.085, 0.085, -0.0175,0.0175, 0.01,-0.01, 0,0}; //GEMS sensors
633 //{ 0.14595, -0.14595, -0.1258, 0.1258, -0.2005, 0.2005, 0.0945, -0.0945, -0.1125, 0.1125, 0.0865, -0.0865 };
634 //for 4 GeV
635 */
636 Double_t xAlign[26] = {3.33029, -2.67835, 2.82585, -3.14615, // I Si
637 2.85175, -3.11705, // II Si
638 8.91918, 2.92088, -3.09752, -9.09492, 8.91198, 2.94528, -2.98182, -8.98402, // III Si
639 -40.718, 40.916, 41.087, -40.6, 42.28, -39.431, -38.452, 43.252,
640 44.078, -37.598, -36.8112, 44.8227}; // GEM
641
642 Double_t yAlign[26] = {
643 4.46297, 5.07597, -1.48803, -0.910531, // I Si
644 1.8115, 1.42, // II Si
645 7.84856, 7.80306, 7.92256, 7.85506, -4.57154, -4.61864, -4.43934, -4.51954, // III Si
646 // 22.5211-0.17, 22.2759-0.17, 22.574-0.15, 22.6124-0.15, 22.8827-0.1, 22.516-0.1, 22.5385-0.05,
647 // 22.9377-0.05, 22.5828+0.025, 22.3481+0.025, 22.7875+0.07, 23.007+0.07};//GEM
648 22.5211, 22.2759, 22.574, 22.6124, 22.8827, 22.516, 22.5385, 22.9377, 22.5828, 22.3481, 22.7875,
649 23.007}; // GEM
650
651 Double_t zAlign3[9] = {14.352475, 17.68975, 26.2729875, // 3 planes of Si
652 42.19995, 67.76995, 115.569, 137.8325, 163.54, 186.103}; // 6 planes of GEM
653
654 Double_t zAlign23[26] = {
655 -0.406375, 0.430925, 0.274025, -0.298575, // I Si
656 -0.43165, 0.43165, // II Si
657 0.540813, -1.06009, 0.653113, -0.915988, -0.719887, 0.761012, -0.386287, 1.12281, // III Si
658 0.08685, -0.08685, -0.07145, 0.07145, -0.128, 0.128, 0.0115, -0.0115,
659 -0.004, 0.004, 0, 0}; // GEMS sensors
660 Double_t driftcorr = 0; //-0.2;
661 Int_t stn = 0, SensCh = 0;
662 // for (Int_t iNode = 1; iNode < sts->GetNdaughters(); iNode++) {
663 for (Int_t iNode = 0; iNode < nStat0; iNode++) {
664 // Go to station node
665 gGeoManager->CdDown(iNode);
666 TGeoNode* stationNode = gGeoManager->GetCurrentNode();
667 TString statName = stationNode->GetName();
668 if (!statName.Contains("station")) {
669 gGeoManager->CdUp();
670 continue;
671 }
672 // AZ Int_t statNr = stationNode->GetNumber();
673 ++stn; // AZ
674 /*
675 if(stn<=3){
676 gGeoManager->CdUp();
677 continue;}
678
679 if(stn<=3) continue;
680 */
681 ++statNr;
682
683 TGeoHMatrix* matr = gGeoManager->GetCurrentMatrix();
684 Double_t* statTrans = matr->GetTranslation();
685 // Double_t* statRot = matr->GetRotationMatrix();
686 TGeoRotation* r2;
687 // if(statNr-1 >1) statTrans[2] = zAlign3[statNr-1]+driftcorr;//zAlign[statNr];
688 cout << "statNr : " << statNr << ", " << statTrans[2] << " zal: " << zAlign3[statNr - 1] << endl;
689 statTrans[2] = zAlign3[statNr - 1]; // zAlign[statNr];
690
691 TGeoTranslation* t2 = new TGeoTranslation(statTrans[0], statTrans[1], statTrans[2]);
692 // if( statNr==2 || statNr==3 || statNr==5)r2 = new TGeoRotation("rot",0,180,180);
693 // if(statNr==1 || statNr==4 || statNr==6) r2 = new TGeoRotation("rot",0,180,180);
694 r2 = new TGeoRotation("rot", 0, 0, 0);
695
696 TGeoCombiTrans* cc1 = new TGeoCombiTrans(*t2, *r2);
697 // statTrans[0] = xAlign[statNr];
698 // statTrans[1] = yAlign[statNr];
699
700 // cout<<"ST TRANS Z: "<<statTrans[2]<<endl;
701
702 // if(statNr-1 >1) statTrans[2]=statTrans[2]+driftcorr;
703 // cout<<"ST TRANS Z: "<<statTrans[2]<<endl;
704 // if(iNode>0) statTrans[2]=statTrans[2]-2.47;
705 // matr->SetTranslation(statTrans);
706 // matr->SetMatrix(statRot);
707
708 TGeoHMatrix* matr0 = new TGeoHMatrix(*cc1);
709 matr0->RegisterYourself();
710
711 // Int_t nModules = stationNode->GetNdaughters();
712 Int_t nModules = GetNofModules(stationNode);
713 // cout<<"nModules: "<<nModules<<endl;
714
715 // sts->GetVolume()->ReplaceNode(stationNode,0,gGeoManager->GetCurrentMatrix()); //AZ
716 // sts->GetVolume()->RemoveNode(stationNode); //AZ
717 sts->GetVolume()->AddNode((TGeoVolumeAssembly*)stationNode->GetVolume(), 0, matr0); // AZ
718
719 // AZ- hard-coded st_node=sts->GetVolume()->GetNode(iNode+6);
720 st_node = (TGeoNode*)sts->GetVolume()->GetNodes()->Last(); // AZ
721
722 // Double_t statZ = statTrans[2];
723 // cout <<"sta: " << statNr << " " << gGeoManager->GetCurrentMatrix()->GetTranslation()[2] << " " <<
724 // sts->GetNdaughters() << endl;
725
726 // gGeoManager->CdUp(); // to sts
727
728 //-----------------------module translate-------------
729 Int_t moduleNr = 0, copy_no = 0;
730 // cout<<"nMODULES: "<< nModules<<endl;
731 /*if(iNode<2){
732 // nMod0sts=nModules;
733 nModules=nModules*2;
734
735 } */
736 // ---> Large sensors
737 for (Int_t iStatD = 0; iStatD < nModules; iStatD++) {
738 gGeoManager->CdDown(iStatD);
739
740 TGeoNode* module = gGeoManager->GetCurrentNode();
741 if (!TString(module->GetName()).Contains("module")) {
742 gGeoManager->CdUp();
743 continue;
744 }
745 cout << iStatD << " " << module->GetName() << " stn: " << statNr << endl;
746
747 //-------------------------add I II III Si----------------------------
748 if (iNode <= 2) {
749
750 Int_t nOfModSi = SensCh; // 10+(moduleNr++);
751 cout << "modSI: " << nOfModSi << " " << zAlign23[nOfModSi] << endl;
752 // fprintf(parFile, "%4d %4d\n", moduleNr, 1);
753 // ---> Sensor number
754 // geoMan->CdUp(); // back to module sensor
755 // SaveSensor(geoMan, parFile, phiStat, sensor);
756 TGeoHMatrix* matrMod = gGeoManager->GetCurrentMatrix();
757 Double_t* modTrans = matrMod->GetTranslation();
758 // cout<<"trans: "<<modTrans[2]<<endl;
759 modTrans[0] = xAlign[nOfModSi];
760 modTrans[1] = yAlign[nOfModSi];
761 modTrans[2] = zAlign23[nOfModSi]; // zAlign2[ipNr];
762 matrMod->SetTranslation(modTrans);
763 TGeoHMatrix* matr0Mod = new TGeoHMatrix(*matrMod);
764 matr0Mod->RegisterYourself();
765
766 // sts->GetVolume()->ReplaceNode(stationNode,0,gGeoManager->GetCurrentMatrix()); //AZ
767 // sts->GetVolume()->RemoveNode(stationNode); //AZ
768 // stationNode->GetVolume()->AddNode((TGeoVolumeAssembly*)module->GetVolume(),0,matr0Mod); //AZ
769
770 // cout<<" 1 st name add: "<< stationNode->GetName()<<" mod name add: "<<module->GetName()<< " "<<
771 // module<< " i: "<<iStatD<<" cols: "<<st_node->GetVolume()->GetNdaughters()<<endl; Double_t* sensTrans
772 // = matrMod->GetTranslation(); cout<<"trans: "<<sensTrans[0]<<" "<<sensTrans[1]<<" "<<sensTrans[2]<< "
773 // Nr mod: "<<moduleNr<<endl;
774 // stationNode->GetVolume()->AddNode((TGeoVolumeAssembly*)module->GetVolume(),copy_no,matr0Mod);
775 st_node->GetVolume()->AddNode((TGeoVolumeAssembly*)module->GetVolume(), copy_no, matr0Mod);
776 // cout<<" 1 st name add: "<< stationNode->GetName()<<" mod name add: "<<module->GetName() << " "<<
777 // module<< " i: "<<iStatD<<" cols: "<<st_node->GetVolume()->GetNdaughters()<<endl; Double_t modZ =
778 // modTrans[2];
779 Double_t* modTrans1 = gGeoManager->GetCurrentMatrix()->GetTranslation();
780 // cout <<"mod: " << nOfModSi << " VEC: " << modTrans1[0] << " "<<modTrans1[1] << " "<<modTrans1[2] <<
781 // endl;
782 copy_no++;
783 // delete matr0Mod;
784 gGeoManager->CdUp(); // back to module
785
786 SensCh++;
787 }
788 //-------------------------add I II III Si----------------------------
789
790 if (iNode > 2) {
791 if (TString(module->GetName()).Contains("Senso")) {
792 // if(iNode>0 && iNode<=3) moduleNr=0;
793 // else{
794 if (iStatD == 0)
795 moduleNr = 0;
796 if (iStatD > 0)
797 moduleNr = 1;
798 // }
799 // moduleNr=iStatD;
800 // fprintf(parFile, "%4d %4d\n", moduleNr, 1);
801 // SaveSensor(geoMan, parFile, phiStat, module);
802
803 Int_t ipNr = 0;
804 // cout<<" iNode: "<< iNode<<endl;
805 // if (iNode>0 && iNode<=3) ipNr = iNode;
806 // else ipNr = iNode*2 - 4 + moduleNr;
807 // Int_t order=-1;
808 // if( (statNr==1 || statNr==4 || statNr==6) && moduleNr==0) order=1;
809
810 TGeoHMatrix* matrMod = gGeoManager->GetCurrentMatrix();
811 Double_t* modTrans = matrMod->GetTranslation();
812 // cout<<" modtransX: "<<modTrans[0]<<" xali: "<<xAlign[ip_r7[statNr-1][moduleNr]]<< "sub:
813 // "<<xAlign[ip_r7[statNr-1][moduleNr]]-modTrans[0] <<endl;
814 cout << " aaaa ^^^^^^^^^^^^^^^^^^^ : " << statNr - 1;
815 cout << "\n" << moduleNr;
816 cout << "\n" << ip_r7[(statNr - 1) - 3][moduleNr] << endl;
817
818 modTrans[0] = xAlign[ip_r7[(statNr - 1) - 3][moduleNr] + 14];
819 modTrans[1] = yAlign[ip_r7[(statNr - 1) - 3][moduleNr]
820 + 14]; // - modTrans[1]; //not need anymore (- mod...)
821 modTrans[2] = zAlign23[ip_r7[(statNr - 1) - 3][moduleNr] + 14]; // zAlign2[ipNr];
822 // cout<<"ST TRANS Z: "<<zAlign23[ipNr]<<endl;
823 matrMod->SetTranslation(modTrans);
824 TGeoHMatrix* matr0Mod = new TGeoHMatrix(*matrMod);
825 matr0Mod->RegisterYourself();
826 ipNr++;
827 // cout<<" ip N: "<< ipNr<<endl;
828 // sts->GetVolume()->ReplaceNode(stationNode,0,gGeoManager->GetCurrentMatrix()); //AZ
829 // sts->GetVolume()->RemoveNode(stationNode); //AZ
830 // stationNode->GetVolume()->AddNode((TGeoVolumeAssembly*)module->GetVolume(),0,matr0Mod); //AZ
831
832 // cout<<" 1 st name add: "<< stationNode->GetName()<<" mod name add: "<<module->GetName()<< " "<<
833 // module<< " i: "<<iStatD<<" cols: "<<st_node->GetVolume()->GetNdaughters()<<endl; Double_t*
834 // sensTrans = matrMod->GetTranslation(); cout<<"trans: "<<sensTrans[0]<<" "<<sensTrans[1]<<"
835 // "<<sensTrans[2]<< " Nr mod: "<<moduleNr<<endl;
836 // stationNode->GetVolume()->AddNode((TGeoVolumeAssembly*)module->GetVolume(),copy_no,matr0Mod);
837 st_node->GetVolume()->AddNode((TGeoVolumeAssembly*)module->GetVolume(), copy_no, matr0Mod);
838 // cout<<" 1 st name add: "<< stationNode->GetName()<<" mod name add: "<<module->GetName() << " "<<
839 // module<< " i: "<<iStatD<<" cols: "<<st_node->GetVolume()->GetNdaughters()<<endl; Double_t modZ =
840 // modTrans[2]; cout <<"mod: " << ipNr << " VEC: " << modTrans[0] << " "<<modTrans[1] << "
841 // "<<modTrans[2] << endl;
842 copy_no++;
843 // delete matr0Mod;
844 }
845
846 gGeoManager->CdUp(); // back to module
847 }
848 } // for (Int_t iStatD = 0; iStatD < nModules;
849 //----------------------end module translate----------
850 // delete matr0;
851 gGeoManager->CdUp(); // to sts
852 } // station loop
853
854 Int_t snr = 1;
855 vector<TGeoNode*> removNodes;
856 // Remove extra nodes
857 for (Int_t iNode = 0; iNode < nStat0; iNode++) {
858 // Go to station node
859 // gGeoManager->CdDown(1);
860 gGeoManager->CdDown(iNode);
861 TGeoNode* stationNode = gGeoManager->GetCurrentNode();
862 cout << " st name del: " << stationNode->GetName() << " iNode: " << iNode << endl;
863 /*if(iNode>2)*/ removNodes.push_back(stationNode); // sts->GetVolume()->RemoveNode(stationNode); //AZ
864 gGeoManager->CdUp(); // to sts
865 }
866
867 for (Int_t o = 0; o < removNodes.size(); o++)
868 sts->GetVolume()->RemoveNode(removNodes[o]);
869 removNodes.clear();
870
871 for (Int_t iNode = 0; iNode < nStat0; iNode++) {
872 // Go to station node
873 // gGeoManager->CdDown(iNode);
874 gGeoManager->CdDown(iNode);
875 TGeoNode* stationNode = gGeoManager->GetCurrentNode();
876
877 TString statName = stationNode->GetName();
878 if (!statName.Contains("station")) {
879 gGeoManager->CdUp();
880 continue;
881 }
882 if (iNode > 2) {
883 Int_t nMod = 2;
884 for (Int_t iStatD = 0; iStatD < nMod; iStatD++) {
885 gGeoManager->CdDown(0); // stationNode->GetNdaughters()-1);
886
887 TGeoNode* module = gGeoManager->GetCurrentNode();
888 // cout<<" 2 st name del: "<< stationNode->GetName()<<" mod name del: "<<module->GetName() << " i:
889 // "<<iStatD<<endl;
890 stationNode->GetVolume()->RemoveNode(module); // AZ
891
892 gGeoManager->CdUp(); // to sts
893 }
894 } else { // cout<<" nModules: "<<GetNofModules(stationNode)<<endl;
895 Int_t currNodes = GetNofModules(stationNode) / 2;
896 for (Int_t iStatD = 0; iStatD < currNodes; iStatD++) {
897 gGeoManager->CdDown(0); // stationNode->GetNdaughters()-1);
898 TString mName = gGeoManager->GetCurrentNode()->GetName();
899 if (!mName.Contains("module")) {
900 gGeoManager->CdUp();
901 continue;
902 }
903 TGeoNode* moduleS = gGeoManager->GetCurrentNode();
904 // cout<<" DELETE VALUME: "<< stationNode->GetName()<<" mod name del: "<<moduleS->GetName() << " i:
905 // "<<iStatD<<endl;
906 stationNode->GetVolume()->RemoveNode(moduleS); // AZ
907
908 gGeoManager->CdUp(); // to sts
909 }
910 }
911
912 gGeoManager->CdUp();
913 }
914
915 for (Int_t iNode = 0; iNode < sts->GetNdaughters(); iNode++) {
916
917 // Go to station node
918 // gGeoManager->CdDown(iNode);
919 gGeoManager->CdDown(iNode);
920 TGeoNode* stationNode = gGeoManager->GetCurrentNode();
921 cout << " Check-in STATION: " << stationNode->GetName()
922 << " zpos : " << gGeoManager->GetCurrentMatrix()->GetTranslation()[2] << endl;
923 for (Int_t iModule = 0; iModule < stationNode->GetNdaughters(); iModule++) {
924 gGeoManager->CdDown(iModule);
925 TGeoNode* moduleNode = gGeoManager->GetCurrentNode();
926
927 TString gName = moduleNode->GetName();
928 if (gName.Contains("frame")) {
929 stationNode->GetVolume()->RemoveNode(moduleNode); // remove frames
930 iModule--;
931 gGeoManager->CdUp();
932 continue;
933 }
934
935 cout << " Check-in st name : " << stationNode->GetName() << " mod name : " << moduleNode->GetName() << endl;
936 Double_t* sensTrans = gGeoManager->GetCurrentMatrix()->GetTranslation();
937 cout << "trans mod: " << sensTrans[0] << " " << sensTrans[1] << " " << sensTrans[2] << endl;
938 for (Int_t iModule2 = 0; iModule2 < moduleNode->GetNdaughters(); iModule2++) {
939 gGeoManager->CdDown(iModule2);
940 TGeoNode* moduleNode2 = gGeoManager->GetCurrentNode();
941
942 TString hName = moduleNode2->GetName();
943 if (hName.Contains("frame")) {
944 moduleNode->GetVolume()->RemoveNode(moduleNode2); // remove frames
945 iModule2--;
946 gGeoManager->CdUp();
947 continue;
948 }
949
950 cout << " Check-in st name : " << stationNode->GetName() << " hot name : " << moduleNode2->GetName()
951 << endl;
952 Double_t* sensTrans = gGeoManager->GetCurrentMatrix()->GetTranslation();
953 cout << "trans hot: " << sensTrans[0] << " " << sensTrans[1] << " " << sensTrans[2] << endl;
954
955 gGeoManager->CdUp();
956 }
957 // stationNode->GetVolume()->RemoveNode(moduleNode); //AZ
958 gGeoManager->CdUp();
959 }
960
961 gGeoManager->CdUp(); // to sts
962 }
963
964 // exit(0);
965 // gGeoManager->SetVisLevel(500);
966 // gGeoManager->GetTopVolume()->SetTransparency(0);
967 // gGeoManager->GetTopVolume()->Draw("ogl");
968 // gGeoManager->Export("SIGEMS_r7_eventDisplay.root","SIGEMS"); //with geomanager
969
970 // Int_t c;
971 // cin>>c;
972}
int i
Definition P4_F32vec4.h:22
bool Bool_t
int Int_t
float Float_t
constexpr Bool_t kFALSE
double Double_t
constexpr Bool_t kTRUE
void Fatal(const char *location, const char *msgfmt,...)
R__EXTERN TGeoManager * gGeoManager
R__EXTERN TSystem * gSystem
void SetScale(Double_t factor)
Set a global field scaling factor.
Definition BmnFieldMap.h:55
virtual void Init()
Double_t GetScale() const
Accessor to global scaling factor.
Definition BmnFieldMap.h:90
BmnFileSource.h.
Field Map does not have any symmetries w.r.t coordinate axises.
Int_t Extrapolate(Double_t z, Double_t *QP0=0, Bool_t line=false)
Access to i-th hit.
Int_t Fit(Bool_t downstream=1, Bool_t line=false)
void GetTrackParam(FairTrackParam &track)
The CbmKF Kalman Filter algorithm.
Definition CbmKF.h:28
CbmStsCluster.h.
Double_t GetQtot() const
virtual Bool_t init(FairParIo *input)
Initialisation from input device.
Bool_t Init()
Initialisation.
virtual Int_t GetStationNr() const
Definition CbmStsHit.h:67
Int_t fDigiB
Definition CbmStsHit.h:102
Int_t GetSectorNr() const
Definition CbmStsHit.h:68
Int_t fDigiF
Definition CbmStsHit.h:101
void SetKFHits(CbmKFTrack &T, CbmStsTrack *track)
void Init()
Initialisation.
FairTrackParam * GetParamLast()
Definition CbmStsTrack.h:70
void SetParamLast(FairTrackParam &par)
Definition CbmStsTrack.h:84
TArrayI * GetStsHits()
Definition CbmStsTrack.h:71
void SetParamFirst(FairTrackParam &par)
Definition CbmStsTrack.h:83
Int_t GetNStsHits() const
Accessors.
Definition CbmStsTrack.h:60
Double_t GetChi2() const
Definition CbmStsTrack.h:66
Int_t GetStsHitIndex(Int_t iHit) const
Definition CbmStsTrack.h:62
FairTrackParam * GetParamFirst()
Definition CbmStsTrack.h:69
Int_t GetNDF() const
Definition CbmStsTrack.h:67
CbmVertex.h.
Definition CbmVertex.h:21
void Set(Int_t n) override
Long64_t GetEntries() const override
TList * GetListOfKeys() const override
TObject * FindObjectAny(const char *name) const override
Int_t Write(const char *name=nullptr, Int_t opt=0, Int_t bufsize=0) const override
virtual Bool_t IsOpen() const
void Close(Option_t *option="") override
virtual Double_t * GetTranslation()
void SetTranslation(const Double_t *vect)
virtual void RegisterYourself()
TGeoVolume * GetVolume() const
Int_t GetNdaughters() const
TGeoNode * GetDaughter(Int_t ind) const
void RemoveNode(TGeoNode *node)
virtual TGeoNode * AddNode(TGeoVolume *vol, Int_t copy_no, TGeoMatrix *mat=nullptr, Option_t *option="")
TObjArray * GetNodes()
void Add(TObject *obj) override
const char * GetName() const override
Int_t GetEntriesFast() const
TObject * Last() const override
TObject * UncheckedAt(Int_t i) const
const char * Data() const
Bool_t BeginsWith(const char *s, ECaseCompare cmp=kExact) const
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
virtual Int_t Fill()
virtual Int_t Branch(const char *folder, Int_t bufsize=32000, Int_t splitlevel=99)
Double_t Z() const
Int_t GetNofModules(TGeoNode *station)
void ApplyAlignment()
Definition ResidOk_r7.C:575
const Int_t ip_r7[6][2]
Definition ResidOk_r7.C:67
Bool_t GetBmnGeom(const char *fileName)
Definition ResidOk_r7.C:539
void ResidOk_r7(TString inFile="", TString histFile="", Double_t Mscale=1.366666667)
Definition ResidOk_r7.C:71
Double_t ATan(Double_t)