SND@LHC Software
Loading...
Searching...
No Matches
ShipMuonShield.cxx
Go to the documentation of this file.
1#include "ShipMuonShield.h"
2
3#include "TGeoManager.h"
4#include "TObjArray.h" // for TObjArray
5#include "TString.h" // for TString
6#include "TGeoBBox.h"
7#include "TGeoTrd1.h"
8#include "TGeoCompositeShape.h"
9#include "TGeoBoolNode.h"
10#include "TGeoTube.h"
11#include "TGeoMaterial.h"
12#include "FairGeoInterface.h"
13#include "FairGeoMedia.h"
14#include "FairGeoBuilder.h"
15#include "FairRuntimeDb.h" // for FairRuntimeDb
16#include "TVectorT.h"
17#include "TFile.h"
18#include <iostream> // for operator<<, basic_ostream, etc
19
20Double_t cm = 1;
21Double_t m = 100 * cm;
22Double_t mm = 0.1 * cm;
23Double_t kilogauss = 1.;
24Double_t tesla = 10 * kilogauss;
25
26ShipMuonShield::~ShipMuonShield() {}
27ShipMuonShield::ShipMuonShield() : FairModule("ShipMuonShield", "") {}
28
29ShipMuonShield::ShipMuonShield(TString geofile,
30 const Int_t withCoMagnet, const Bool_t StepGeo,
31 const Bool_t WithConstAbsorberField, const Bool_t WithConstShieldField)
32 : FairModule("MuonShield", "ShipMuonShield")
33{
34 fWithConstAbsorberField = WithConstAbsorberField;
35 fWithConstShieldField = WithConstShieldField;
36 fStepGeo = StepGeo;
37 fWithCoMagnet = withCoMagnet;
38 fGeofile = geofile;
39 auto f = TFile::Open(geofile, "read");
40 TVectorT<Double_t> params;
41 params.Read("params");
42 Double_t LE = 10. * m, floor = 5. * m;
43 fDesign = 8;
44 fField = 1.7;
45 dZ0 = 1 * m;
46 dZ1 = 0.4 * m;
47 dZ2 = 2.31 * m;
48 dZ3 = params[2];
49 dZ4 = params[3];
50 dZ5 = params[4];
51 dZ6 = params[5];
52 dZ7 = params[6];
53 dZ8 = params[7];
54 fMuonShieldLength = 2 * (dZ1 + dZ2 + dZ3 + dZ4 + dZ5 + dZ6 + dZ7 + dZ8) + LE;
55
56 fFloor = floor;
57 fSupport = true;
58
59 Double_t Z = -25 * m - fMuonShieldLength / 2.;
60
61 zEndOfAbsorb = Z + dZ0 - fMuonShieldLength / 2.;
62 zEndOfAbsorb -= dZ0;
63}
64
65ShipMuonShield::ShipMuonShield(const char* name, const Int_t Design, const char* Title,
66 Double_t Z, Double_t L0, Double_t L1, Double_t L2, Double_t L3, Double_t L4, Double_t L5, Double_t L6,
67 Double_t L7, Double_t L8, Double_t gap, Double_t LE, Double_t, Double_t floor, Double_t field,
68 const Int_t withCoMagnet, const Bool_t StepGeo,
69 const Bool_t WithConstAbsorberField, const Bool_t WithConstShieldField)
70 : FairModule(name ,Title)
71{
72 fDesign = Design;
73 fField = field;
74 fGeofile = "";
75 fWithConstAbsorberField = WithConstAbsorberField;
76 fWithConstShieldField = WithConstShieldField;
77 fStepGeo = StepGeo;
78 fWithCoMagnet = withCoMagnet;
79 if (fDesign==1){
80 fMuonShieldLength = L1;
81 }
82 if (fDesign==2 || fDesign==3 || fDesign==4 ){
83 Fatal("ShipMuonShield","Design %i not anymore supported",fDesign);
84 }
85 if (fDesign==5 || fDesign==6){
86 dZ0 = L0;
87 dZ1 = L1;
88 dZ2 = L2;
89 dZ3 = L3;
90 dZ4 = L4;
91 dZ5 = L5;
92 dZ6 = L6;
93 dZ7 = L7;
94 dZ8 = L8;
95 dXgap= gap;
96 fMuonShieldLength = 2*(dZ1+dZ2+dZ3+dZ4+dZ5+dZ6+dZ7+dZ8) + LE ; //leave some space for nu-tau detector
97 }
98
99 if (fDesign>=7){
100 dZ1 = L1;
101 dZ2 = L2;
102 dZ3 = L3;
103 dZ4 = L4;
104 dZ5 = L5;
105 dZ6 = L6;
106 dZ7 = L7;
107 dZ8 = L8;
108 fMuonShieldLength =
109 2 * (dZ1 + dZ2 + dZ3 + dZ4 + dZ5 + dZ6 + dZ7 + dZ8) + LE;
110 }
111
112 fFloor = (fDesign >= 7) ? floor : 0;
113
114 zEndOfAbsorb = Z + dZ0 - fMuonShieldLength/2.;
115 if(fDesign>=6){zEndOfAbsorb = Z - fMuonShieldLength/2.;}
116 fSupport = true;
117}
118
119// ----- Private method InitMedium
120Int_t ShipMuonShield::InitMedium(TString name)
121{
122 static FairGeoLoader *geoLoad=FairGeoLoader::Instance();
123 static FairGeoInterface *geoFace=geoLoad->getGeoInterface();
124 static FairGeoMedia *media=geoFace->getMedia();
125 static FairGeoBuilder *geoBuild=geoLoad->getGeoBuilder();
126
127 FairGeoMedium *ShipMedium=media->getMedium(name);
128
129 if (!ShipMedium)
130 Fatal("InitMedium","Material %s not defined in media file.", name.Data());
131 TGeoMedium* medium=gGeoManager->GetMedium(name);
132 if (medium)
133 return ShipMedium->getMediumIndex();
134 return geoBuild->createMedium(ShipMedium);
135}
136
137void ShipMuonShield::CreateTube(TString tubeName, TGeoMedium *medium,
138 Double_t dX, Double_t dY, Double_t dZ,
139 Int_t color, TGeoVolume *tShield,
140 Double_t x_translation, Double_t y_translation,
141 Double_t z_translation) {
142 TGeoVolume *absorber = gGeoManager->MakeTube(tubeName, medium, dX, dY, dZ);
143 absorber->SetLineColor(color);
144 tShield->AddNode(
145 absorber, 1,
146 new TGeoTranslation(x_translation, y_translation, z_translation));
147}
148
149void ShipMuonShield::CreateArb8(TString arbName, TGeoMedium *medium,
150 Double_t dZ, std::array<Double_t, 16> corners,
151 Int_t color, TGeoUniformMagField *magField,
152 TGeoVolume *tShield, Double_t x_translation,
153 Double_t y_translation,
154 Double_t z_translation) {
155 TGeoVolume *magF =
156 gGeoManager->MakeArb8(arbName, medium, dZ, corners.data());
157 magF->SetLineColor(color);
158 if (arbName.Contains("Absorb")) {
159 if (fWithConstAbsorberField) {
160 magF->SetField(magField);
161 }
162 } else if (fWithConstShieldField) {
163 magF->SetField(magField);
164 }
165 tShield->AddNode(magF, 1, new TGeoTranslation(x_translation, y_translation,
166 z_translation));
167}
168
169void ShipMuonShield::CreateArb8(TString arbName, TGeoMedium *medium,
170 Double_t dZ, std::array<Double_t, 16> corners,
171 Int_t color, TGeoUniformMagField *magField,
172 TGeoVolume *tShield, Double_t x_translation,
173 Double_t y_translation,
174 Double_t z_translation,
175 Bool_t stepGeo) {
176 if (!stepGeo)
177 {
178 CreateArb8 (arbName, medium, dZ, corners, color, magField, tShield, x_translation, y_translation, z_translation);
179 return;
180 }
181 Double_t partLength = 0.5;
182 Int_t zParts = std::ceil(2.0*dZ/m/partLength);
183 Double_t finalCorners[zParts][16];
184 Double_t dxdy[4][2];
185 Double_t dZp = dZ/Double_t(zParts);
186
187 for (int i = 0; i < 4; ++i)
188 {
189 dxdy[i][0] = (corners[8+2*i] - corners[0+2*i])/Double_t(zParts);
190 dxdy[i][1] = (corners[9+2*i] - corners[1+2*i])/Double_t(zParts);
191 }
192
193 std::copy(corners.data() + 0, corners.data() + 8, finalCorners[0]);
194
195 for (int i = 0; i < zParts; ++i)
196 {
197 for (int k = 0; k < 4; ++k)
198 {
199 finalCorners[i][8+2*k] = finalCorners[i][0+2*k] + dxdy[k][0];
200 finalCorners[i][9+2*k] = finalCorners[i][1+2*k] + dxdy[k][1];
201 }
202 if (i != zParts-1)
203 {
204 std::copy(finalCorners[i] + 8, finalCorners[i] + 16, finalCorners[i+1]);
205 }
206 }
207
208 for (int i = 0; i < zParts; ++i)
209 {
210 for (int k = 0; k < 4; ++k)
211 {
212 finalCorners[i][8+2*k] = finalCorners[i][0+2*k] = (finalCorners[i][0+2*k] + finalCorners[i][8+2*k]) / 2.0;
213 finalCorners[i][9+2*k] = finalCorners[i][1+2*k] = (finalCorners[i][9+2*k] + finalCorners[i][1+2*k]) / 2.0;
214 }
215 }
216
217 std::vector<TGeoVolume*> magF;
218
219 for (int i = 0; i < zParts; ++i)
220 {
221 magF.push_back(gGeoManager->MakeArb8(arbName + '_' + std::to_string(i), medium, dZp - 0.00001*m, finalCorners[i]));
222 magF[i]->SetLineColor(color);
223 if (fWithConstShieldField) {
224 magF[i]->SetField(magField);
225 }
226 }
227
228 for (int i = 0; i < zParts; ++i)
229 {
230 Double_t true_z_translation = z_translation + 2.0 * Double_t(i) * dZp - dZ + dZp;
231 tShield->AddNode(magF[i], 1, new TGeoTranslation(x_translation, y_translation, true_z_translation));
232 }
233}
234
235void ShipMuonShield::CreateMagnet(TString magnetName,TGeoMedium* medium,TGeoVolume *tShield,TGeoUniformMagField *fields[4],FieldDirection fieldDirection,
236 Double_t dX, Double_t dY, Double_t dX2, Double_t dY2, Double_t dZ,
237 Double_t middleGap,Double_t middleGap2,
238 Double_t HmainSideMag, Double_t HmainSideMag2,
239 Double_t gap,Double_t gap2, Double_t Z, Bool_t NotMagnet,
240 Bool_t stepGeo)
241 {
242 Double_t coil_gap,coil_gap2;
243 Int_t color[4] = {45,31,30,38};
244
245 if (NotMagnet) {
246 coil_gap = gap;
247 coil_gap2 = gap2;
248 } else if (fDesign > 7) {
249 // Assuming 0.5A/mm^2 and 10000At needed, about 200cm^2 gaps are necessary
250 // Current design safely above this. Will consult with MISiS to get a better minimum.
251 gap = std::ceil(std::max(100. / dY, gap));
252 gap2 = std::ceil(std::max(100. / dY2, gap2));
253 coil_gap = gap;
254 coil_gap2 = gap2;
255 } else {
256 coil_gap = std::max(20., gap);
257 coil_gap2 = std::max(20., gap2);
258 gap = std::max(2., gap);
259 gap2 = std::max(2., gap2);
260 }
261
262 Double_t anti_overlap = (fDesign == 5) ? 0.0 : 0.1; // gap between fields in the
263 // corners for mitred joints
264 // (Geant goes crazy when
265 // they touch each other)
266
267 std::array<Double_t, 16> cornersMainL = {
268 middleGap, -(dY + dX - anti_overlap),
269 middleGap, dY + dX - anti_overlap,
270 dX + middleGap, dY - anti_overlap,
271 dX + middleGap, -(dY - anti_overlap),
272 middleGap2, -(dY2 + dX2 - anti_overlap),
273 middleGap2, dY2 + dX2 - anti_overlap,
274 dX2 + middleGap2, dY2 - anti_overlap,
275 dX2 + middleGap2, -(dY2 - anti_overlap)
276 };
277
278 std::array<Double_t, 16> cornersTL = {middleGap + dX,
279 dY,
280 middleGap,
281 dY + dX,
282 2 * dX + middleGap + coil_gap,
283 dY + dX,
284 dX + middleGap + coil_gap,
285 dY,
286 middleGap2 + dX2,
287 dY2,
288 middleGap2,
289 dY2 + dX2,
290 2 * dX2 + middleGap2 + coil_gap2,
291 dY2 + dX2,
292 dX2 + middleGap2 + coil_gap2,
293 dY2};
294
295 std::array<Double_t, 16> cornersMainSideL =
296 fDesign <= 7 ?
297 std::array<Double_t, 16>{
298 dX + middleGap + gap, -HmainSideMag,
299 dX + middleGap + gap, HmainSideMag,
300 2 * dX + middleGap + gap, HmainSideMag,
301 2 * dX + middleGap + gap, -HmainSideMag,
302 dX2 + middleGap2 + gap2, -HmainSideMag2,
303 dX2 + middleGap2 + gap2, HmainSideMag2,
304 2 * dX2 + middleGap2 + gap2, HmainSideMag2,
305 2 * dX2 + middleGap2 + gap2, -HmainSideMag2
306 } :
307 std::array<Double_t, 16>{
308 dX + middleGap + gap, -(dY - anti_overlap),
309 dX + middleGap + gap, dY - anti_overlap,
310 2 * dX + middleGap + gap, dY + dX - anti_overlap,
311 2 * dX + middleGap + gap, -(dY + dX - anti_overlap),
312 dX2 + middleGap2 + gap2, -(dY2 - anti_overlap),
313 dX2 + middleGap2 + gap2, dY2 - anti_overlap,
314 2 * dX2 + middleGap2 + gap2, dY2 + dX2 - anti_overlap,
315 2 * dX2 + middleGap2 + gap2, -(dY2 + dX2 - anti_overlap)
316 };
317
318 std::array<Double_t, 16> cornersMainR, cornersCLBA,
319 cornersMainSideR, cornersCLTA, cornersCRBA,
320 cornersCRTA, cornersTR, cornersBL, cornersBR;
321
322 if (fDesign <= 7) {
323 cornersCLBA = {dX + middleGap + gap,
324 -HmainSideMag,
325 2 * dX + middleGap + gap,
326 -HmainSideMag,
327 2 * dX + middleGap + coil_gap,
328 -(dY + dX - anti_overlap),
329 dX + middleGap + coil_gap,
330 -(dY - anti_overlap),
331 dX2 + middleGap2 + gap2,
332 -HmainSideMag2,
333 2 * dX2 + middleGap2 + gap2,
334 -HmainSideMag2,
335 2 * dX2 + middleGap2 + coil_gap2,
336 -(dY2 + dX2 - anti_overlap),
337 dX2 + middleGap2 + coil_gap2,
338 -(dY2 - anti_overlap)};
339 }
340
341 // Use symmetries to define remaining magnets
342 for (int i = 0; i < 16; ++i) {
343 cornersMainR[i] = -cornersMainL[i];
344 cornersMainSideR[i] = -cornersMainSideL[i];
345 cornersCRTA[i] = -cornersCLBA[i];
346 cornersBR[i] = -cornersTL[i];
347 }
348 // Need to change order as corners need to be defined clockwise
349 for (int i = 0, j = 4; i < 8; ++i) {
350 j = (11 - i) % 8;
351 cornersCLTA[2 * j] = cornersCLBA[2 * i];
352 cornersCLTA[2 * j + 1] = -cornersCLBA[2 * i + 1];
353 cornersTR[2 * j] = -cornersTL[2 * i];
354 cornersTR[2 * j + 1] = cornersTL[2 * i + 1];
355 }
356 for (int i = 0; i < 16; ++i) {
357 cornersCRBA[i] = -cornersCLTA[i];
358 cornersBL[i] = -cornersTR[i];
359 }
360
361 TString str1L = "_MiddleMagL";
362 TString str1R = "_MiddleMagR";
363 TString str2 = "_MagRetL";
364 TString str3 = "_MagRetR";
365 TString str4 = "_MagCLB";
366 TString str5 = "_MagCLT";
367 TString str6 = "_MagCRT";
368 TString str7 = "_MagCRB";
369 TString str8 = "_MagTopLeft";
370 TString str9 = "_MagTopRight";
371 TString str10 = "_MagBotLeft";
372 TString str11 = "_MagBotRight";
373
374 switch (fieldDirection){
375
376 case FieldDirection::up:
377 CreateArb8(magnetName + str1L, medium, dZ, cornersMainL, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
378 CreateArb8(magnetName + str1R, medium, dZ, cornersMainR, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
379 CreateArb8(magnetName + str2, medium, dZ, cornersMainSideL, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
380 CreateArb8(magnetName + str3, medium, dZ, cornersMainSideR, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
381 if (fDesign <= 7) {
382 CreateArb8(magnetName + str4, medium, dZ, cornersCLBA, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
383 CreateArb8(magnetName + str5, medium, dZ, cornersCLTA, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
384 CreateArb8(magnetName + str6, medium, dZ, cornersCRTA, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
385 CreateArb8(magnetName + str7, medium, dZ, cornersCRBA, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
386 }
387 CreateArb8(magnetName + str8, medium, dZ, cornersTL, color[3], fields[3], tShield, 0, 0, Z, stepGeo);
388 CreateArb8(magnetName + str9, medium, dZ, cornersTR, color[2], fields[2], tShield, 0, 0, Z, stepGeo);
389 CreateArb8(magnetName + str10, medium, dZ, cornersBL, color[2], fields[2], tShield, 0, 0, Z, stepGeo);
390 CreateArb8(magnetName + str11, medium, dZ, cornersBR, color[3], fields[3], tShield, 0, 0, Z, stepGeo);
391 break;
392 case FieldDirection::down:
393 CreateArb8(magnetName + str1L, medium, dZ, cornersMainL, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
394 CreateArb8(magnetName + str1R, medium, dZ, cornersMainR, color[1], fields[1], tShield, 0, 0, Z, stepGeo);
395 CreateArb8(magnetName + str2, medium, dZ, cornersMainSideL, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
396 CreateArb8(magnetName + str3, medium, dZ, cornersMainSideR, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
397 if (fDesign <= 7) {
398 CreateArb8(magnetName + str4, medium, dZ, cornersCLBA, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
399 CreateArb8(magnetName + str5, medium, dZ, cornersCLTA, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
400 CreateArb8(magnetName + str6, medium, dZ, cornersCRTA, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
401 CreateArb8(magnetName + str7, medium, dZ, cornersCRBA, color[0], fields[0], tShield, 0, 0, Z, stepGeo);
402 }
403 CreateArb8(magnetName + str8, medium, dZ, cornersTL, color[2], fields[2], tShield, 0, 0, Z, stepGeo);
404 CreateArb8(magnetName + str9, medium, dZ, cornersTR, color[3], fields[3], tShield, 0, 0, Z, stepGeo);
405 CreateArb8(magnetName + str10, medium, dZ, cornersBL, color[3], fields[3], tShield, 0, 0, Z, stepGeo);
406 CreateArb8(magnetName + str11, medium, dZ, cornersBR, color[2], fields[2], tShield, 0, 0, Z, stepGeo);
407 break;
408 }
409 }
410
411Int_t ShipMuonShield::Initialize(std::vector<TString> &magnetName,
412 std::vector<FieldDirection> &fieldDirection,
413 std::vector<Double_t> &dXIn, std::vector<Double_t> &dYIn,
414 std::vector<Double_t> &dXOut, std::vector<Double_t> &dYOut,
415 std::vector<Double_t> &dZ, std::vector<Double_t> &midGapIn,
416 std::vector<Double_t> &midGapOut,
417 std::vector<Double_t> &HmainSideMagIn,
418 std::vector<Double_t> &HmainSideMagOut,
419 std::vector<Double_t> &gapIn, std::vector<Double_t> &gapOut,
420 std::vector<Double_t> &Z) {
421
422 const Int_t nMagnets = (fDesign >= 7) ? 9 : 8;
423 magnetName.reserve(nMagnets);
424 fieldDirection.reserve(nMagnets);
425 for (auto i :
426 {&dXIn, &dXOut, &dYIn, &dYOut, &dZ, &midGapIn, &midGapOut,
427 &HmainSideMagIn, &HmainSideMagOut, &gapIn, &gapOut, &Z}) {
428 i->reserve(nMagnets);
429 }
430
431 Double_t zgap = (fDesign > 6) ? 10 : 0; // fixed distance between magnets in Z-axis
432
433 if (fDesign == 8) {
434
435 magnetName = {"MagnAbsorb1", "MagnAbsorb2", "Magn1", "Magn2", "Magn3",
436 "Magn4", "Magn5", "Magn6", "Magn7"};
437
438 fieldDirection = {
439 FieldDirection::up, FieldDirection::up, FieldDirection::up,
440 FieldDirection::up, FieldDirection::up, FieldDirection::down,
441 FieldDirection::down, FieldDirection::down, FieldDirection::down,
442 };
443
444 auto f = TFile::Open(fGeofile, "read");
445 TVectorT<Double_t> params;
446 params.Read("params");
447
448 const int offset = 7;
449
450 dXIn[0] = 0.4 * m;
451 dXOut[0] = 0.40 * m;
452 gapIn[0] = 0.1 * mm;
453 dYIn[0] = 1.5 * m;
454 dYOut[0] = 1.5 * m;
455 gapOut[0] = 0.1 * mm;
456 dXIn[1] = 0.5 * m;
457 dXOut[1] = 0.5 * m;
458 gapIn[1] = 0.02 * m;
459 dYIn[1] = 1.3 * m;
460 dYOut[1] = 1.3 * m;
461 gapOut[1] = 0.02 * m;
462
463 for (Int_t i = 2; i < nMagnets - 1; ++i) {
464 dXIn[i] = params[offset + i * 6 + 1];
465 dXOut[i] = params[offset + i * 6 + 2];
466 dYIn[i] = params[offset + i * 6 + 3];
467 dYOut[i] = params[offset + i * 6 + 4];
468 gapIn[i] = params[offset + i * 6 + 5];
469 gapOut[i] = params[offset + i * 6 + 6];
470 }
471
472 dZ[0] = dZ1 - zgap / 2;
473 Z[0] = zEndOfAbsorb + dZ[0] + zgap;
474 dZ[1] = dZ2 - zgap / 2;
475 Z[1] = Z[0] + dZ[0] + dZ[1] + zgap;
476 dZ[2] = dZ3 - zgap / 2;
477 Z[2] = Z[1] + dZ[1] + dZ[2] + 2 * zgap;
478 dZ[3] = dZ4 - zgap / 2;
479 Z[3] = Z[2] + dZ[2] + dZ[3] + zgap;
480 dZ[4] = dZ5 - zgap / 2;
481 Z[4] = Z[3] + dZ[3] + dZ[4] + zgap;
482 dZ[5] = dZ6 - zgap / 2;
483 Z[5] = Z[4] + dZ[4] + dZ[5] + zgap;
484 dZ[6] = dZ7 - zgap / 2;
485 Z[6] = Z[5] + dZ[5] + dZ[6] + zgap;
486 dZ[7] = dZ8 - zgap / 2;
487 Z[7] = Z[6] + dZ[6] + dZ[7] + zgap;
488
489 dXIn[8] = dXOut[7];
490 dYIn[8] = dYOut[7];
491 dXOut[8] = dXIn[8];
492 dYOut[8] = dYIn[8];
493 gapIn[8] = gapOut[7];
494 gapOut[8] = gapIn[8];
495 dZ[8] = 0.1 * m;
496 Z[8] = Z[7] + dZ[7] + dZ[8];
497
498 for (int i = 0; i < nMagnets; ++i) {
499 midGapIn[i] = 0.;
500 midGapOut[i] = 0.;
501 HmainSideMagIn[i] = dYIn[i] / 2;
502 HmainSideMagOut[i] = dYOut[i] / 2;
503 }
504
505 } else if (fDesign == 9) {
506 magnetName = {"MagnAbsorb1", "MagnAbsorb2", "Magn1", "Magn2", "Magn3",
507 "Magn4", "Magn5", "Magn6", "Magn7"
508 };
509
510 fieldDirection = {
511 FieldDirection::up, FieldDirection::up, FieldDirection::up,
512 FieldDirection::up, FieldDirection::up, FieldDirection::down,
513 FieldDirection::down, FieldDirection::down, FieldDirection::down,
514 };
515
516 dXIn[0] = 0.4 * m;
517 dXOut[0] = 0.40 * m;
518 dYIn[0] = 1.5 * m;
519 dYOut[0] = 1.5 * m;
520 gapIn[0] = 0.1 * mm;
521 gapOut[0] = 0.1 * mm;
522 dZ[0] = dZ1 - zgap / 2;
523 Z[0] = zEndOfAbsorb + dZ[0] + zgap;
524
525 dXIn[1] = 0.5 * m;
526 dXOut[1] = 0.5 * m;
527 dYIn[1] = 1.3 * m;
528 dYOut[1] = 1.3 * m;
529 gapIn[1] = 0.02 * m;
530 gapOut[1] = 0.02 * m;
531 dZ[1] = dZ2 - zgap / 2;
532 Z[1] = Z[0] + dZ[0] + dZ[1] + zgap;
533
534 dXIn[2] = 0.72 * m;
535 dXOut[2] = 0.51 * m;
536 dYIn[2] = 0.29 * m;
537 dYOut[2] = 0.46 * m;
538 gapIn[2] = 0.10 * m;
539 gapOut[2] = 0.07 * m;
540 dZ[2] = dZ3 - zgap / 2;
541 Z[2] = Z[1] + dZ[1] + dZ[2] + 2 * zgap;
542
543 dXIn[3] = 0.54 * m;
544 dXOut[3] = 0.38 * m;
545 dYIn[3] = 0.46 * m;
546 dYOut[3] = 1.92 * m;
547 gapIn[3] = 0.14 * m;
548 gapOut[3] = 0.09 * m;
549 dZ[3] = dZ4 - zgap / 2;
550 Z[3] = Z[2] + dZ[2] + dZ[3] + zgap;
551
552 dXIn[4] = 0.10 * m;
553 dXOut[4] = 0.31 * m;
554 dYIn[4] = 0.35 * m;
555 dYOut[4] = 0.31 * m;
556 gapIn[4] = 0.51 * m;
557 gapOut[4] = 0.11 * m;
558 dZ[4] = dZ5 - zgap / 2;
559 Z[4] = Z[3] + dZ[3] + dZ[4] + zgap;
560
561 dXIn[5] = 0.03 * m;
562 dXOut[5] = 0.32 * m;
563 dYIn[5] = 0.54 * m;
564 dYOut[5] = 0.24 * m;
565 gapIn[5] = 0.08 * m;
566 gapOut[5] = 0.08 * m;
567 dZ[5] = dZ6 - zgap / 2;
568 Z[5] = Z[4] + dZ[4] + dZ[5] + zgap;
569
570 dXIn[6] = 0.22 * m;
571 dXOut[6] = 0.32 * m;
572 dYIn[6] = 2.09 * m;
573 dYOut[6] = 0.35 * m;
574 gapIn[6] = 0.08 * m;
575 gapOut[6] = 0.13 * m;
576 dZ[6] = dZ7 - zgap / 2;
577 Z[6] = Z[5] + dZ[5] + dZ[6] + zgap;
578
579 dXIn[7] = 0.33 * m;
580 dXOut[7] = 0.77 * m;
581 dYIn[7] = 0.85 * m;
582 dYOut[7] = 2.41 * m;
583 gapIn[7] = 0.09 * m;
584 gapOut[7] = 0.26 * m;
585 dZ[7] = dZ8 - zgap / 2;
586 Z[7] = Z[6] + dZ[6] + dZ[7] + zgap;
587
588 dXIn[8] = dXOut[7];
589 dYIn[8] = dYOut[7];
590 dXOut[8] = dXIn[8];
591 dYOut[8] = dYIn[8];
592 gapIn[8] = gapOut[7];
593 gapOut[8] = gapIn[8];
594 dZ[8] = 0.1 * m;
595 Z[8] = Z[7] + dZ[7] + dZ[8];
596
597 for (int i = 0; i < nMagnets; ++i) {
598 midGapIn[i] = 0.;
599 midGapOut[i] = 0.;
600 HmainSideMagIn[i] = dYIn[i] / 2;
601 HmainSideMagOut[i] = dYOut[i] / 2;
602 }
603
604 } else if (fDesign == 7) {
605 magnetName = {"MagnAbsorb1", "MagnAbsorb2", "Magn1", "Magn2", "Magn3",
606 "Magn4", "Magn5", "Magn6", "Magn7"};
607
608 fieldDirection[0] = FieldDirection::up;
609 dXIn[0] = 0.4*m; dYIn[0] = 1.5*m;
610 dXOut[0] = 0.40*m; dYOut[0]= 1.5*m;
611 gapIn[0] = 0.02 * m; gapOut[0] = 0.02 * m;
612 dZ[0] = dZ1-zgap/2; Z[0] = zEndOfAbsorb + dZ[0]+zgap;
613
614 fieldDirection[1] = FieldDirection::up;
615 dXIn[1] = 0.8*m; dYIn[1] = 1.5*m;
616 dXOut[1] = 0.8*m; dYOut[1]= 1.5*m;
617 gapIn[1] = 0.02*m; gapOut[1] = 0.02*m;
618 dZ[1] = dZ2-zgap/2; Z[1] = Z[0] + dZ[0] + dZ[1]+zgap;
619
620 fieldDirection[2] = FieldDirection::up;
621 dXIn[2] = 0.87*m; dYIn[2] = 0.35*m;
622 dXOut[2] = 0.65*m; dYOut[2]= 1.21*m;
623 gapIn[2] = 0.11 * m; gapOut[2] = 0.065 * m;
624 dZ[2] = dZ3-zgap/2; Z[2] = Z[1] + dZ[1] + dZ[2] + zgap;
625
626 fieldDirection[3] = FieldDirection::up;
627 dXIn[3] = 0.65*m; dYIn[3] = 1.21*m;
628 dXOut[3] = 0.43*m; dYOut[3]= 2.07*m;
629 gapIn[3] = 0.065 * m; gapOut[3] = 0.02 * m;
630 dZ[3] = dZ4-zgap/2; Z[3] = Z[2] + dZ[2] + dZ[3]+zgap;
631
632 fieldDirection[4] = FieldDirection::up;
633 dXIn[4] = 0.06*m; dYIn[4] = 0.32*m;
634 dXOut[4] = 0.33*m; dYOut[4]= 0.13*m;
635 gapIn[4] = 0.7*m; gapOut[4] = 0.11*m;
636 dZ[4] = dZ5-zgap/2; Z[4] = Z[3] + dZ[3] + dZ[4]+zgap;
637
638 fieldDirection[5] = FieldDirection::down;
639 dXIn[5] = 0.05*m; dYIn[5] = 1.12*m;
640 dXOut[5] =0.16*m; dYOut[5]= 0.05*m;
641 gapIn[5] = 0.04*m; gapOut[5] = 0.02*m;
642 dZ[5] = dZ6-zgap/2; Z[5] = Z[4] + dZ[4] + dZ[5]+zgap;
643
644 fieldDirection[6] = FieldDirection::down;
645 dXIn[6] = 0.15*m; dYIn[6] = 2.35*m;
646 dXOut[6] = 0.34*m; dYOut[6]= 0.32*m;
647 gapIn[6] = 0.05*m; gapOut[6] = 0.08*m;
648 dZ[6] = dZ7-zgap/2; Z[6] = Z[5] + dZ[5] + dZ[6]+zgap;
649
650 Double_t clip_width = 0.1*m; // clip field width by this width
651 fieldDirection[7] = FieldDirection::down;
652 dXIn[7] = 0.31*m; dYIn[7] = 1.86*m;
653 dXOut[7] = 0.9*m - clip_width; dYOut[7]= 3.1*m;
654 Double_t clip_len =
655 (dZ8-zgap/2) * (1 - (dXOut[7] - dXIn[7]) / (dXOut[7] + clip_width - dXIn[7]));
656 gapIn[7] = 0.02*m; gapOut[7] = 0.55*m;
657 dZ[7] = dZ8 - clip_len - zgap / 2; Z[7] = Z[6] + dZ[6] + dZ[7] + zgap;
658
659 fieldDirection[8] = FieldDirection::down;
660 dXIn[8] = dXOut[7]; dYIn[8] = dYOut[7];
661 dXOut[8] = dXOut[7]; dYOut[8]= dYOut[7];
662 gapIn[8] = 0.55*m; gapOut[8] = 0.55*m;
663 dZ[8] = clip_len; Z[8] = Z[7] + dZ[7] + dZ[8];
664
665 for (int i = 0; i < nMagnets; ++i) {
666 midGapIn[i] = 0.;
667 midGapOut[i] = 0.;
668 HmainSideMagIn[i] = dYIn[i] / 2;
669 HmainSideMagOut[i] = dYOut[i] / 2;
670 }
671
672 } else {
673
674 magnetName = {"1", "2", "3", "4", "5", "6", "7"};
675
676 fieldDirection[0] = FieldDirection::up;
677 dXIn[0] = 0.7*m; dYIn[0] = 1.*m;
678 dXOut[0] = 0.7*m; dYOut[0]= 0.8158*m;
679 midGapIn[0] = 0; midGapOut[0] = 0;
680 HmainSideMagIn[0] = dYIn[0]; HmainSideMagOut[0] = dYOut[0];
681 gapIn[0] = 20; gapOut[0] = 20;
682 dZ[0] = dZ1-zgap; Z[0] = zEndOfAbsorb + dZ[0]+zgap;
683
684 fieldDirection[1] = FieldDirection::up;
685 dXIn[1] = 0.36*m; dYIn[1] = 0.8158*m;
686 dXOut[1] = 0.19*m; dYOut[1]= 0.499*m;
687 midGapIn[1] = 0; midGapOut[1] = 0;
688 HmainSideMagIn[1] = dYIn[1]/2; HmainSideMagOut[1] = dYOut[1]/2;
689 gapIn[1] = 88; gapOut[1] = 122;
690 dZ[1] = dZ2-zgap/2; Z[1] = Z[0] + dZ[0] + dZ[1]+zgap;
691
692 fieldDirection[2] = FieldDirection::down;
693 dXIn[2] = 0.075*m; dYIn[2] = 0.499*m;
694 dXOut[2] = 0.25*m; dYOut[2]= 1.10162*m;
695 midGapIn[2] = 0; midGapOut[2] = 0;
696 HmainSideMagIn[2] = dYIn[2]/2; HmainSideMagOut[2] = dYOut[2]/2;
697 gapIn[2] = 0; gapOut[2] = 0;
698 dZ[2] = dZ3-zgap/2; Z[2] = Z[1] + dZ[1] + dZ[2]+zgap;
699
700 fieldDirection[3] = FieldDirection::down;
701 dXIn[3] = 0.25*m; dYIn[3] = 1.10262*m;
702 dXOut[3] = 0.3*m; dYOut[3]= 1.82697*m;
703 midGapIn[3] = 0; midGapOut[3] = 0;
704 HmainSideMagIn[3] = dXIn[3]; HmainSideMagOut[3] = dXOut[3];
705 gapIn[3] = 0; gapOut[3] = 25;
706 dZ[3] = dZ4-zgap/2; Z[3] = Z[2] + dZ[2] + dZ[3]+zgap;
707
708 fieldDirection[4] = FieldDirection::down;
709 dXIn[4] = 0.3*m; dYIn[4] = 1.82697*m;
710 dXOut[4] = 0.4*m; dYOut[4]= 2.55131*m;
711 midGapIn[4] = 5; midGapOut[4] = 25;
712 HmainSideMagIn[4] = dXIn[4]; HmainSideMagOut[4] = dXOut[4];
713 gapIn[4] = 20; gapOut[4] = 20;
714 dZ[4] = dZ6-zgap/2; Z[4] = Z[3] + dZ[3] + dZ[4]+zgap;
715
716 fieldDirection[5] = FieldDirection::down;
717 dXIn[5] = 0.4*m; dYIn[5] = 2.55131*m;
718 dXOut[5] =0.4*m; dYOut[5]= 3.27566*m;
719 midGapIn[5] = 25; midGapOut[5] = 65;
720 HmainSideMagIn[5] = dXIn[5]; HmainSideMagOut[5] = dXOut[5];
721 gapIn[5] = 20; gapOut[5] = 20;
722 dZ[5] = dZ7-zgap/2; Z[5] = Z[4] + dZ[4] + dZ[5]+zgap;
723
724 fieldDirection[6] = FieldDirection::down;
725 dXIn[6] = 0.4*m; dYIn[6] = 3.27566*m;
726 dXOut[6] = 0.75*m; dYOut[6]= 4*m;
727 midGapIn[6] = 65; midGapOut[6] = 75;
728 HmainSideMagIn[6] = dXIn[6]; HmainSideMagOut[6] = dXOut[6];
729 gapIn[6] = 20; gapOut[6] = 20;
730 dZ[6] = dZ8-zgap/2; Z[6] = Z[5] + dZ[5] + dZ[6]+zgap;
731 }
732 return nMagnets;
733}
734void ShipMuonShield::ConstructGeometry()
735{
736 TGeoVolume *top=gGeoManager->GetTopVolume();
737 TGeoVolume *tShield = new TGeoVolumeAssembly("MuonShieldArea");
738 InitMedium("steel");
739 TGeoMedium *steel =gGeoManager->GetMedium("steel");
740 InitMedium("iron");
741 TGeoMedium *iron =gGeoManager->GetMedium("iron");
742 InitMedium("Concrete");
743 TGeoMedium *concrete =gGeoManager->GetMedium("Concrete");
744
745 if (fDesign >= 5 && fDesign <= 9) {
746 Double_t ironField = fField*tesla;
747 TGeoUniformMagField *magFieldIron = new TGeoUniformMagField(0.,ironField,0.);
748 TGeoUniformMagField *RetField = new TGeoUniformMagField(0.,-ironField,0.);
749 TGeoUniformMagField *ConRField = new TGeoUniformMagField(-ironField,0.,0.);
750 TGeoUniformMagField *ConLField = new TGeoUniformMagField(ironField,0.,0.);
751 TGeoUniformMagField *fields[4] = {magFieldIron,RetField,ConRField,ConLField};
752
753 std::vector<TString> magnetName;
754 std::vector<FieldDirection> fieldDirection;
755 std::vector<Double_t> dXIn, dYIn, dXOut, dYOut, dZf, midGapIn, midGapOut,
756 HmainSideMagIn, HmainSideMagOut, gapIn, gapOut, Z;
757 const Int_t nMagnets = Initialize(magnetName, fieldDirection, dXIn, dYIn, dXOut, dYOut, dZf,
758 midGapIn, midGapOut, HmainSideMagIn, HmainSideMagOut, gapIn,
759 gapOut, Z);
760
761
762 if (fDesign == 6){
763 Double_t dA = 3*m;
764 CreateMagnet("AbsorberStop-1",iron,tShield,fields,FieldDirection::up,
765 dA/6.,dA/6.,dA/6.,dA/6.,dZ0/3.,0,0,dA/12.,dA/12.,0,0,zEndOfAbsorb - 5.*dZ0/3.,0, false);
766 CreateMagnet("AbsorberStop-2",iron,tShield,fields,FieldDirection::up,
767 dA/2.,dA/2.,dA/2.,dA/2.,dZ0*2./3.,0,0,dA/4.,dA/4.,0,0,zEndOfAbsorb - 2.*dZ0/3.,0, false);
768 TGeoBBox* fullAbsorber = new TGeoBBox("fullAbsorber", dA, dA, dZ0/3.);
769 TGeoBBox* cutOut = new TGeoBBox("cutout", dA/3.+20*cm, dA/3.+20*cm, dZ0/3.+0.1*mm); //no idea why to add 20cm
770 TGeoSubtraction *subtraction = new TGeoSubtraction("fullAbsorber","cutout");
771 TGeoCompositeShape *Tc = new TGeoCompositeShape("passiveAbsorberStopSubtr", subtraction);
772 TGeoVolume* passivAbsorber = new TGeoVolume("passiveAbsorberStop-1",Tc, iron);
773 tShield->AddNode(passivAbsorber, 1, new TGeoTranslation(0,0,zEndOfAbsorb - 5.*dZ0/3.));
774 } else if (fDesign >= 7) {
775 float mField = 1.6 * tesla;
776 TGeoUniformMagField *fieldsAbsorber[4] = {
777 new TGeoUniformMagField(0., mField, 0.),
778 new TGeoUniformMagField(0., -mField, 0.),
779 new TGeoUniformMagField(-mField, 0., 0.),
780 new TGeoUniformMagField(mField, 0., 0.)
781 };
782
783 for (Int_t nM = (fDesign == 7) ? 0 : 1; nM < 2; nM++) {
784 CreateMagnet(magnetName[nM], iron, tShield, fieldsAbsorber,
785 fieldDirection[nM], dXIn[nM], dYIn[nM], dXOut[nM],
786 dYOut[nM], dZf[nM], midGapIn[nM], midGapOut[nM],
787 HmainSideMagIn[nM], HmainSideMagOut[nM], gapIn[nM],
788 gapOut[nM], Z[nM], true, false);
789 }
790
791 std::vector<TGeoTranslation*> mag_trans;
792
793 if (fDesign == 7) {
794 auto mag1 = new TGeoTranslation("mag1", 0, 0, -dZ2);
795 mag1->RegisterYourself();
796 mag_trans.push_back(mag1);
797 }
798 auto mag2 = new TGeoTranslation("mag2", 0, 0, +dZ1);
799 mag2->RegisterYourself();
800 mag_trans.push_back(mag2);
801
802 Double_t zgap = 10;
803 Double_t absorber_offset = zgap;
804 Double_t absorber_half_length = (dZf[0] + dZf[1]) + zgap / 2.;
805 auto abs = new TGeoBBox("absorber", 3.95 * m, 3.4 * m, absorber_half_length);
806 const std::vector<TString> absorber_magnets =
807 (fDesign == 7) ? std::vector<TString>{"MagnAbsorb1", "MagnAbsorb2"} : std::vector<TString>{"MagnAbsorb2"};
808 const std::vector<TString> magnet_components = fDesign == 7 ? std::vector<TString>{
809 "_MiddleMagL", "_MiddleMagR", "_MagRetL", "_MagRetR",
810 "_MagCLB", "_MagCLT", "_MagCRT", "_MagCRB",
811 "_MagTopLeft", "_MagTopRight", "_MagBotLeft", "_MagBotRight",
812 }: std::vector<TString>{
813 "_MiddleMagL", "_MiddleMagR", "_MagRetL", "_MagRetR",
814 "_MagTopLeft", "_MagTopRight", "_MagBotLeft", "_MagBotRight",
815 };
816 TString absorber_magnet_components;
817 for (auto &&magnet_component : magnet_components) {
818 // format: "-<magnetName>_<magnet_component>:<translation>"
819 absorber_magnet_components +=
820 ("-" + absorber_magnets[0] + magnet_component + ":" +
821 mag_trans[0]->GetName());
822 if (fDesign == 7) {
823 absorber_magnet_components +=
824 ("-" + absorber_magnets[1] + magnet_component + ":" +
825 mag_trans[1]->GetName());
826 }
827 }
828 TGeoCompositeShape *absorberShape = new TGeoCompositeShape(
829 "Absorber", "absorber" + absorber_magnet_components); // cutting out
830 // magnet parts
831 // from absorber
832 TGeoVolume *absorber = new TGeoVolume("AbsorberVol", absorberShape, iron);
833 absorber->SetLineColor(42); // brown / light red
834 tShield->AddNode(absorber, 1, new TGeoTranslation(0, 0, zEndOfAbsorb + absorber_half_length + absorber_offset));
835
836 if (fDesign > 7) {
837 auto coatBox = new TGeoBBox("coat", 10 * m - 1 * mm, 10 * m - 1 * mm, absorber_half_length);
838 auto coatShape = new TGeoCompositeShape("CoatShape", "coat-absorber");
839 auto coat = new TGeoVolume("CoatVol", coatShape, concrete);
840 tShield->AddNode(coat, 1, new TGeoTranslation(0, 0, zEndOfAbsorb + absorber_half_length + absorber_offset ));
841 TGeoVolume *coatWall = gGeoManager->MakeBox("CoatWall",concrete,10 * m - 1 * mm, 10 * m - 1 * mm, 7 * cm - 1 * mm);
842 coatWall->SetLineColor(kRed);
843 tShield->AddNode(coatWall, 1, new TGeoTranslation(0, 0, zEndOfAbsorb + 2*absorber_half_length + absorber_offset+7 * cm));
844
845 }
846 std::array<double, 9> fieldScale = {{1., 1., 1., 1., 1., 1., 1., 1., 1.}};
847 if (fWithCoMagnet > 0)
848 {
849 Double_t lengthSum = 0.;
850 for (int i = 2; i < 9; ++i)
851 {
852 lengthSum += dZf[i];
853 }
854 fieldScale.fill((fField * lengthSum - 2.2 * dZf[fWithCoMagnet])/fField/(lengthSum - dZf[fWithCoMagnet]));
855 fieldScale[0] = 1.;
856 fieldScale[1] = 1.;
857 try
858 {
859 fieldScale.at(fWithCoMagnet) = 2.2 / fField;
860 }
861 catch(const std::out_of_range& e)
862 {
863 Fatal( "ShipMuonShield", "Exception out of range for --coMuonShield occurred \n");
864 }
865 }
866 for (Int_t nM = 2; nM <= (nMagnets - 1); nM++) {
867 Double_t ironField_s = fField * fieldScale[nM] * tesla;
868 TGeoUniformMagField *magFieldIron_s = new TGeoUniformMagField(0.,ironField_s,0.);
869 TGeoUniformMagField *RetField_s = new TGeoUniformMagField(0.,-ironField_s,0.);
870 TGeoUniformMagField *ConRField_s = new TGeoUniformMagField(-ironField_s,0.,0.);
871 TGeoUniformMagField *ConLField_s = new TGeoUniformMagField(ironField_s,0.,0.);
872 TGeoUniformMagField *fields_s[4] = {magFieldIron_s,RetField_s,ConRField_s,ConLField_s};
873 CreateMagnet(magnetName[nM], iron, tShield, fields_s, fieldDirection[nM],
874 dXIn[nM], dYIn[nM], dXOut[nM], dYOut[nM], dZf[nM],
875 midGapIn[nM], midGapOut[nM], HmainSideMagIn[nM],
876 HmainSideMagOut[nM], gapIn[nM], gapOut[nM], Z[nM], nM==8, fStepGeo);
877
878 if (nM==8 || !fSupport) continue;
879 Double_t dymax = std::max(dYIn[nM] + dXIn[nM], dYOut[nM] + dXOut[nM]);
880 Double_t dymin = std::min(dYIn[nM] + dXIn[nM], dYOut[nM] + dXOut[nM]);
881 Double_t slope =
882 (dYIn[nM] + dXIn[nM] - dYOut[nM] - dXOut[nM]) / (2 * dZf[nM]);
883 Double_t w1 = 2 * dXIn[nM] + std::max(20., gapIn[nM]);
884 Double_t w2 = 2 * dXOut[nM] + std::max(20., gapOut[nM]);
885 Double_t anti_overlap = 0.1;
886 Double_t h1 = 0.5 * (dYIn[nM] + dXIn[nM] + anti_overlap - 10 * m + fFloor);
887 Double_t h2 = 0.5 * (dYOut[nM] + dXOut[nM] + anti_overlap - 10 * m + fFloor);
888 Double_t length = std::min(0.5 * m, std::abs(dZf[nM]/2. - 5 * cm));
889 std::array<Double_t, 16> verticesIn = {
890 -w1, -h1,
891 +w1, -h1,
892 +w1, +h1,
893 -w1, +h1,
894 -w1, -h1 + slope * 2. * length,
895 +w1, -h1 + slope * 2. * length,
896 +w1, +h1,
897 -w1, +h1,
898 };
899 std::array<Double_t, 16> verticesOut = {
900 -w2, -h2 - slope * 2. * length,
901 +w2, -h2 - slope * 2. * length,
902 +w2, +h2,
903 -w2, +h2,
904 -w2, -h2,
905 +w2, -h2,
906 +w2, +h2,
907 -w2, +h2,
908 };
909 if (!fStepGeo)
910 {
911
912
913 TGeoVolume *pillar1 =
914 gGeoManager->MakeArb8(TString::Format("pillar_%d", 2 * nM - 1),
915 steel, length, verticesIn.data());
916 TGeoVolume *pillar2 =
917 gGeoManager->MakeArb8(TString::Format("pillar_%d", 2 * nM), steel,
918 length, verticesOut.data());
919 pillar1->SetLineColor(kGreen-5);
920 pillar2->SetLineColor(kGreen-5);
921 tShield->AddNode(pillar1, 1, new TGeoTranslation(
922 0, -0.5 * (dYIn[nM] + dXIn[nM] + 10 * m - fFloor),
923 Z[nM] - dZf[nM] + length));
924 tShield->AddNode(pillar2, 1, new TGeoTranslation(
925 0, -0.5 * (dYOut[nM] + dXOut[nM] + 10 * m - fFloor),
926 Z[nM] + dZf[nM] - length));
927 }
928 }
929
930 } else {
931 CreateTube("AbsorberAdd", iron, 15, 400, dZ0, 43, tShield, 0, 0, zEndOfAbsorb - dZ0);
932 CreateTube("AbsorberAddCore", iron, 0, 15, dZ0, 38, tShield, 0, 0, zEndOfAbsorb - dZ0);
933
934 for (Int_t nM = 0; nM < (nMagnets - 1); nM++) {
935 CreateMagnet(magnetName[nM],iron,tShield,fields,fieldDirection[nM],
936 dXIn[nM],dYIn[nM],dXOut[nM],dYOut[nM],dZf[nM],
937 midGapIn[nM],midGapOut[nM],HmainSideMagIn[nM],HmainSideMagOut[nM],
938 gapIn[nM],gapOut[nM],Z[nM],0, fStepGeo);
939 }
940 }
941 Double_t dX1 = dXIn[0];
942 Double_t dY = dYIn[0];
943
944 // Place in origin of SHiP coordinate system as subnodes placed correctly
945 top->AddNode(tShield, 1);
946
947// Concrete around first magnets. i.e. Tunnel
948 Double_t dZ = dZ1 + dZ2;
949 Double_t ZT = zEndOfAbsorb + dZ;
950 TGeoBBox *box1 = new TGeoBBox("box1", 10*m,10*m,dZ);
951 TGeoBBox *box2 = new TGeoBBox("box2", 15*m,15*m,dZ);
952 TGeoCompositeShape *compRockS = new TGeoCompositeShape("compRockS", "box2-box1");
953 TGeoVolume *rockS = new TGeoVolume("rockS", compRockS, concrete);
954 rockS->SetLineColor(11); // grey
955 rockS->SetTransparency(50);
956 top->AddNode(rockS, 1, new TGeoTranslation(0, 0, ZT ));
957// Concrete around decay tunnel
958 Double_t dZD = 100*m + fMuonShieldLength;
959 TGeoBBox *box3 = new TGeoBBox("box3", 15*m, 15*m,dZD/2.);
960 TGeoBBox *box4 = new TGeoBBox("box4", 10*m, 10*m,dZD/2.);
961
962 if (fDesign >= 7 && fFloor > 0) {
963 // Only add floor for new shield
964 TGeoBBox *box5 = new TGeoBBox("shield_floor", 10 * m, fFloor / 2.,
965 fMuonShieldLength / 2. - dZ - 14.2 * cm); // substract CoatWall
966 TGeoVolume *floor = new TGeoVolume("floorM", box5, concrete);
967 floor->SetLineColor(11); // grey
968 top->AddNode(floor, 1, new TGeoTranslation(0, -10 * m + fFloor / 2.,
969 zEndOfAbsorb +
970 fMuonShieldLength / 2. +dZ + 14.2 * cm)); // avoiding overlap with CoatWall
971 }
972 TGeoCompositeShape *compRockD =
973 new TGeoCompositeShape("compRockD", "(box3-box4)");
974 TGeoVolume *rockD = new TGeoVolume("rockD", compRockD, concrete);
975 rockD->SetLineColor(11); // grey
976 rockD->SetTransparency(50);
977 top->AddNode(rockD, 1, new TGeoTranslation(0, 0, zEndOfAbsorb + 2*dZ + dZD/2.));
978//
979 } else {
980 Fatal("ShipMuonShield","Design %i does not match implemented designs",fDesign);
981 }
982}
983ClassImp(ShipMuonShield)
tesla
Double_t tesla
Definition ShipBellField.cxx:17
cm
Double_t cm
Definition ShipBellField.cxx:19
m
Double_t m
Definition ShipBellField.cxx:20
mm
Double_t mm
Definition ShipBellField.cxx:21
tesla
Double_t tesla
Definition ShipMuonShield.cxx:24
cm
Double_t cm
Definition ShipMuonShield.cxx:20
m
Double_t m
Definition ShipMuonShield.cxx:21
kilogauss
Double_t kilogauss
Definition ShipMuonShield.cxx:23
mm
Double_t mm
Definition ShipMuonShield.cxx:22
FieldDirection
FieldDirection
Definition ShipMuonShield.h:16
ShipMuonShield::~ShipMuonShield
virtual ~ShipMuonShield()
Definition ShipMuonShield.cxx:26
ShipMuonShield::Initialize
Int_t Initialize(std::vector< TString > &magnetName, std::vector< FieldDirection > &fieldDirection, std::vector< Double_t > &dXIn, std::vector< Double_t > &dYIn, std::vector< Double_t > &dXOut, std::vector< Double_t > &dYOut, std::vector< Double_t > &dZ, std::vector< Double_t > &midGapIn, std::vector< Double_t > &midGapOut, std::vector< Double_t > &HmainSideMagIn, std::vector< Double_t > &HmainSideMagOut, std::vector< Double_t > &gapIn, std::vector< Double_t > &gapOut, std::vector< Double_t > &Z)
Definition ShipMuonShield.cxx:411
ShipMuonShield::fMuonShieldLength
Double_t fMuonShieldLength
Definition ShipMuonShield.h:44
ShipMuonShield::zEndOfAbsorb
Double_t zEndOfAbsorb
Definition ShipMuonShield.h:47
ShipMuonShield::CreateMagnet
void CreateMagnet(TString magnetName, TGeoMedium *medium, TGeoVolume *tShield, TGeoUniformMagField *fields[4], FieldDirection fieldDirection, Double_t dX, Double_t dY, Double_t dX2, Double_t dY2, Double_t dZ, Double_t middleGap, Double_t middleGap2, Double_t HmainSideMag, Double_t HmainSideMag2, Double_t gap, Double_t gap2, Double_t Z, Bool_t NotMagnet, Bool_t stepGeo)
Definition ShipMuonShield.cxx:235
ShipMuonShield::fWithConstShieldField
Bool_t fWithConstShieldField
Definition ShipMuonShield.h:53
ShipMuonShield::CreateArb8
void CreateArb8(TString arbName, TGeoMedium *medium, Double_t dZ, std::array< Double_t, 16 > corners, Int_t color, TGeoUniformMagField *magField, TGeoVolume *top, Double_t x_translation, Double_t y_translation, Double_t z_translation)
Definition ShipMuonShield.cxx:149
ShipMuonShield::InitMedium
Int_t InitMedium(TString name)
Definition ShipMuonShield.cxx:120
ShipMuonShield::CreateTube
void CreateTube(TString tubeName, TGeoMedium *medium, Double_t dX, Double_t dY, Double_t dZ, Int_t color, TGeoVolume *top, Double_t x_translation, Double_t y_translation, Double_t z_translation)
Definition ShipMuonShield.cxx:137
ShipMuonShield::fWithConstAbsorberField
Bool_t fWithConstAbsorberField
Definition ShipMuonShield.h:52
ClassImp
ClassImp(ecalContFact) ecalContFact
Definition ecalContFact.cxx:8