FairRoot/PandaRoot
createRootGeometry_DIRC_sepEV_06_2013.C
Go to the documentation of this file.
1 // macro that creates the Barrel DIRC geometry
2 // to be able to run this macro please be sure you have the following configuration in gconfig/g4Config.C file:
3 // TG4RunConfiguration* runConfiguration
4 // = new TG4RunConfiguration("geomRoot", "QGSP_BERT_EMV+optical", "stepLimiter+specialCuts+specialControls");
5 
6 // input parameters are:
7 // fFocusingSystem = 0 - no focusing is used
8 // fFocusingSystem = 1 - lens
9 // fFocusingSystem = 2 - forward mirror is used
10 // fFocusingSystem = 3 - new lens (no air gap) is used
11 // fprizm = kFALSE - no prism
12 // fprizm = kTRUE - prism
13 
14 // allowed configurations:
15 // fFocusingSystem = 0 && fprizm = kFALSE output file is in geometry/dirc_l0_p0.root
16 // fFocusingSystem = 0 && fprizm = kTRUE output file is in geometry/dirc_l0_p1.root
17 // fFocusingSystem = 1 && fprizm = kFALSE output file is in geometry/dirc_l1_p0.root
18 // fFocusingSystem = 2 && fprizm = kFALSE output file is in geometry/dirc_l2_p0.root
19 // fFocusingSystem = 2 && fprizm = kTRUE output file is in geometry/dirc_l2_p1.root
20 // corresponding .root files with geometry please fone in the geometry directory
21 
22 // B E F O R E R U N N I N G T H E P R O G R A M
23 // C H E C K T H E N A M E O F T H E O U T P U T F I L E ! ! !
24 const Double_t pi = 3.1415926535;
25 
26 void getXrange(double dphi, double radiusCornerOut, double radiusCornerIn, double MCPactiveArea, int cas1, int cas2, double ylow, double step, double &x1, double &x2);
27 double getX(TVector3 corner1, TVector3 corner2, double yy);
28 int findSectorIn(double y, double dphi, double radius, double EVdrop, double hthick);
29 int findSectorOut(double y, double dphi_rad, double radiusCornerOut);
30 
31 int createRootGeometry_DIRC_sepEV_06_2013(Int_t fFocusingSystem = 0, Bool_t fprizm = kFALSE){
32 
33  const Double_t pi = 3.1415926535;
34 
35  gROOT->Macro("$VMCWORKDIR/gconfig/rootlogon.C");
36 
37  TString vmcWorkdir = getenv("VMCWORKDIR");
38 
39  // Load this libraries
40  gSystem->Load("libGeoBase");
41  gSystem->Load("libParBase");
42  gSystem->Load("libBase");
43  gSystem->Load("libPndData");
44  gSystem->Load("libPassive");
45  gSystem->Load("libEve");
46 
47  // variable to achieve the DIRC basic parameters
48  PndGeoDrc* fGeo = new PndGeoDrc();
49 
50  // units = cm
51 
52  Double_t eps = 0.01; // epsilon
53  Double_t mirr_hthick = 0.02;
54  Double_t PDbaseLayer = 5.; // [cm] thickness of the carbon at the back of the EV
55 
56  Double_t radius = fGeo->radius(); // 50. radius in middle of the barbox (x and y)
57  Double_t hthick = fGeo->barHalfThick(); // 1.7/2. half thickness of the bars
58  Double_t barwidth = fGeo->BarWidth(); // 3.2 width of the radiator bar
59  Double_t barnum = fGeo->barNum(); // 6 number of bars per barbox
60  Double_t bbnum = fGeo->BBoxNum(); //16. total number of sides = barboxes
61  Double_t bbGapAngle = fGeo->BBoxGap(); //1.5 gap btw the neighboring barboxes (at the middle height)
62  Double_t pipehAngle = fGeo->PipehAngle(); //3.6 [degrees] half of the angular space needed for the target pipe
63  Double_t barWin_hthick = 0.1/2.; // [cm]=15um thickness of the 'glas' window at the readout end of the barbox;
64  Double_t phi0 = (180.-2.*pipehAngle)/bbnum + pipehAngle;
65  Double_t dphi = (180.-2.*pipehAngle)/bbnum*2.;
66  Double_t dphi_rad = dphi/180.*pi;
67 
68  Double_t bbox_zdown = fGeo->barBoxZDown(); // 130. bar box z downstream
69  Double_t bbox_zup = fGeo->barBoxZUp(); // -120. bar box z upstream
70  Double_t bbox_hlen = 0.5*(bbox_zdown - bbox_zup); // 125. bar box half length
71  Double_t bbox_shift = bbox_zup + bbox_hlen; // 5. bar box shift
72  Double_t barhgap = fGeo->barhGap(); // 0.01 half gap between bars
73  Double_t boxgap = fGeo->boxGap(); // 0.1 gap between bars and the bar box
74  Double_t boxthick = fGeo->boxThick(); // 0.05 thickness of the bar box
75  Double_t gluehthick = 0.0005;//[cm] // glue layer, connects two halfs into one long radiator bar
76  Double_t len = 0.; // length of the lenses block. see further
77  Double_t fSlabEnd = 0.; // [cm] position at which bar in front of EV ends
78  Double_t fdz_mirr1 = 0.;
79  Double_t fdz_mirr2 = 0.;
80  Double_t fdz_lens3 = 0.;
81  Double_t fdz_lens2 = 0.;
82  Double_t fdz_lens1 = 0.;
83 
84  //parameters for MCPs:
85  Double_t MCPsize = 5.78;//5.9; //[cm] width=height of an MCP
86  Double_t MCPactiveArea = 5.3; //[cm] width=height of the active area of an MCP
87  Double_t MCPgap = 0.01; //[cm] gap between MCPs
88  Double_t PixelSize = 0.65; //[cm] size of a pixel
89  Int_t Npix = 8; // number of mcps in a row/column
90  const Int_t NpixTotal = Npix*Npix; // total number of pixels for 1 mcp
91  Double_t PixelGap = (MCPactiveArea - Double_t(Npix)*PixelSize)/(Double_t(Npix)-1.);
92  Double_t PDwindowThick = 0.1; //[cm]
93  Double_t PhCathodeThick = 0.0001; //[cm] = 1 um
94  Double_t PDsensitiveThick = 0.01; //[cm]
95  Double_t PDgreaseLayer = 0.1; //[cm]
96  Double_t hgap = 0.5*(MCPsize - MCPactiveArea + MCPgap); // half gap btw active areas of MCPs
97  Double_t step = MCPactiveArea + 2.*hgap; // step in which to locate MCPs = 6 for Mcp1, Mcp2; = 5.88 for Mcp2a. This is total size of one MCP with gaps in between
98 
99  // EV parameters
100  Double_t sob_len = fGeo->EVlen(); // 30. in current version
101  Double_t sob_shift = -bbox_hlen + bbox_shift - sob_len; // -150.
102  Double_t sob_angleB = fGeo->EVbackAngle(); //90. [degrees] angle of the EV (usually it is 90)
103  Double_t EVdrop = 0.;//fGeo->EVdrop(); //0.7 [cm] drop of the EV - inner radius
104  Double_t EVoffset = fGeo->EVoffset(); //0. [cm] offset of the EV - outer radius
105  Double_t EVgreaseLayer = 0.0015; //[cm] grease layer thickness btw the bar window and the EV
106  //Double_t sob_angle = 45.;//fGeo->EVangle(); //80. [degrees] opening angle of the EV
107  Double_t sob_angle = atan((5.*step-2.*hthick-EVoffset-EVdrop)/sob_len)/pi*180.;// [degrees] opening angle of the EV
108  Double_t sob_Rout = (radius + hthick + EVoffset +(MCPsize-MCPactiveArea)+ 2.*boxgap+2.*boxthick +(sob_len + EVgreaseLayer + PDbaseLayer)*tan(sob_angle/180.*pi));// [cm] radius in the middle of the section
109 
110  Double_t bbSideGap = 0.5*( ((barwidth+2.*barhgap)*barnum) - (2.*radius*sin((dphi-bbGapAngle)/180.*pi/2.)+2.*barhgap) );
111  Double_t bbX = 2.*radius*sin((dphi-bbGapAngle)/180.*pi/2.)+2.*barhgap;
112 
113  Double_t radiusCornerIn = (radius-hthick-EVdrop)/cos(dphi/2./180.*pi);
114  Double_t radiusCornerOut = sob_Rout/cos(dphi_rad/2.);
115  Double_t radiusMiddleSmall = radius-hthick-EVdrop;
116 
117  cout<<"radius = "<<radius<<"rad corner in = "<< radiusCornerIn<<", rad corner out = "<< radiusCornerOut<<endl; cout<<"pixel gap = "<< PixelGap<<endl;
118  cout<<"sob_angle = "<<sob_angle<<endl;
119  cout<<"sob_Rout = "<<sob_Rout<<endl;
120  cout<<"bbSideGap = "<<bbSideGap<<", bbX+2.sidegap = "<<bbX+2.*bbSideGap<<", 5*barwidth+5*2*bargap = "<<5.*barwidth+5.*2.*barhgap<<endl;
121  cout<<", bbX = "<<bbX<<endl;
122  cout<<"dphi = "<<dphi<<", phi0 = "<<phi0<<endl;
123  //----------------------------------------------------------------------
124 
125  // Rotations:
126  TGeoRotation rot1;
127  rot1.RotateZ(90.);
128 
129  TString fGeoFile= Form("../../geometry/dirc_l%d_p%d_mirrorGap_sepEV_06_2013.root",fFocusingSystem, fprizm);
130  TFile* fi = new TFile(fGeoFile,"RECREATE");
131  cout<<"Output file = "<<fGeoFile<<endl;
132 
133  FairGeoLoader* geoLoad = new FairGeoLoader("TGeo","FairGeoLoader");
134  FairGeoInterface *geoFace = geoLoad->getGeoInterface();
135  cout << "geoface setmediafile" << endl;
136  geoFace->setMediaFile("../../geometry/media_pnd.geo");
137  //cout << "geoface readmedia" << endl;
138  geoFace->readMedia();
139  //cout << "geoface print" << endl;
140  //geoFace->print();
141  //cout << "geoface done" << endl;
142 
143  FairGeoMedia *Media = geoFace->getMedia();
144  FairGeoBuilder *geobuild=geoLoad->getGeoBuilder();
145 
146  FairGeoMedium *DrcAir = Media->getMedium("air");
147  FairGeoMedium *DrcAirNoSens = Media->getMedium("DIRCairNoSens");
148  FairGeoMedium *DrcEpotek301_2 = Media->getMedium("Epotek301_2");
149  FairGeoMedium *DrcOpticalGrease = Media->getMedium("OpticalGrease");
150  FairGeoMedium *DrcCarbonFiber = Media->getMedium("DIRCcarbonFiber");
151  FairGeoMedium *DrcFusedSil = Media->getMedium("FusedSil");
152  FairGeoMedium *DrcMirror = Media->getMedium("Mirror");
153  FairGeoMedium *DrcMarcol82_7 = Media->getMedium("Marcol82_7");
154  FairGeoMedium *DrcNLAK33A = Media->getMedium("NLAK33A");
155  FairGeoMedium *DrcPhotocathode= Media->getMedium("Photocathode");
156 
157  Int_t nmed=geobuild->createMedium(DrcAir);
158  nmed=geobuild->createMedium(DrcAirNoSens);
159  nmed=geobuild->createMedium(DrcEpotek301_2);
160  nmed=geobuild->createMedium(DrcOpticalGrease);
161  nmed=geobuild->createMedium(DrcCarbonFiber);
162  nmed=geobuild->createMedium(DrcFusedSil);
163  nmed=geobuild->createMedium(DrcMirror);
164  nmed=geobuild->createMedium(DrcMarcol82_7);
165  nmed=geobuild->createMedium(DrcNLAK33A);
166  nmed=geobuild->createMedium(DrcPhotocathode);
167  //-------------------------------------------------------------------------------------------------------------
168  // focusing systems:
169  // 0 no ficusing:
170  if(fFocusingSystem == 0){
171  Double_t flen = 0.;
172  fSlabEnd = -bbox_hlen + bbox_shift + flen;
173  cout<<"bar ends at = "<<fSlabEnd<<endl;
174 
175  //fAtBarEnd = "DrcBar";
176  }
177 // 1 lens
178  if (fFocusingSystem == 1){ // L E N S E S (no airgaps, thin nlak33)
179  // some notes to lens operations (revision 9649) is in
180  // ~carsten/work/documents/software/PandaRoot/lens_definitions_2.pdf
181 
182  Double_t r = 75.18; // first lens radius (cm)
183  Double_t alpha = TMath::ASin(barwidth/2./r); // lside -> bbX, assuming that there are 5 bars in the bar box??
184  Double_t a = r - r*TMath::Cos(alpha);
185  Double_t b = a + 0.5; // box dimension .6 instead of .5 due to strong curvature
186 
187  cout<<" DIRC a,b = "<<a<<" "<<b<<endl;
188 
189 
190  Double_t r2 = 75.18; // radius second lens (cm)
191  Double_t b2 = 0.2; // + a2;
192 
193  Double_t r3 = 17.95; // third lens radius (cm)
194  Double_t alpha3 = TMath::ASin(barwidth/2./r3); // lside -> bbX, assuming that there are 5 bars in the bar box??
195  Double_t a3 = r3 - r3*TMath::Cos(alpha3);
196  Double_t b3 = a3 + 0.0; // box dimension .6 instead of .5 due to strong curvature
197 
198  // lens1 (b) + lens2 (0.6) + lens3 (b3) + gap (0.5)
199 
200  len = b + 0.2 + b3 + 0.5;//+ a2; // dimension of the box containing both lenses
201 
202  cout<<"DIRC len= "<<len<<endl;
203 
204  fSlabEnd = -bbox_hlen + bbox_shift + len; // used in processHits
205  cout<<"bar ends at = "<<fSlabEnd<<endl;
206 
207  // Lenses
208 
209  // Lens 1
210  Double_t t = -r +b/2;
211 
212  TGeoSphere* logicSphere= new TGeoSphere("S",0.,r, 0. ,180.,0.,360.);
213  TGeoBBox* lBox = new TGeoBBox("B", barwidth/2., hthick, b/2.);
214  TGeoTranslation *tr1 = new TGeoTranslation("tr1", 0.,0., t);
215  tr1->RegisterYourself();
216  TGeoCompositeShape *cs = new TGeoCompositeShape("cs","S*(B:tr1)");
217  TGeoVolume *lens1 = new TGeoVolume("DrcLENS1Sensor",cs, gGeoManager->GetMedium("FusedSil"));
218  lens1->SetLineColor(kRed-8);
219  lens1->SetTransparency(40);
220 
221  // position lens within already shifted bar container at -(bbox_hlen-eps)+len, the lens base is -r + b
222  // with -0.01 one can make a gap visible (.1mm) for orientation
223  //barContainer->(lens1, 1,new TGeoCombiTrans(0., 0., -(bbox_hlen-eps)+len -(-r+b) /*-0.01*/, new TGeoRotation (0)));
224 /* barContainer->AddNode(lens1, 1,new TGeoCombiTrans(0., 0., -(bbox_hlen)+len -(-r+b), new TGeoRotation (0)));
225 */ fdz_lens1 = -(bbox_hlen)+len -(-r+b);
226 
227 
228  //Lens 2
229  Double_t t2 = -r2;// +b2/2 r2 is the reference point (concave lens)
230  TGeoSphere* logicSphere2 = new TGeoSphere("S2",0 ,r2, 0. ,180.,0.,360.);
231  TGeoBBox* lBox2 = new TGeoBBox("B2", barwidth/2., hthick, b2/2.); // lside -> bbX
232  TGeoTranslation *tr2 = new TGeoTranslation("tr2", 0.,0., t2);
233  tr2->RegisterYourself();
234  TGeoCompositeShape *cs2 = new TGeoCompositeShape("cs2","(B2:tr2)-S2");
235  TGeoVolume *lens2 = new TGeoVolume("DrcLENS2Sensor",cs2, gGeoManager->GetMedium("NLAK33A"));
236  lens2->SetLineColor(kRed+2);
237  lens2->SetTransparency(40);
238 
239  // place the lens exactly on lens1
240  // position lens within already shifted bar container at -(bbox_hlen-eps)+len, the lens base is -r2
241  // the tip of lens1 is at b
242  // with -0.02 one can make a gap visible (.1mm due to lens1) for orientation
243  //barContainer->AddNode(lens2, 1,new TGeoCombiTrans(0., 0., -(bbox_hlen-eps)+len -(-r2) -b /*-0.02*/, new TGeoRotation (0)));
244 /* barContainer->AddNode(lens2, 1,new TGeoCombiTrans(0., 0., -(bbox_hlen)+len -(-r2) -b , new TGeoRotation (0)));
245 */ fdz_lens2 = -(bbox_hlen)+len -(-r2) -b;
246 
247  //Lens3 (like lens1, same treatment)
248  Double_t t3 = -r3+b3/2;
249  TGeoSphere* logicSphere3= new TGeoSphere("S3",0.,r3, 0. ,180.,0.,360.);
250  TGeoBBox* lBox3 = new TGeoBBox("B3", barwidth/2., hthick, b3/2.); // lside -> bbX
251  TGeoTranslation *tr3 = new TGeoTranslation("tr3", 0.,0., t3);
252  tr3->RegisterYourself();
253  TGeoCompositeShape *cs3 = new TGeoCompositeShape("cs3","S3*(B3:tr3)");
254  TGeoVolume *lens3 = new TGeoVolume("DrcLENS3Sensor",cs3, gGeoManager->GetMedium("NLAK33A"));
255  lens3->SetLineColor(kRed-6);
256  lens3->SetTransparency(40);
257 
258  // place the lens exactly on lens2 plane side
259  // position lens within already shifted bar container at -(bbox_hlen-eps)+len, the lens base is -r3 + b3
260  // with -0.03 one can make a gap visible (.1mm due to lens 1&2) for orientation
261  // b2/2 is the thickness of lens2 in the middle
262  //barContainer->AddNode(lens3, 1,new TGeoCombiTrans(0., 0., -(bbox_hlen-eps)+len -(-r3+b3) -b - b2/2 /*-0.03*/, new TGeoRotation (0)));
263 /* barContainer->AddNode(lens3, 1,new TGeoCombiTrans(0., 0., -(bbox_hlen)+len -(-r3+b3) -b - b2/2 , new TGeoRotation (0)));
264  //cout<<" DIRC r,r2,r3 = "<<r<<" "<<r2<<" "<<r3<<endl;
265 */ fdz_lens3 = -(bbox_hlen)+len -(-r3+b3) -b - b2/2. ;
266 
267  //fAtBarEnd = "DrcLENS3";
268  } // E N D O F N E W L E N S E S
269 
270  // 2 - mirror:
271  if (fFocusingSystem == 2){ // Mirrors at front
272  // Put some mirrors at the downstream end with a focal plane at the PD.
273 
274  double zpos=120.;
275 
276  // The bar is produced with "bbox_hlen-fabs(len)/2.-mirr_hthick"
277  len = -130. + zpos + 1.;
278  //len = -247.0; // negative to make space at downstream end of bar.
279 
280  Double_t len1 = 1.5; // 1st block
281 
282  Double_t mirror_angle = 0.0;// 0.0 is pointing upstream, -90.0 is pointing to the beam axis
283  Double_t focal_length = 2. * bbox_hlen + sob_len - fabs(len) + len1;
284  Double_t mirror_radius = 2. * focal_length;
285 
286  cout<<" mirror radius: "<<mirror_radius<<endl;
287 
288 
289  // no angle for the moment
290  // block
291  TGeoSphere* logicSphere = new TGeoSphere("S",0.,mirror_radius, 0. ,180.,0.,360.);
292  TGeoBBox* lBox = new TGeoBBox("B", barwidth/2., hthick, fabs(len1)/2.); // lside -> bbX
293 
294  Double_t t = mirror_radius - len1/2.;
295 
296  TGeoTranslation *tr1 = new TGeoTranslation("tr1", 0.,0., t);
297  tr1->RegisterYourself();
298  TGeoCompositeShape *cs = new TGeoCompositeShape("cs","S*(B:tr1)");
299 
300  TGeoVolume *block1 = new TGeoVolume("DrcBlock1",cs, gGeoManager->GetMedium("FusedSil"));
301  block1->SetLineColor(kRed);
302  block1->SetTransparency(40);
303 
304 
305  Double_t shift1 = len1-mirror_radius; // Now the start of block1 is at zero
306  shift1 += bbox_hlen-fabs(len)-2.*mirr_hthick;
307 
308 /*
309  barContainer->AddNode(block1, 1, new TGeoCombiTrans(0., 0., shift1, new TGeoRotation (0)));
310 */ fdz_mirr1 = shift1;
311 
312 
313  Double_t gap = 0;
314  Double_t len2 = fabs(len) - len1 - gap;
315 
316  // no angle for the moment
317  // block
318  TGeoSphere* logicSphere2 = new TGeoSphere("S2",0,mirror_radius, 0. ,180.,0.,360.);
319  TGeoBBox* lBox2 = new TGeoBBox("B2", barwidth/2., hthick, fabs(len2)/2.); // lside -> bbX
320 
321  // make radius part of block
322  Double_t t2 = mirror_radius + len2/2 - 1;
323 
324  TGeoTranslation *tr2 = new TGeoTranslation("tr2", 0.,0., t2);
325  tr2->RegisterYourself();
326  TGeoCompositeShape *cs2 = new TGeoCompositeShape("cs2","(B2:tr2)-S2");
327 
328  TGeoVolume *block2 = new TGeoVolume("DrcBlock2",cs2, gGeoManager->GetMedium("Mirror"));
329  block2->SetLineColor(kGreen);
330  block2->SetTransparency(40);
331 
332 
333  Double_t shift2 = -mirror_radius; // Now the start of the block is at zero
334  shift2 += bbox_hlen-fabs(len)-2*mirr_hthick; // at end of bar
335  shift2 += len1 + gap;
336 
337  fdz_mirr2 = shift2;
338 
339 /*
340  barContainer->AddNode(block2,
341  1,
342  new TGeoCombiTrans(0.,
343  0.,
344  shift2,
345  // bbox_hlen-mirror_radius-2*mirr_hthick-len2+len1+1,
346  //bbox_hlen-mirror_radius-len2+0.1,
347  new TGeoRotation (0)
348  )
349  );
350 */
351 
352 
353  Double_t flen = 0.;
354  fSlabEnd = -bbox_hlen + bbox_shift + flen;
355  cout<<"bar ends at = "<<fSlabEnd<<endl;
356 
357  //fAtBarEnd = "DrcBar";
358  } // E N D O F MIRRORS
359 
360  if(fFocusingSystem == 3){ // cylindrical lens w/o airgap
361 
362  // main parameters:
363 
364  Double_t Hcyl2 = 0.1; //[cm] thickness in the middle of the second (FusedSil) lens
365  Double_t Rcyl = 7.35; // [cm]
366  Double_t Lcyl = barwidth/2.; //cylinder length/2.
367  Double_t Acyl = TMath::ASin(hthick/Rcyl); // angle
368  Double_t Tcyl1 = 0.5*Rcyl*(1. + TMath::Cos(Acyl)); // distance between centers of the cylinder and the box
369 
370  len = Hcyl2 + Rcyl*(1. - TMath::Cos(Acyl));
371  cout<<"DIRC: len = "<<len<<endl;
372 
373  cout<<"length = "<< Lcyl<<", angle = "<<Acyl/3.1415*180.<<", translation = "<<Tcyl1<<endl;
374  cout<<"lens 0.5*width = "<<hthick<<", 0.5*thickness = "<< 0.5*Rcyl*(1.-TMath::Cos(Acyl))<<", 0.5*heigth = "<<Lcyl <<endl;
375 
376  //Lens1
377  TGeoEltu* lCylinder = new TGeoEltu("Cyl",Rcyl,Rcyl,Lcyl);
378  TGeoBBox* lCylBox = new TGeoBBox("CylBox", hthick, hthick, Lcyl);
379  //TGeoBBox* lCylBox = new TGeoBBox("CylBox",36.75, 36.75, Lcyl);
380  //TGeoTranslation *trCyl = new TGeoTranslation("trCyl", 0., 100.4, 0.);
381  TGeoTranslation *trCyl1 = new TGeoTranslation("trCyl", 0., Rcyl*TMath::Cos(Acyl)+hthick, 0.);// !!!
382  trCyl1->RegisterYourself();
383  TGeoCompositeShape *llens1 = new TGeoCompositeShape("llens1","Cyl*(CylBox:trCyl)");
384  TGeoVolume *CylLens1 = new TGeoVolume("DrcLENS1Sensor",llens1, gGeoManager->GetMedium("NLAK33A"));
385  CylLens1->SetLineColor(kRed-8);
386  CylLens1->SetTransparency(40);
387 
388  fdz_lens1 = -(bbox_hlen - barWin_hthick) + len - (-Rcyl + Rcyl*(1.-TMath::Cos(Acyl)));
389 
390  //Lens2
391  TGeoTranslation *trCyl2 = new TGeoTranslation("trCyl2", 0., Rcyl + Hcyl2 - hthick, 0.);
392  trCyl2->RegisterYourself();
393  TGeoCompositeShape *llens2 = new TGeoCompositeShape("llens2", "(CylBox:trCyl2) - Cyl");
394  TGeoVolume* CylLens2 = new TGeoVolume("DrcLENS2Sensor", llens2, gGeoManager->GetMedium("FusedSil"));
395  CylLens2->SetLineColor(kRed+2);
396  CylLens2->SetTransparency(40);
397 
398  // TGeoRotation rot_lens3;
399  // rot_lens3.RotateZ(-90.);
400  // rot_lens3.RotateX( 90.);
401 
402 
403 
404  fdz_lens2 = -(bbox_hlen - barWin_hthick) + len + (hthick - Rcyl*(1.-TMath::Cos(Acyl)) - Hcyl2);
405  fSlabEnd = -bbox_hlen + bbox_shift;
406  cout<<"bar ends at = "<<fSlabEnd<<endl;
407 
408  }
409  //---------------------------------------------------------------------------------------------
410 
411 
412  // create top volume:
413  TGeoManager* gGeoMan = (TGeoManager*)gROOT->FindObject("FAIRGeom");
414 
415  TGeoBBox* lTop = new TGeoBBox(500,500,300);
416  TGeoVolume* top = new TGeoVolume("DIRC", lTop, gGeoManager->GetMedium("air"));
417  gGeoManager->SetTopVolume(top);
418 
419  // create pre-top volume:
420  TGeoVolume* vLocalMother;
421  TGeoPcon* shape = new TGeoPcon("BarrelDIRCShape", 0, 360., 4);
422  shape->DefineSection(0, bbox_zdown, 35., 60.);
423  shape->DefineSection(1, bbox_zup, 35., 60.);
424  shape->DefineSection(2, bbox_zup - sob_len, 35., sob_Rout+EVoffset+1.);
425  shape->DefineSection(3, bbox_zup - sob_len - PDbaseLayer - EVgreaseLayer, 35., sob_Rout+EVoffset+1.);
426  vLocalMother = new TGeoVolume("BarrelDIRC", shape, gGeoManager->GetMedium("DIRCairNoSens")); //DIRCairNoSens_m); //################
427  top->AddNode(vLocalMother, 0,0);
428 
429  cout<<"bbox length = "<<2.*(bbox_hlen-0.5*(boxgap+boxthick))<<endl;
430 
431  // create BarBoxes:
432  TGeoBBox* logicbbL;
433  TGeoVolume *bbox;
434  logicbbL = new TGeoBBox("logicbbL", 0.5*barnum*(barwidth+2.*barhgap)+bbSideGap+boxthick, hthick+boxgap+boxthick, bbox_hlen);
435 
436  TGeoTrap* logicPrizmBox = new TGeoTrap("logicPrizmBox",sob_len/2. + PDbaseLayer/2. + EVgreaseLayer/2.,
437  atan(tan(sob_angle/180.*pi)/2.)/pi*180., //2
438  270., //3
439  (2.*hthick+EVdrop+EVoffset+(MCPsize-MCPactiveArea)+ 2.*boxgap+2.*boxthick+(sob_len + PDbaseLayer + EVgreaseLayer)*tan(sob_angle*pi/180.))/2.,
440  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea))+ boxgap + boxthick, //4
441  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea))+ boxgap + boxthick, //4
442  0, //7
443  (2.*hthick+EVdrop+EVoffset+(MCPsize-MCPactiveArea) +2.*boxgap+2.*boxthick)/2.,
444  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea))+ boxgap + boxthick, //8
445  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea))+ boxgap + boxthick, //8
446  0); //11
447 
448 
449  TGeoTranslation* trprb = new TGeoTranslation("trprb", 0., 0.5*(logicPrizmBox->GetH1()+ logicPrizmBox->GetH2())-hthick-EVdrop-boxthick-boxgap-0.5*(MCPsize-MCPactiveArea), -bbox_hlen-sob_len/2.-PDbaseLayer/2.-EVgreaseLayer/2.);
450  trprb->RegisterYourself();
451  TGeoCompositeShape *cspb = new TGeoCompositeShape("cspb","logicbbL + logicPrizmBox:trprb");
452  //bbox = new TGeoVolume("DrcBarBox", cspb,gGeoManager->GetMedium("DIRCcarbonFiber"));
453  bbox = new TGeoVolume("DrcBarBox", logicbbL,gGeoManager->GetMedium("DIRCcarbonFiber"));
454  bbox->SetLineColor(30);
455 
456 
457  TGeoBBox* logicbbS;
458  TGeoVolume *abox;
459  logicbbS = new TGeoBBox("logicbbS", 0.5*barnum*(barwidth+2.*barhgap)+bbSideGap, hthick+boxgap, bbox_hlen-barWin_hthick);
460  abox = new TGeoVolume("DrcAirBox", logicbbS, gGeoManager->GetMedium("DIRCairNoSens"));
461  bbox->AddNode(abox, 1, new TGeoCombiTrans(0., 0., barWin_hthick, new TGeoRotation(0)));
462  abox->SetLineColor(19); // gray
463 
464  //create windows at the readout end of the bar boxes:
465  TGeoBBox* logicBarWin = new TGeoBBox("logicBarWin", 0.5*barnum*(barwidth+2.*barhgap)+bbSideGap, hthick+boxgap, barWin_hthick);
466  TGeoVolume* barwin = new TGeoVolume("DrcBarboxWindowSensor", logicBarWin, gGeoManager->GetMedium("FusedSil"));
467  barwin->SetLineColor(kBlue-4);
468  bbox->AddNode(barwin, 1, new TGeoCombiTrans(0.,0.,-bbox_hlen+barWin_hthick,new TGeoRotation(0)));
469 
470 
471  /*
472  if(fFocusingSystem == 3){
473  TGeoRotation rot_lens3;
474  rot_lens3.RotateZ(-90.);
475  rot_lens3.RotateY( 90.);
476  bbox->AddNode(CylLens1, 1, new TGeoCombiTrans(0., 0., fdz_lens1+2.*barWin_hthick, new TGeoRotation(rot_lens3)));
477  bbox->AddNode(CylLens2, 1, new TGeoCombiTrans(0., 0., fdz_lens1+2.*barWin_hthick, new TGeoRotation(rot_lens3)));
478  }
479  */
480  //create layers of grease at the readout end of the bar boxes, between the windows and the EV:
481  TGeoBBox* logicEVgrease = new TGeoBBox("logicEVgrease", 0.5*barnum*(barwidth+2.*barhgap)+bbSideGap, hthick+boxgap, EVgreaseLayer/2.);
482  TGeoVolume* evgrease = new TGeoVolume("DrcEVgrease", logicEVgrease, gGeoManager->GetMedium("OpticalGrease"));
483  evgrease->SetLineColor(kSpring);
484 
485  // Expansion volume:
486  cout<<"width of the expansion volume = "<<2.*(0.5*barnum*(barwidth+2.*barhgap)-barhgap)<<", step * 3 = "<<step*3.<<endl;
487  TGeoTrap *Trd = new TGeoTrap("Trd",sob_len/2., //1
488  atan(((2.*hthick+EVdrop+EVoffset+sob_len*tan(sob_angle*pi/180.))-(2.*hthick+EVdrop+EVoffset))/(2.*sob_len))*180./pi, //2
489  270., //3
490  (2.*hthick+EVdrop+EVoffset+sob_len*tan(sob_angle*pi/180.))/2.,
491  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea)), //4
492  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea)), //4
493  0, //7
494  (2.*hthick+EVdrop+EVoffset)/2.,
495  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea)), //8
496  0.5*(3.*step - MCPgap - (MCPsize - MCPactiveArea)), //8
497  0); //11
498  TGeoVolume *EVsep = new TGeoVolume("SepEVSensor",Trd,gGeoManager->GetMedium("FusedSil"));
499  EVsep->SetLineColor(kCyan-9);
500  // EVsep->SetTransparency(40);
501 
502  //PhotoDetector Basic material - Carbon. It is placed on the back side of the EV to simulate the support structure of the MCPs, the photons are to hit it in gaps between MCPs:
503  TGeoBBox* logicPDbase = new TGeoBBox("logicPDbase", 0.5*3.*step, (2.*hthick+EVdrop+EVoffset+sob_len*tan(sob_angle*pi/180.))/2.+ 0.5*(MCPsize-MCPactiveArea), PDbaseLayer/2.);
504  TGeoVolume *pdbase = new TGeoVolume("DrcPDbase", logicPDbase, gGeoManager->GetMedium("DIRCcarbonFiber"));
505  pdbase->SetLineColor(kGreen-6);
506  pdbase->SetTransparency(40);
507 
508  //create pixel plates: one for each MCP:
509  TGeoBBox* logicPD = new TGeoBBox("logicPD", MCPactiveArea/2., MCPactiveArea/2., PDsensitiveThick/2.);
510  TGeoVolume *pixelholder = new TGeoVolume("DrcPDSensor", logicPD, gGeoManager->GetMedium("FusedSil"));
511  pixelholder->SetLineColor(kGreen+1);
512 
513  // create photo cathodes for each MCP
514  TGeoBBox* logicPhCathode = new TGeoBBox("logicPhCathode", MCPactiveArea/2., MCPactiveArea/2., PhCathodeThick/2.);
515  TGeoVolume *phcathode = new TGeoVolume("DrcPhCathodeSensor", logicPhCathode, gGeoManager->GetMedium("Photocathode"));//("FusedSil"));
516  phcathode->SetLineColor(kGray+1);
517 
518  // create Windows for each MCP
519  TGeoBBox* logicWindow = new TGeoBBox("logicWindow", MCPsize/2., MCPsize/2., PDwindowThick/2.);
520  TGeoVolume *window = new TGeoVolume("DrcPDwindowSensor", logicWindow, gGeoManager->GetMedium("FusedSil"));
521  window->SetLineColor(kBlue-4);
522 
523  // create grease layers between MCP window and the EV back side
524  TGeoBBox* logicMCPgrease = new TGeoBBox("logicMCPgrease", MCPsize/2., MCPsize/2., PDgreaseLayer/2.);
525  TGeoVolume *mcpgrease = new TGeoVolume("DrcMcpGreaseSensor", logicMCPgrease, gGeoManager->GetMedium("OpticalGrease"));
526  mcpgrease->SetLineColor(kSpring);
527 
528  //put MCPs in the PDbase:
529  Double_t xcurr = 0.;
530  Double_t ycurr = 0.;
531  Int_t TotalNmcp = 0;
532  for(Int_t i=0; i<5; i++){ // loop over y
533  for(Int_t j=0; j<3; j++){ // loop over x
534  xcurr = -0.5*(3.*step) + 0.5*step+j*step;
535  ycurr = -(2.*hthick+EVdrop+EVoffset+sob_len*tan(sob_angle*pi/180.))/2. + 0.5*step+i*step - 0.5*(MCPsize-MCPactiveArea);
536  TotalNmcp = TotalNmcp+1;
537  // put Photo Cathodes
538  pdbase->AddNode(phcathode, TotalNmcp, new TGeoCombiTrans(xcurr, ycurr, PDbaseLayer/2.-PDgreaseLayer-PhCathodeThick/2.-PDwindowThick, new TGeoRotation(0)));
539  //put windows
540  pdbase->AddNode(window, TotalNmcp, new TGeoCombiTrans(xcurr, ycurr, PDbaseLayer/2.-PDgreaseLayer-PDwindowThick/2., new TGeoRotation(0)));
541  //put grease layers
542  pdbase->AddNode(mcpgrease, TotalNmcp, new TGeoCombiTrans(xcurr, ycurr, PDbaseLayer/2.-PDgreaseLayer/2., new TGeoRotation(0)));
543  //put pixel holders:
544  pdbase->AddNode(pixelholder, TotalNmcp, new TGeoCombiTrans(xcurr, ycurr, PDbaseLayer/2.-PDgreaseLayer-PDwindowThick-PhCathodeThick-PDsensitiveThick/2., new TGeoRotation(0)));
545  //TotalNmcp = TotalNmcp + 1;
546  //xpos = xpos + step;
547  }
548  //ypos = ypos + step;
549  }
550 
551  cout<<"bar width = "<<barwidth << ", bar with gaps = "<<(bbX/barnum)<<endl;
552  cout<<"barboxL width = "<<bbX/2.+boxgap+boxthick<<", barboxS width = "<<bbX/2.+boxgap<<endl;
553 
554  // create logic bar:
555  TGeoBBox* logicBar = new TGeoBBox("logicBar", barwidth/2., hthick, bbox_hlen-fabs(len)/2.-mirr_hthick-barWin_hthick);
556 
557  TGeoVolume* bar;
558  bar = new TGeoVolume("DrcBarSensor",logicBar, gGeoManager->GetMedium("FusedSil"));
559  bar->SetLineColor(kCyan-9);
560  bar->SetTransparency(60);
561 
562  // create logic mirror:
563  TGeoBBox* logicMirror = new TGeoBBox("logicMirror", barwidth/2., hthick, mirr_hthick /2.); //$$$$$$$
564  TGeoVolume *mirr = new TGeoVolume("DrcMirr", logicMirror, gGeoManager->GetMedium("Mirror"));
565  mirr->SetLineColor(5);
566 
567  // create glue layer inside the bar (connects two halves):
568  TGeoBBox* logicBarGlue = new TGeoBBox("logicBarGlue", barwidth/2., hthick, gluehthick); //$$$$$$$
569  TGeoVolume *barglue = new TGeoVolume("DrcBarGlueSensor", logicBarGlue, gGeoManager->GetMedium("Epotek301_2"));
570  barglue->SetLineColor(kSpring-5);
571  bar->AddNode(barglue, 1, new TGeoCombiTrans(0., 0., 0., new TGeoRotation (0)));
572 
573  // put barboxes into right positions:
574  Double_t dx_bbox, dy_bbox, dz_bbox, phi_curr;
575  Double_t evLocShift = (Trd->GetH1()+Trd->GetH2())/2.-hthick - EVdrop;
576  Double_t pdLocShift = logicPDbase->GetDY()-hthick-0.5*(MCPsize-MCPactiveArea);
577  for(Int_t m = 0; m < bbnum; m ++){
578  phi_curr = (90. - phi0 - dphi*m)/180.*pi;
579  if(m > bbnum/2-1){ phi_curr = (90. - phi0 - dphi*m - 2.*pipehAngle)/180.*pi; }
580  dx_bbox = radius * cos(phi_curr);
581  dy_bbox = radius * sin(phi_curr);
582  dz_bbox = bbox_shift;
583 
584  TGeoRotation rot_bbox;
585  rot_bbox.RotateZ( -phi0 - m*dphi - (TMath::Floor(2.*m/bbnum))*(2.*pipehAngle));
586  vLocalMother->AddNode(bbox, m+1, new TGeoCombiTrans(dx_bbox, dy_bbox, dz_bbox, new TGeoRotation(rot_bbox)));
587  // vLocalMother->AddNode(evgrease, m+1, new TGeoCombiTrans(dx_bbox, dy_bbox, bbox_zup-EVgreaseLayer/2.,new TGeoRotation(rot_bbox)));
588  vLocalMother->AddNode(EVsep, m+1, new TGeoCombiTrans(dx_bbox+evLocShift*cos(phi_curr), dy_bbox+evLocShift*sin(phi_curr), -bbox_hlen+bbox_shift-sob_len/2.-EVgreaseLayer, new TGeoRotation(rot_bbox)));
589  vLocalMother->AddNode(pdbase, m+1, new TGeoCombiTrans(dx_bbox+pdLocShift*cos(phi_curr), dy_bbox+pdLocShift*sin(phi_curr), -bbox_hlen+bbox_shift-sob_len-EVgreaseLayer-PDbaseLayer/2., new TGeoRotation(rot_bbox)));
590  }
591 
592 
593  Double_t dx, dy, dz_bar, dz_mirr;
594 
595  for(Int_t j=0; j<barnum; j++){
596  dx = - (bbX/2.) - bbSideGap + (barwidth+2.*barhgap)/2. + j * (barwidth+2.*barhgap);
597  dy = 0.;//len/2. - mirr_hthick;
598  if(fFocusingSystem == 0 || fFocusingSystem == 1 || fFocusingSystem == 2){
599  dz_mirr = bbox_hlen - barWin_hthick - mirr_hthick;
600  dz_bar = len/2.-mirr_hthick;
601  }
602  if(fFocusingSystem == 3){
603  dz_mirr = bbox_hlen - barWin_hthick - mirr_hthick;
604  dz_bar = -mirr_hthick + len/2.;
605  }
606  if(fFocusingSystem == 1){ // lens
607  abox->AddNode(lens1, 1+j, new TGeoCombiTrans(dx, dy, fdz_lens1, new TGeoRotation (0)));
608  abox->AddNode(lens2, 1+j, new TGeoCombiTrans(dx, dy, fdz_lens2, new TGeoRotation (0)));
609  abox->AddNode(lens3, 1+j, new TGeoCombiTrans(dx, dy, fdz_lens3, new TGeoRotation (0)));
610  }
611  if(fFocusingSystem == 2){ // mirror
612  abox->AddNode(block1, 1+j, new TGeoCombiTrans(dx, dy, fdz_mirr1, new TGeoRotation (0)));
613  abox->AddNode(block2, 1+j, new TGeoCombiTrans(dx, dy, fdz_mirr2, new TGeoRotation (0)));
614  }
615  if(fFocusingSystem == 3){
616  TGeoRotation rot_lens3;
617  rot_lens3.RotateZ(-90.);
618  rot_lens3.RotateY( 90.);
619  abox->AddNode(CylLens1, 1+j, new TGeoCombiTrans(dx, dy, fdz_lens1, new TGeoRotation(rot_lens3)));
620  abox->AddNode(CylLens2, 1+j, new TGeoCombiTrans(dx, dy, fdz_lens1, new TGeoRotation(rot_lens3)));
621  }
622  abox->AddNode(bar, 1+j, new TGeoCombiTrans(dx, dy, dz_bar, new TGeoRotation(0)));
623  if(fFocusingSystem != 2){ // not a forward mirror
624  abox->AddNode(mirr, 1+j, new TGeoCombiTrans(dx, dy, dz_mirr, new TGeoRotation(0)));
625  }
626  }
627 
628  gGeoManager->CloseGeometry();
629 
630  top->CheckOverlaps(0.0001, "");
631  gGeoManager->CheckOverlaps(0.00001,""); // [cm]
632  gGeoManager->SetVisLevel(2);
633  //gGeoManager->CheckGeometryFull();
634 
635  top->Write();
636  fi->Close();
637  top->Draw("ogl");
638  //bbox->Draw("ogl");
639  //evgrease->Draw("ogl same");
640 
641 
642  TObjArray *listOfOverlaps = gGeoManager->GetListOfOverlaps();
643  cout<<listOfOverlaps->GetEntries()<<endl;
644  listOfOverlaps->Print();
645  return 0;
646 }
CAMath::ASin
static T ASin(const T &x)
cos
friend F32vec4 cos(const F32vec4 &a)
Definition: P4_F32vec4.h:112
PndGeoDrc::EVbackAngle
Double_t EVbackAngle()
Definition: PndGeoDrc.h:148
dy
double dy
Definition: runGeaneTestSimBox.cxx:45
t2
Int_t t2
Definition: hist-t7.C:106
r
double r
Definition: RiemannTest.C:14
PndGeoDrc::BBoxNum
Double_t BBoxNum()
Definition: PndGeoDrc.h:136
geoLoad
FairGeoLoader * geoLoad
Definition: createRootGeoFile.C:16
i
Int_t i
Definition: run_full.C:25
Media
FairGeoMedia * Media
Definition: createRootFscGeometryFile.C:68
m
__m128 m
Definition: P4_F32vec4.h:28
b
TTree * b
Definition: GammaSpectraAnalysis.C:24
gGeoMan
TGeoManager * gGeoMan
Definition: createRootFscGeometryFile.C:85
sin
friend F32vec4 sin(const F32vec4 &a)
Definition: P4_F32vec4.h:111
PndGeoDrc::EVlen
Double_t EVlen()
Definition: PndGeoDrc.h:127
PndGeoDrc::BarWidth
Double_t BarWidth()
Definition: PndGeoDrc.h:100
createRootGeometry_DIRC_sepEV_06_2013
int createRootGeometry_DIRC_sepEV_06_2013(Int_t fFocusingSystem=0, Bool_t fprizm=kFALSE)
Definition: createRootGeometry_DIRC_sepEV_06_2013.C:31
pi
#define pi
Definition: createSTT.C:60
PndGeoDrc::BBoxGap
Double_t BBoxGap()
Definition: PndGeoDrc.h:130
PndGeoDrc::boxGap
Double_t boxGap()
Definition: PndGeoDrc.h:116
gGeoManager
TGeoManager * gGeoManager
Definition: createRootGeoFile.C:97
CAMath::Cos
static T Cos(const T &x)
Definition: PndCAMath.h:43
PndGeoDrc::barBoxZDown
Double_t barBoxZDown()
Definition: PndGeoDrc.h:104
PndGeoDrc::EVoffset
Double_t EVoffset()
Definition: PndGeoDrc.h:145
top
TGeoVolume * top
Definition: createRootFscGeometryFile.C:89
rot1
TGeoRotation * rot1
Definition: createDesySetupAli.C:115
PndGeoDrc::PipehAngle
Double_t PipehAngle()
Definition: PndGeoDrc.h:139
b2
Double_t b2
Definition: createRootGeoFileBwEnd_FwEndCapGeo.C:22
a
Int_t a
Definition: anaLmdDigi.C:126
PndGeoDrc::barBoxZUp
Double_t barBoxZUp()
Definition: PndGeoDrc.h:108
geobuild
FairGeoBuilder * geobuild
Definition: createRootFscGeometryFile.C:69
shape
TGeoShape * shape
Definition: createDummySensor.C:18
fi
TFile * fi
Definition: createRootGeoFile.C:238
Double_t
Double_t
Definition: SimCompleteLinkDef.h:6
phi0
Double_t phi0
Definition: checkhelixhit.C:60
eps
double eps(TVector3 v1, TVector3 v2)
Definition: softtrigger_kin5.C:163
zpos
const Double_t zpos
Definition: createRootFscGeometryFile.C:11
vLocalMother
vLocalMother
Definition: createRootGeoFile.C:121
findSectorOut
int findSectorOut(double y, double dphi_rad, double radiusCornerOut)
Definition: createRootGeometry_DIRC.C:852
t3
Int_t t3
Definition: hist-t7.C:106
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: P4_F32vec4.h:47
PndGeoDrc::boxThick
Double_t boxThick()
Definition: PndGeoDrc.h:120
PndGeoDrc::barhGap
Double_t barhGap()
Definition: PndGeoDrc.h:112
dx
double dx
Definition: runGeaneTestSimBox.cxx:44
PndGeoDrc::barHalfThick
Double_t barHalfThick()
Definition: PndGeoDrc.h:96
t
TTree * t
Definition: bump_analys.C:13
geoFace
FairGeoInterface * geoFace
Definition: createRootGeoFile.C:17
getXrange
void getXrange(double dphi, double radiusCornerOut, double radiusCornerIn, double MCPactiveArea, int cas1, int cas2, double ylow, double step, double &x1, double &x2)
Definition: createRootGeometry_DIRC.C:869
findSectorIn
int findSectorIn(double y, double dphi, double radius, double EVdrop, double hthick)
Definition: createRootGeometry_DIRC.C:836
y
Double_t y
Definition: createRootFscGeometryFile.C:159
alpha
double alpha
Definition: f_Init.h:9
fGeoFile
TString fGeoFile
Definition: createRootGeoFile.C:26
getX
double getX(TVector3 corner1, TVector3 corner2, double yy)
Definition: createRootGeometry_DIRC.C:944
PndGeoDrc::barNum
Double_t barNum()
Definition: PndGeoDrc.h:124
r2
double r2
Definition: RiemannTrackFinderTest.C:19
PndGeoDrc::radius
Double_t radius()
Definition: PndGeoDrc.h:92
Generated on Wed Apr 3 2019 10:02:45 for FairRoot/PandaRoot by   doxygen 1.8.5