generate_dunevd10kt_2view_v2_refactored.pl
Go to the documentation of this file.
1 #!/usr/bin/perl
2 
3 #
4 #
5 # First attempt to make a GDML fragment generator for the DUNE vertical drift
6 # 10kt detector geometry with 2 orthogonal view readout
7 # The lower chamber is not added yet.
8 # !!!NOTE!!!: the readout is on a positive Y plane (drift along horizontal X)
9 # due to current reco limitations)
10 # No photon detectors declared
11 # Simplified treatment of inter-module dead spaces
12 #
13 # Created: Thu Oct 1 16:45:27 CEST 2020
14 # Vyacheslav Galymov <vgalymov@ipnl.in2p3.fr>
15 #
16 # Modified:
17 # VG: Added defs to enable use in the refactored sim framework
18 # Laura Paulucci (lpaulucc@fnal.gov): Sept 2021 PDS added.
19 # Use option -pds=1 for backup design (membrane only coverage).
20 # Default (pds=0) is the reference design (~4-pi).
21 # This is linked with a larger geometry to account for photon propagation, generate it with -k=4.
22 # Field Cage is turned on with reference and backup designs to match PDS option.
23 # For not including the pds, please use option -pds=-1
24 #
25 #
26 #################################################################################
27 
28 # Each subroutine generates a fragment GDML file, and the last subroutine
29 # creates an XML file that make_gdml.pl will use to appropriately arrange
30 # the fragment GDML files to create the final desired DUNE GDML file,
31 # to be named by make_gdml output command
32 
33 ##################################################################################
34 
35 
36 #use warnings;
37 use gdmlMaterials;
38 use Math::Trig;
39 use Getopt::Long;
40 use Math::BigFloat;
41 Math::BigFloat->precision(-16);
42 
43 GetOptions( "help|h" => \$help,
44  "suffix|s:s" => \$suffix,
45  "output|o:s" => \$output,
46  "wires|w:s" => \$wires,
47  "workspace|k:s" => \$wkspc,
48  "pdsconfig|pds:s" => \$pdsconfig);
49 
50 my $FieldCage_switch="off";
51 my $Cathode_switch="off";
52 
53 if ( defined $help )
54 {
55  # If the user requested help, print the usage notes and exit.
56  usage();
57  exit;
58 }
59 
60 if ( ! defined $suffix )
61 {
62  # The user didn't supply a suffix, so append nothing to the file
63  # names.
64  $suffix = "";
65 }
66 else
67 {
68  # Otherwise, stick a "-" before the suffix, so that a suffix of
69  # "test" applied to filename.gdml becomes "filename-test.gdml".
70  $suffix = "-" . $suffix;
71 }
72 
73 
74 $workspace = 0;
75 if(defined $wkspc )
76 {
77  $workspace = $wkspc;
78 }
79 elsif ( $workspace != 0 )
80 {
81  print "\t\tCreating smaller workspace geometry.\n";
82 }
83 
84 if ( ! defined $pdsconfig )
85 {
86  $pdsconfig = 0;
87  print "\t\tCreating reference design: 4-pi PDS converage.\n";
88 }
89 elsif ( $pdsconfig == 1 )
90 {
91  print "\t\tCreating backup design: membrane-only PDS coverage.\n";
92 }
93 
94 # set wires on to be the default, unless given an input by the user
95 $wires_on = 1; # 1=on, 0=off
96 if (defined $wires)
97 {
98  $wires_on = $wires;
99 }
100 
101 $tpc_on = 1;
102 $basename="_";
103 
104 
105 ##################################################################
106 ############## Parameters for One Readout Panel ##################
107 
108 # parameters for 1.5 x 1.7 sub-unit Charge Readout Module / Unit
109 #$widthPCBActive = 169.0; # cm
110 #$lengthPCBActive = 150.0; # cm
111 $nChannelsViewInd = 320;
112 $nChannelsViewCol = 288;
113 
114 $wirePitchY = 0.527; # $widthPCBActive / $nChannelsViewInd;
115 $wirePitchZ = 0.518; # $lengthPCBActive / $nChannelsViewCol;
116 
117 $widthPCBActive = $wirePitchY * $nChannelsViewInd;
118 $lengthPCBActive = $wirePitchZ * $nChannelsViewCol;
119 
120 $borderCRM = 0.05; # dead space at the border of each CRM
121 
122 $widthCRM_active = $widthPCBActive;
123 $lengthCRM_active = $lengthPCBActive;
124 
125 $widthCRM = $widthPCBActive + 2 * $borderCRM;
126 $lengthCRM = $lengthPCBActive + 2 * $borderCRM;
127 
128 # number of CRMs in y and z
129 $nCRM_x = 4 * 2;
130 $nCRM_z = 20 * 2;
131 
132 # create a smaller geometry
133 if( $workspace == 1 )
134 {
135  $nCRM_x = 1 * 2;
136  $nCRM_z = 1 * 2;
137 }
138 
139 # create a smaller geometry
140 if( $workspace == 2 )
141 {
142  $nCRM_x = 2 * 2;
143  $nCRM_z = 2 * 2;
144 }
145 
146 # create a smaller geometry
147 if( $workspace == 3 )
148 {
149  $nCRM_x = 3 * 2;
150  $nCRM_z = 3 * 2;
151 }
152 
153 # create pds geometry
154 if( $workspace == 4 )
155 {
156  $nCRM_x = 4 * 2;
157  $nCRM_z = 7 * 2;
158 }
159 
160 
161 # calculate tpc area based on number of CRMs and their dimensions
162 $widthTPCActive = $nCRM_x * $widthCRM; # around 1200 for full module
163 $lengthTPCActive = $nCRM_z * $lengthCRM; # around 6000 for full module
164 
165 # active volume dimensions
166 $driftTPCActive = 650.0;
167 
168 # model anode strips as wires of some diameter
169 $padWidth = 0.02;
170 $ReadoutPlane = 2 * $padWidth; ## + 0.5; # 5 mm thick PCB?
171 
172 ##################################################################
173 ############## Parameters for TPC and inner volume ###############
174 
175 # inner volume dimensions of the cryostat
176 $Argon_x = 1510;
177 $Argon_y = 1510;
178 $Argon_z = 6200;
179 
180 # width of gas argon layer on top
181 $HeightGaseousAr = 100;
182 
183 if( $workspace != 0 )
184 {
185  #active tpc + 1.0 m buffer on each side
186  $Argon_x = $driftTPCActive + $HeightGaseousAr + $ReadoutPlane + 100;
187  $Argon_y = $widthTPCActive + 200;
188  $Argon_z = $lengthTPCActive + 200;
189 }
190 
191 
192 # size of liquid argon buffer
193 $xLArBuffer = $Argon_x - $driftTPCActive - $HeightGaseousAr - $ReadoutPlane;
194 $yLArBuffer = 0.5 * ($Argon_y - $widthTPCActive);
195 $zLArBuffer = 0.5 * ($Argon_z - $lengthTPCActive);
196 
197 # cryostat
198 $SteelThickness = 0.12; # membrane
199 
200 $Cryostat_x = $Argon_x + 2*$SteelThickness;
201 $Cryostat_y = $Argon_y + 2*$SteelThickness;
202 $Cryostat_z = $Argon_z + 2*$SteelThickness;
203 
204 ##################################################################
205 ############## DetEnc and World relevant parameters #############
206 
207 $SteelSupport_x = 100;
208 $SteelSupport_y = 100;
209 $SteelSupport_z = 100;
210 $FoamPadding = 80;
211 $FracMassOfSteel = 0.5; #The steel support is not a solid block, but a mixture of air and steel
212 $FracMassOfAir = 1 - $FracMassOfSteel;
213 
214 
215 $SpaceSteelSupportToWall = 100;
216 $SpaceSteelSupportToCeiling = 100;
217 
218 $DetEncX = $Cryostat_x
219  + 2*($SteelSupport_x + $FoamPadding) + $SpaceSteelSupportToCeiling;
220 
221 $DetEncY = $Cryostat_y
222  + 2*($SteelSupport_y + $FoamPadding) + 2*$SpaceSteelSupportToWall;
223 
224 $DetEncZ = $Cryostat_z
225  + 2*($SteelSupport_z + $FoamPadding) + 2*$SpaceSteelSupportToWall;
226 
227 $posCryoInDetEnc_x = - $DetEncX/2 + $SteelSupport_x + $FoamPadding + $Cryostat_x/2;
228 
229 
230 $RockThickness = 4000;
231 
232  # We want the world origin to be vertically centered on active TPC
233  # This is to be added to the x and y position of every volume in volWorld
234 
235 $OriginXSet = $DetEncX/2.0
236  - $SteelSupport_x
237  - $FoamPadding
238  - $SteelThickness
239  - $xLArBuffer
240  - $driftTPCActive/2.0;
241 
242 $OriginYSet = $DetEncY/2.0
243  - $SpaceSteelSupportToWall
244  - $SteelSupport_y
245  - $FoamPadding
246  - $SteelThickness
247  - $yLArBuffer
248  - $widthTPCActive/2.0;
249 
250  # We want the world origin to be at the very front of the fiducial volume.
251  # move it to the front of the enclosure, then back it up through the concrete/foam,
252  # then through the Cryostat shell, then through the upstream dead LAr (including the
253  # dead LAr on the edge of the TPC)
254  # This is to be added to the z position of every volume in volWorld
255 
256 $OriginZSet = $DetEncZ/2.0
257  - $SpaceSteelSupportToWall
258  - $SteelSupport_z
259  - $FoamPadding
260  - $SteelThickness
261  - $zLArBuffer
262  - $borderCRM;
263 
264 ##################################################################
265 ############## Field Cage Parameters ###############
266 $FieldShaperLongTubeLength = $lengthTPCActive;
267 $FieldShaperShortTubeLength = $widthTPCActive;
268 #$FieldShaperInnerRadius = 1.485;
269 #$FieldShaperOuterRadius = 1.685;
270 #$FieldShaperTorRad = 1.69;
271 $FieldShaperInnerRadius = 0.5; #cm
272 $FieldShaperOuterRadius = 2.285; #cm
273 $FieldShaperOuterRadiusSlim = 0.75; #cm
274 $FieldShaperTorRad = 2.3; #cm
275 
276 $FieldShaperLength = $FieldShaperLongTubeLength + 2*$FieldShaperOuterRadius+ 2*$FieldShaperTorRad;
277 $FieldShaperWidth = $FieldShaperShortTubeLength + 2*$FieldShaperOuterRadius+ 2*$FieldShaperTorRad;
278 
279 $FieldShaperSeparation = 6.0; #cm
280 $NFieldShapers = ($driftTPCActive/$FieldShaperSeparation) - 1;
281 
282 $FieldCageSizeX = $FieldShaperSeparation*$NFieldShapers+2;
283 $FieldCageSizeY = $FieldShaperWidth+2;
284 $FieldCageSizeZ = $FieldShaperLength+2;
285 
286 
287 ####################################################################
288 ######################## ARAPUCA Dimensions ########################
289 ## in cm
290 
291 $ArapucaOut_x = 65.0;
292 $ArapucaOut_y = 2.5;
293 $ArapucaOut_z = 65.0;
294 $ArapucaIn_x = 60.0;
295 $ArapucaIn_y = 2.0;
296 $ArapucaIn_z = 60.0;
297 $ArapucaAcceptanceWindow_x = 60.0;
298 $ArapucaAcceptanceWindow_y = 1.0;
299 $ArapucaAcceptanceWindow_z = 60.0;
300 $BlockPD_x = 0.5*$widthCRM_active - 6.35; #Sub-division of Frame (4 n x, 4 in z)
301 $BlockPD_z = 0.5*$lengthCRM_active - 6.25;
302 $GapPD = 0.5; #Size of supporting bars in Frame
303 $FrameToArapucaSpace = 1.0; #At this moment, should cover the thickness of Frame + small gap to prevent overlap. VALUE NEEDS TO BE CHECKED!!!
304 $FrameToArapucaSpaceLat = $yLArBuffer - 60.0; #Arapucas 60 cm behind FC. At this moment, should cover the thickness of Frame + small gap to prevent overlap. VALUE NEEDS TO BE CHECKED!!!
305 $VerticalPDdist = 75.0; #distance of arapucas (center to center) in the y direction
306 $HorizontalPDdist = 150.0; #distance of arapucas (center to center) in the x direction
307 
308 #Positions of the 4 arapucas with respect to the Frame center --> arapucas over the cathode
309 $list_posx_bot[0]=-2*$BlockPD_x - 1.5*$GapPD + 0.5*$ArapucaOut_x;
310 $list_posz_bot[0]= 0.5*$BlockPD_z + $GapPD;
311 $list_posx_bot[1]=- $GapPD - 0.5*$ArapucaOut_x;
312 $list_posz_bot[1]=-1.5*$BlockPD_z - 1.5*$GapPD;
313 $list_posx_bot[2]=-$list_posx_bot[1];
314 $list_posz_bot[2]=-$list_posz_bot[1];
315 $list_posx_bot[3]=-$list_posx_bot[0];
316 $list_posz_bot[3]=-$list_posz_bot[0];
317 
318 
319 
320 #+++++++++++++++++++++++++ End defining variables ++++++++++++++++++++++++++
321 
322 
323 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
324 #+++++++++++++++++++++++++++++++++++++++++ usage +++++++++++++++++++++++++++++++++++++++++
325 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
326 
327 sub usage()
328 {
329  print "Usage: $0 [-h|--help] [-o|--output <fragments-file>] [-s|--suffix <string>]\n";
330  print " if -o is omitted, output goes to STDOUT; <fragments-file> is input to make_gdml.pl\n";
331  print " -s <string> appends the string to the file names; useful for multiple detector versions\n";
332  print " -h prints this message, then quits\n";
333 }
334 
335 
336 sub gen_Extend()
337 {
338 
339 # Create the <define> fragment file name,
340 # add file to list of fragments,
341 # and open it
342  $DEF = $basename."_Ext" . $suffix . ".gdml";
343  push (@gdmlFiles, $DEF);
344  $DEF = ">" . $DEF;
345  open(DEF) or die("Could not open file $DEF for writing");
346 
347 print DEF <<EOF;
348 <?xml version='1.0'?>
349 <gdml>
350 <extension>
351  <color name="magenta" R="0.0" G="1.0" B="0.0" A="1.0" />
352  <color name="green" R="0.0" G="1.0" B="0.0" A="1.0" />
353  <color name="red" R="1.0" G="0.0" B="0.0" A="1.0" />
354  <color name="blue" R="0.0" G="0.0" B="1.0" A="1.0" />
355  <color name="yellow" R="1.0" G="1.0" B="0.0" A="1.0" />
356 </extension>
357 </gdml>
358 EOF
359  close (DEF);
360 }
361 
362 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
363 #++++++++++++++++++++++++++++++++++++++ gen_Define +++++++++++++++++++++++++++++++++++++++
364 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
365 
366 sub gen_Define()
367 {
368 
369 # Create the <define> fragment file name,
370 # add file to list of fragments,
371 # and open it
372  $DEF = $basename."_Def" . $suffix . ".gdml";
373  push (@gdmlFiles, $DEF);
374  $DEF = ">" . $DEF;
375  open(DEF) or die("Could not open file $DEF for writing");
376 
377 
378 print DEF <<EOF;
379 <?xml version='1.0'?>
380 <gdml>
381 <define>
382 
383 <!--
384 
385 
386 
387 -->
388 
389  <position name="posCryoInDetEnc" unit="cm" x="$posCryoInDetEnc_x" y="0" z="0"/>
390  <position name="posCenter" unit="cm" x="0" y="0" z="0"/>
391  <rotation name="rPlus90AboutX" unit="deg" x="90" y="0" z="0"/>
392  <rotation name="rPlus90AboutY" unit="deg" x="0" y="90" z="0"/>
393  <rotation name="rPlus90AboutXPlus90AboutY" unit="deg" x="90" y="90" z="0"/>
394  <rotation name="rMinus90AboutX" unit="deg" x="270" y="0" z="0"/>
395  <rotation name="rMinus90AboutY" unit="deg" x="0" y="270" z="0"/>
396  <rotation name="rMinus90AboutYMinus90AboutX" unit="deg" x="270" y="270" z="0"/>
397  <rotation name="rPlus180AboutX" unit="deg" x="180" y="0" z="0"/>
398  <rotation name="rPlus180AboutY" unit="deg" x="0" y="180" z="0"/>
399  <rotation name="rPlus180AboutXPlus180AboutY" unit="deg" x="180" y="180" z="0"/>
400  <rotation name="rIdentity" unit="deg" x="0" y="0" z="0"/>
401 </define>
402 </gdml>
403 EOF
404  close (DEF);
405 }
406 
407 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
408 #+++++++++++++++++++++++++++++++++++++ gen_Materials +++++++++++++++++++++++++++++++++++++
409 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
410 
411 sub gen_Materials()
412 {
413 
414 # Create the <materials> fragment file name,
415 # add file to list of output GDML fragments,
416 # and open it
417  $MAT = $basename."_Materials" . $suffix . ".gdml";
418  push (@gdmlFiles, $MAT);
419  $MAT = ">" . $MAT;
420 
421  open(MAT) or die("Could not open file $MAT for writing");
422 
423  # Add any materials special to this geometry by defining a mulitline string
424  # and passing it to the gdmlMaterials::gen_Materials() function.
425 my $asmix = <<EOF;
426  <!-- preliminary values -->
427  <material name="AirSteelMixture" formula="AirSteelMixture">
428  <D value=" 0.001205*(1-$FracMassOfSteel) + 7.9300*$FracMassOfSteel " unit="g/cm3"/>
429  <fraction n="$FracMassOfSteel" ref="STEEL_STAINLESS_Fe7Cr2Ni"/>
430  <fraction n="$FracMassOfAir" ref="Air"/>
431  </material>
432  <material name="vm2000" formula="vm2000">
433  <D value="1.2" unit="g/cm3"/>
434  <composite n="2" ref="carbon"/>
435  <composite n="4" ref="hydrogen"/>
436  </material>
437 EOF
438 
439  # add the general materials used anywere
440  print MAT gdmlMaterials::gen_Materials( $asmix );
441 
442  close(MAT);
443 }
444 
445 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
446 #++++++++++++++++++++++++++++++++++++++++ gen_TPC ++++++++++++++++++++++++++++++++++++++++
447 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
448 sub gen_TPC()
449 {
450  # CRM active volume
451  my $TPCActive_x = $driftTPCActive;
452  my $TPCActive_y = $widthCRM_active;
453  my $TPCActive_z = $lengthCRM_active;
454 
455  # CRM total volume
456  my $TPC_x = $TPCActive_x + $ReadoutPlane;
457  my $TPC_y = $widthCRM;
458  my $TPC_z = $lengthCRM;
459 
460 
461  $TPC = $basename."_TPC" . $suffix . ".gdml";
462  push (@gdmlFiles, $TPC);
463  $TPC = ">" . $TPC;
464  open(TPC) or die("Could not open file $TPC for writing");
465 
466 # The standard XML prefix and starting the gdml
467  print TPC <<EOF;
468 <?xml version='1.0'?>
469 <gdml>
470 EOF
471 
472 
473  # All the TPC solids save the wires.
474  print TPC <<EOF;
475 <solids>
476  <box name="CRM" lunit="cm"
477  x="$TPC_x"
478  y="$TPC_y"
479  z="$TPC_z"/>
480  <box name="CRMVPlane" lunit="cm"
481  x="$padWidth"
482  y="$TPCActive_y"
483  z="$TPCActive_z"/>
484  <box name="CRMZPlane" lunit="cm"
485  x="$padWidth"
486  y="$TPCActive_y"
487  z="$TPCActive_z"/>
488  <box name="CRMActive" lunit="cm"
489  x="$TPCActive_x"
490  y="$TPCActive_y"
491  z="$TPCActive_z"/>
492 EOF
493 
494 
495 #++++++++++++++++++++++++++++ Wire Solids ++++++++++++++++++++++++++++++
496 
497 print TPC <<EOF;
498 
499  <tube name="CRMWireV"
500  rmax="0.5*$padWidth"
501  z="$TPCActive_z"
502  deltaphi="360"
503  aunit="deg"
504  lunit="cm"/>
505  <tube name="CRMWireZ"
506  rmax="0.5*$padWidth"
507  z="$TPCActive_y"
508  deltaphi="360"
509  aunit="deg"
510  lunit="cm"/>
511 </solids>
512 
513 EOF
514 
515 
516 # Begin structure and create wire logical volumes
517 print TPC <<EOF;
518 <structure>
519  <volume name="volTPCActive">
520  <materialref ref="LAr"/>
521  <solidref ref="CRMActive"/>
522  <auxiliary auxtype="SensDet" auxvalue="SimEnergyDeposit"/>
523  <auxiliary auxtype="StepLimit" auxunit="cm" auxvalue="0.5*cm"/>
524  <auxiliary auxtype="Efield" auxunit="V/cm" auxvalue="500*V/cm"/>
525  <colorref ref="blue"/>
526  </volume>
527 EOF
528 
529 if ($wires_on==1)
530 {
531  print TPC <<EOF;
532  <volume name="volTPCWireV">
533  <materialref ref="Copper_Beryllium_alloy25"/>
534  <solidref ref="CRMWireV"/>
535  </volume>
536 
537  <volume name="volTPCWireZ">
538  <materialref ref="Copper_Beryllium_alloy25"/>
539  <solidref ref="CRMWireZ"/>
540  </volume>
541 EOF
542 }
543 
544 print TPC <<EOF;
545 
546  <volume name="volTPCPlaneV">
547  <materialref ref="LAr"/>
548  <solidref ref="CRMVPlane"/>
549 EOF
550 
551 if ($wires_on==1) # add wires to Y plane (plane with wires reading y position)
552 {
553 for($i=0;$i<$nChannelsViewInd;++$i)
554 {
555 my $ypos = -0.5 * $TPCActive_y + ($i+0.5)*$wirePitchY + 0.5*$padWidth;
556 
557 print TPC <<EOF;
558  <physvol>
559  <volumeref ref="volTPCWireV"/>
560  <position name="posWireV$i" unit="cm" x="0" y="$ypos" z="0"/>
561  <rotationref ref="rIdentity"/>
562  </physvol>
563 EOF
564 }
565 }
566 
567 print TPC <<EOF;
568  </volume>
569 
570  <volume name="volTPCPlaneZ">
571  <materialref ref="LAr"/>
572  <solidref ref="CRMZPlane"/>
573 EOF
574 
575 
576 if ($wires_on==1) # add wires to Z plane (plane with wires reading z position)
577 {
578 for($i=0;$i<$nChannelsViewCol;++$i)
579 {
580 
581 my $zpos = -0.5 * $TPCActive_z + ($i+0.5)*$wirePitchZ + 0.5*$padWidth;
582 print TPC <<EOF;
583  <physvol>
584  <volumeref ref="volTPCWireZ"/>
585  <position name="posWireZ$i" unit="cm" x="0" y="0" z="$zpos"/>
586  <rotationref ref="rPlus90AboutX"/>
587  </physvol>
588 EOF
589 }
590 }
591 
592 
593 print TPC <<EOF;
594  </volume>
595 EOF
596 
597 
598 $posVplane[0] = 0.5*$TPC_x - 1.5*$padWidth;
599 $posVplane[1] = 0;
600 $posVplane[2] = 0;
601 
602 $posZplane[0] = 0.5*$TPC_x - 0.5*$padWidth;
603 $posZplane[1] = 0;
604 $posZplane[2] = 0;
605 
606 $posTPCActive[0] = -$ReadoutPlane/2;
607 $posTPCActive[1] = 0;
608 $posTPCActive[2] = 0;
609 
610 
611 #wrap up the TPC file
612 print TPC <<EOF;
613 
614  <volume name="volTPC">
615  <materialref ref="LAr"/>
616  <solidref ref="CRM"/>
617  <physvol>
618  <volumeref ref="volTPCPlaneV"/>
619  <position name="posPlaneV" unit="cm"
620  x="$posVplane[0]" y="$posVplane[1]" z="$posVplane[2]"/>
621  <rotationref ref="rIdentity"/>
622  </physvol>
623  <physvol>
624  <volumeref ref="volTPCPlaneZ"/>
625  <position name="posPlaneZ" unit="cm"
626  x="$posZplane[0]" y="$posZplane[1]" z="$posZplane[2]"/>
627  <rotationref ref="rIdentity"/>
628  </physvol>
629  <physvol>
630  <volumeref ref="volTPCActive"/>
631  <position name="posActive" unit="cm"
632  x="$posTPCActive[0]" y="$posTPCActive[1]" z="$posTPCActive[2]"/>
633  <rotationref ref="rIdentity"/>
634  </physvol>
635  </volume>
636 EOF
637 
638 
639 print TPC <<EOF;
640 </structure>
641 </gdml>
642 EOF
643 
644 close(TPC);
645 }
646 
647 
648 
649 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
650 #++++++++++++++++++++++++++++++++++++++ gen_FieldCage ++++++++++++++++++++++++++++++++++++
651 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
652 
653 sub gen_FieldCage {
654 
655  $FieldCage = $basename."_FieldCage" . $suffix . ".gdml";
656  push (@gdmlFiles, $FieldCage);
657  $FieldCage = ">" . $FieldCage;
658  open(FieldCage) or die("Could not open file $FieldCage for writing");
659 
660 # The standard XML prefix and starting the gdml
661 print FieldCage <<EOF;
662  <?xml version='1.0'?>
663  <gdml>
664 EOF
665 # The printing solids used in the Field Cage
666 #print "lengthTPCActive : $lengthTPCActive \n";
667 #print "widthTPCActive : $widthTPCActive \n";
668 
669 
670 print FieldCage <<EOF;
671 <solids>
672  <torus name="FieldShaperCorner" rmin="$FieldShaperInnerRadius" rmax="$FieldShaperOuterRadius" rtor="$FieldShaperTorRad" deltaphi="90" startphi="0" aunit="deg" lunit="cm"/>
673  <tube name="FieldShaperLongtube" rmin="$FieldShaperInnerRadius" rmax="$FieldShaperOuterRadius" z="$FieldShaperLongTubeLength" deltaphi="360" startphi="0" aunit="deg" lunit="cm"/>
674  <tube name="FieldShaperLongtubeSlim" rmin="$FieldShaperInnerRadius" rmax="$FieldShaperOuterRadiusSlim" z="$FieldShaperLongTubeLength" deltaphi="360" startphi="0" aunit="deg" lunit="cm"/>
675  <tube name="FieldShaperShorttube" rmin="$FieldShaperInnerRadius" rmax="$FieldShaperOuterRadius" z="$FieldShaperShortTubeLength" deltaphi="360" startphi="0" aunit="deg" lunit="cm"/>
676 
677  <union name="FSunion1">
678  <first ref="FieldShaperLongtube"/>
679  <second ref="FieldShaperCorner"/>
680  <position name="esquinapos1" unit="cm" x="@{[-$FieldShaperTorRad]}" y="0" z="@{[0.5*$FieldShaperLongTubeLength]}"/>
681  <rotation name="rot1" unit="deg" x="90" y="0" z="0" />
682  </union>
683 
684  <union name="FSunion2">
685  <first ref="FSunion1"/>
686  <second ref="FieldShaperShorttube"/>
687  <position name="esquinapos2" unit="cm" x="@{[-0.5*$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[+0.5*$FieldShaperLongTubeLength+$FieldShaperTorRad]}"/>
688  <rotation name="rot2" unit="deg" x="0" y="90" z="0" />
689  </union>
690 
691  <union name="FSunion3">
692  <first ref="FSunion2"/>
693  <second ref="FieldShaperCorner"/>
694  <position name="esquinapos3" unit="cm" x="@{[-$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[0.5*$FieldShaperLongTubeLength]}"/>
695  <rotation name="rot3" unit="deg" x="90" y="270" z="0" />
696  </union>
697 
698  <union name="FSunion4">
699  <first ref="FSunion3"/>
700  <second ref="FieldShaperLongtube"/>
701  <position name="esquinapos4" unit="cm" x="@{[-$FieldShaperShortTubeLength-2*$FieldShaperTorRad]}" y="0" z="0"/>
702  </union>
703 
704  <union name="FSunion5">
705  <first ref="FSunion4"/>
706  <second ref="FieldShaperCorner"/>
707  <position name="esquinapos5" unit="cm" x="@{[-$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[-0.5*$FieldShaperLongTubeLength]}"/>
708  <rotation name="rot5" unit="deg" x="90" y="180" z="0" />
709  </union>
710 
711  <union name="FSunion6">
712  <first ref="FSunion5"/>
713  <second ref="FieldShaperShorttube"/>
714  <position name="esquinapos6" unit="cm" x="@{[-0.5*$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[-0.5*$FieldShaperLongTubeLength-$FieldShaperTorRad]}"/>
715  <rotation name="rot6" unit="deg" x="0" y="90" z="0" />
716  </union>
717 
718  <union name="FieldShaperSolid">
719  <first ref="FSunion6"/>
720  <second ref="FieldShaperCorner"/>
721  <position name="esquinapos7" unit="cm" x="@{[-$FieldShaperTorRad]}" y="0" z="@{[-0.5*$FieldShaperLongTubeLength]}"/>
722  <rotation name="rot7" unit="deg" x="90" y="90" z="0" />
723  </union>
724 
725  <union name="FSunionSlim1">
726  <first ref="FieldShaperLongtubeSlim"/>
727  <second ref="FieldShaperCorner"/>
728  <position name="esquinapos1" unit="cm" x="@{[-$FieldShaperTorRad]}" y="0" z="@{[0.5*$FieldShaperLongTubeLength]}"/>
729  <rotation name="rot1" unit="deg" x="90" y="0" z="0" />
730  </union>
731 
732  <union name="FSunionSlim2">
733  <first ref="FSunionSlim1"/>
734  <second ref="FieldShaperShorttube"/>
735  <position name="esquinapos2" unit="cm" x="@{[-0.5*$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[+0.5*$FieldShaperLongTubeLength+$FieldShaperTorRad]}"/>
736  <rotation name="rot2" unit="deg" x="0" y="90" z="0" />
737  </union>
738 
739  <union name="FSunionSlim3">
740  <first ref="FSunionSlim2"/>
741  <second ref="FieldShaperCorner"/>
742  <position name="esquinapos3" unit="cm" x="@{[-$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[0.5*$FieldShaperLongTubeLength]}"/>
743  <rotation name="rot3" unit="deg" x="90" y="270" z="0" />
744  </union>
745 
746  <union name="FSunionSlim4">
747  <first ref="FSunionSlim3"/>
748  <second ref="FieldShaperLongtubeSlim"/>
749  <position name="esquinapos4" unit="cm" x="@{[-$FieldShaperShortTubeLength-2*$FieldShaperTorRad]}" y="0" z="0"/>
750  </union>
751 
752  <union name="FSunionSlim5">
753  <first ref="FSunionSlim4"/>
754  <second ref="FieldShaperCorner"/>
755  <position name="esquinapos5" unit="cm" x="@{[-$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[-0.5*$FieldShaperLongTubeLength]}"/>
756  <rotation name="rot5" unit="deg" x="90" y="180" z="0" />
757  </union>
758 
759  <union name="FSunionSlim6">
760  <first ref="FSunionSlim5"/>
761  <second ref="FieldShaperShorttube"/>
762  <position name="esquinapos6" unit="cm" x="@{[-0.5*$FieldShaperShortTubeLength-$FieldShaperTorRad]}" y="0" z="@{[-0.5*$FieldShaperLongTubeLength-$FieldShaperTorRad]}"/>
763  <rotation name="rot6" unit="deg" x="0" y="90" z="0" />
764  </union>
765 
766  <union name="FieldShaperSolidSlim">
767  <first ref="FSunionSlim6"/>
768  <second ref="FieldShaperCorner"/>
769  <position name="esquinapos7" unit="cm" x="@{[-$FieldShaperTorRad]}" y="0" z="@{[-0.5*$FieldShaperLongTubeLength]}"/>
770  <rotation name="rot7" unit="deg" x="90" y="90" z="0" />
771  </union>
772 
773 </solids>
774 
775 EOF
776 
777 print FieldCage <<EOF;
778 
779 <structure>
780 <volume name="volFieldShaper">
781  <materialref ref="Al2O3"/>
782  <solidref ref="FieldShaperSolid"/>
783 </volume>
784 <volume name="volFieldShaperSlim">
785  <materialref ref="Al2O3"/>
786  <solidref ref="FieldShaperSolidSlim"/>
787 </volume>
788 
789 </structure>
790 
791 EOF
792 
793 print FieldCage <<EOF;
794 
795 </gdml>
796 EOF
797 close(FieldCage);
798 }
799 
800 
801 
802 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
803 #++++++++++++++++++++++++++++++++++++++ gen_Cryostat +++++++++++++++++++++++++++++++++++++
804 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
805 
806 sub gen_Cryostat()
807 {
808 
809 # Create the cryostat fragment file name,
810 # add file to list of output GDML fragments,
811 # and open it
812  $CRYO = $basename."_Cryostat" . $suffix . ".gdml";
813  push (@gdmlFiles, $CRYO);
814  $CRYO = ">" . $CRYO;
815  open(CRYO) or die("Could not open file $CRYO for writing");
816 
817 
818 # The standard XML prefix and starting the gdml
819  print CRYO <<EOF;
820 <?xml version='1.0'?>
821 <gdml>
822 EOF
823 
824 # All the cryostat solids.
825 # External active are two side volumes for generating light outside the field cage (no top or bottom buffers included)
826 print CRYO <<EOF;
827 <solids>
828  <box name="Cryostat" lunit="cm"
829  x="$Cryostat_x"
830  y="$Cryostat_y"
831  z="$Cryostat_z"/>
832 
833  <box name="ArgonInterior" lunit="cm"
834  x="$Argon_x"
835  y="$Argon_y"
836  z="$Argon_z"/>
837 
838  <box name="GaseousArgon" lunit="cm"
839  x="$HeightGaseousAr"
840  y="$Argon_y"
841  z="$Argon_z"/>
842 
843  <box name="ExternalAuxOut" lunit="cm"
844  x="$Argon_x - $xLArBuffer"
845  y="$Argon_y - 2*$ArapucaOut_y - 2*$FrameToArapucaSpaceLat"
846  z="$Argon_z"/>
847 
848  <box name="ExternalAuxIn" lunit="cm"
849  x="$Argon_x"
850  y="$FieldCageSizeY"
851  z="$Argon_z + 1"/>
852 
853  <subtraction name="ExternalActive">
854  <first ref="ExternalAuxOut"/>
855  <second ref="ExternalAuxIn"/>
856  </subtraction>
857 
858  <subtraction name="SteelShell">
859  <first ref="Cryostat"/>
860  <second ref="ArgonInterior"/>
861  </subtraction>
862 
863 </solids>
864 EOF
865 
866 #PDS
867 #Double sided detectors should only be included when both top and bottom volumes become available
868 #Optical sensitive volumes cannot be rotated because Larsoft cannot pick up the rotation when obtinaing the lengths needed for the semi-analytic model --> two acceptance windows for single sided lateral and cathode
869 print CRYO <<EOF;
870 <solids>
871  <box name="ArapucaOut" lunit="cm"
872  x="@{[$ArapucaOut_x]}"
873  y="@{[$ArapucaOut_y]}"
874  z="@{[$ArapucaOut_z]}"/>
875 
876  <box name="ArapucaIn" lunit="cm"
877  x="@{[$ArapucaIn_x]}"
878  y="@{[$ArapucaOut_y]}"
879  z="@{[$ArapucaIn_z]}"/>
880 
881  <subtraction name="ArapucaWalls">
882  <first ref="ArapucaOut"/>
883  <second ref="ArapucaIn"/>
884  <position name="posArapucaSub" x="0" y="@{[$ArapucaOut_y/2.0]}" z="0." unit="cm"/>
885  </subtraction>
886 
887  <box name="ArapucaAcceptanceWindow" lunit="cm"
888  x="@{[$ArapucaAcceptanceWindow_x]}"
889  y="@{[$ArapucaAcceptanceWindow_y]}"
890  z="@{[$ArapucaAcceptanceWindow_z]}"/>
891 
892  <box name="ArapucaDoubleIn" lunit="cm"
893  x="@{[$ArapucaIn_x]}"
894  y="@{[$ArapucaOut_y+1.0]}"
895  z="@{[$ArapucaIn_z]}"/>
896 
897  <subtraction name="ArapucaDoubleWalls">
898  <first ref="ArapucaOut"/>
899  <second ref="ArapucaDoubleIn"/>
900  <position name="posArapucaDoubleSub" x="0" y="0" z="0" unit="cm"/>
901  </subtraction>
902 
903  <box name="ArapucaDoubleAcceptanceWindow" lunit="cm"
904  x="@{[$ArapucaOut_y-0.02]}"
905  y="@{[$ArapucaAcceptanceWindow_x]}"
906  z="@{[$ArapucaAcceptanceWindow_z]}"/>
907 
908  <box name="ArapucaCathodeAcceptanceWindow" lunit="cm"
909  x="@{[$ArapucaAcceptanceWindow_y]}"
910  y="@{[$ArapucaAcceptanceWindow_x]}"
911  z="@{[$ArapucaAcceptanceWindow_z]}"/>
912 
913 </solids>
914 EOF
915 
916 # Cryostat structure
917 print CRYO <<EOF;
918 <structure>
919  <volume name="volSteelShell">
920  <materialref ref="STEEL_STAINLESS_Fe7Cr2Ni" />
921  <solidref ref="SteelShell" />
922  </volume>
923  <volume name="volGaseousArgon">
924  <materialref ref="ArGas"/>
925  <solidref ref="GaseousArgon"/>
926  </volume>
927  <volume name="volExternalActive">
928  <materialref ref="LAr"/>
929  <solidref ref="ExternalActive"/>
930  <auxiliary auxtype="SensDet" auxvalue="SimEnergyDeposit"/>
931  <auxiliary auxtype="StepLimit" auxunit="cm" auxvalue="0.5208*cm"/>
932  <auxiliary auxtype="Efield" auxunit="V/cm" auxvalue="0*V/cm"/>
933  <colorref ref="green"/>
934  </volume>
935 EOF
936 #including single sided arapucas over the cathode while there is only the top volume
937 #if double sided, use
938 # <volume name="volArapucaDouble_$i\-$j\-$p">
939 # <materialref ref="G10" />
940 # <solidref ref="ArapucaDoubleWalls" />
941 # </volume>
942 # <volume name="volOpDetSensitive_ArapucaDouble_$i\-$j\-$p">
943 # <materialref ref="LAr"/>
944 # <solidref ref="ArapucaDoubleAcceptanceWindow"/>
945 # </volume>
946 if ($pdsconfig == 0){ #4-pi PDS converage
947 for($i=0 ; $i<$nCRM_x/2 ; $i++){ #arapucas over the cathode
948 for($j=0 ; $j<$nCRM_z/2 ; $j++){
949 for($p=0 ; $p<4 ; $p++){
950  print CRYO <<EOF;
951  <volume name="volArapucaDouble_$i\-$j\-$p">
952  <materialref ref="G10" />
953  <solidref ref="ArapucaWalls" />
954  </volume>
955  <volume name="volOpDetSensitive_ArapucaDouble_$i\-$j\-$p">
956  <materialref ref="LAr"/>
957  <solidref ref="ArapucaCathodeAcceptanceWindow"/>
958  </volume>
959 EOF
960 }
961 }
962 }
963 }
964 
965 if ($nCRM_x==8){ #arapucas on the laterals
966 for($j=0 ; $j<$nCRM_z/2 ; $j++){
967 for($p=0 ; $p<8 ; $p++){
968  print CRYO <<EOF;
969  <volume name="volArapucaLat_$j\-$p">
970  <materialref ref="G10" />
971  <solidref ref="ArapucaWalls" />
972  </volume>
973  <volume name="volOpDetSensitive_ArapucaLat_$j\-$p">
974  <materialref ref="LAr"/>
975  <solidref ref="ArapucaAcceptanceWindow"/>
976  </volume>
977 EOF
978 }
979 }
980 }
981 
982 if ($pdsconfig == 1){ #Membrane PDS converage
983 if ($nCRM_x==8) { #arapucas on the laterals
984 for($j=0 ; $j<$nCRM_z/2 ; $j++){
985 for($p=8 ; $p<18 ; $p++){
986  print CRYO <<EOF;
987  <volume name="volArapucaLat_$j\-$p">
988  <materialref ref="G10" />
989  <solidref ref="ArapucaWalls" />
990  </volume>
991  <volume name="volOpDetSensitive_ArapucaLat_$j\-$p">
992  <materialref ref="LAr"/>
993  <solidref ref="ArapucaAcceptanceWindow"/>
994  </volume>
995 EOF
996 }
997 }
998 }
999 }
1000 
1001  print CRYO <<EOF;
1002 
1003  <volume name="volCryostat">
1004  <materialref ref="LAr" />
1005  <solidref ref="Cryostat" />
1006  <physvol>
1007  <volumeref ref="volGaseousArgon"/>
1008  <position name="posGaseousArgon" unit="cm" x="@{[$Argon_x/2-$HeightGaseousAr/2]}" y="0" z="0"/>
1009  </physvol>
1010  <physvol>
1011  <volumeref ref="volSteelShell"/>
1012  <position name="posSteelShell" unit="cm" x="0" y="0" z="0"/>
1013  </physvol>
1014  <physvol>
1015  <volumeref ref="volExternalActive"/>
1016  <position name="posExternalActive" unit="cm" x="-$xLArBuffer/2" y="0" z="0"/>
1017  </physvol>
1018 EOF
1019 
1020 
1021 if ($tpc_on==1) # place TPC inside croysotat
1022 {
1023  $posX = $Argon_x/2 - $HeightGaseousAr - 0.5*($driftTPCActive + $ReadoutPlane);
1024  $idx = 0;
1025  for($ii=0;$ii<$nCRM_z;$ii++)
1026  {
1027  $posZ = -0.5*$Argon_z + $zLArBuffer + ($ii+0.5)*$lengthCRM;
1028 
1029  for($jj=0;$jj<$nCRM_x;$jj++)
1030  {
1031  $posY = -0.5*$Argon_y + $yLArBuffer + ($jj+0.5)*$widthCRM;
1032  print CRYO <<EOF;
1033  <physvol>
1034  <volumeref ref="volTPC"/>
1035  <position name="posTPC\-$idx" unit="cm"
1036  x="$posX" y="$posY" z="$posZ"/>
1037  </physvol>
1038 EOF
1039 $idx++
1040  }
1041  }
1042 
1043 }
1044 
1045 #The +50 in the x positions must depend on some other parameter
1046  if ( $FieldCage_switch eq "on" ) {
1047  for ( $i=0; $i<$NFieldShapers; $i=$i+1 ) {
1048  $dist=$i*$FieldShaperSeparation;
1049 $posX = $Argon_x/2 - $HeightGaseousAr - 0.5*($driftTPCActive + $ReadoutPlane);
1050  if ($pdsconfig==0){
1051  if ($dist>250){
1052  print CRYO <<EOF;
1053  <physvol>
1054  <volumeref ref="volFieldShaperSlim"/>
1055  <position name="posFieldShaper$i" unit="cm" x="@{[-$OriginXSet+50+($i-$NFieldShapers*0.5)*$FieldShaperSeparation]}" y="@{[-0.5*$FieldShaperShortTubeLength-$FieldShaperTorRad]}" z="0" />
1056  <rotation name="rotFS$i" unit="deg" x="0" y="0" z="90" />
1057  </physvol>
1058 EOF
1059  }else{
1060  print CRYO <<EOF;
1061  <physvol>
1062  <volumeref ref="volFieldShaper"/>
1063  <position name="posFieldShaper$i" unit="cm" x="@{[-$OriginXSet+50+($i-$NFieldShapers*0.5)*$FieldShaperSeparation]}" y="@{[-0.5*$FieldShaperShortTubeLength-$FieldShaperTorRad]}" z="0" />
1064  <rotation name="rotFS$i" unit="deg" x="0" y="0" z="90" />
1065  </physvol>
1066 EOF
1067  }
1068  }else{
1069  print CRYO <<EOF;
1070  <physvol>
1071  <volumeref ref="volFieldShaperSlim"/>
1072  <position name="posFieldShaper$i" unit="cm" x="@{[-$OriginXSet+50+($i-$NFieldShapers*0.5)*$FieldShaperSeparation]}" y="@{[-0.5*$FieldShaperShortTubeLength-$FieldShaperTorRad]}" z="0" />
1073  <rotation name="rotFS$i" unit="deg" x="0" y="0" z="90" />
1074  </physvol>
1075 EOF
1076  }
1077  }
1078  }
1079 
1080 
1081 $CathodePosX =-$OriginXSet+50+(-1-$NFieldShapers*0.5)*$FieldShaperSeparation;
1082 $CathodePosY = 0;
1083 $CathodePosZ = 0;
1084  if ( $Cathode_switch eq "on" )
1085  {
1086  print CRYO <<EOF;
1087  <physvol>
1088  <volumeref ref="volGroundGrid"/>
1089  <position name="posGroundGrid01" unit="cm" x="$CathodePosX" y="@{[-$CathodePosY]}" z="@{[$CathodePosZ]}"/>
1090  <rotation name="rotGG01" unit="deg" x="0" y="0" z="90" />
1091  </physvol>
1092 
1093 EOF
1094  }
1095 
1096 $FrameLenght_x= 2.0*$widthCRM_active;
1097 $FrameLenght_z= 2.0*$lengthCRM_active;
1098 
1099 if ($pdsconfig == 0) { #4-pi PDS converage
1100 
1101 #for placing the Arapucas over the cathode
1102  $FrameCenter_x=-1.5*$FrameLenght_x+(4-$nCRM_x/2)/2*$FrameLenght_x;
1103  $FrameCenter_y=$CathodePosX;
1104  $FrameCenter_z=-9.5*$FrameLenght_z+(20-$nCRM_z/2)/2*$FrameLenght_z;
1105 for($i=0;$i<$nCRM_x/2;$i++){
1106 for($j=0;$j<$nCRM_z/2;$j++){
1107  place_OpDetsCathode($FrameCenter_x, $FrameCenter_y, $FrameCenter_z, $i, $j);
1108  $FrameCenter_z+=$FrameLenght_z;
1109 }
1110  $FrameCenter_x+=$FrameLenght_x;
1111  $FrameCenter_z=-9.5*$FrameLenght_z+(20-$nCRM_z/2)/2*$FrameLenght_z;
1112 }
1113 }
1114 
1115 if ($pdsconfig == 0) { #4-pi PDS converage
1116 #for placing the Arapucas on laterals
1117 if ($nCRM_x==8) {
1118  $FrameCenter_y=$posZplane[0]; #anode position
1119  $FrameCenter_z=-19*$FrameLenght_z/2+(40-$nCRM_z)/2*$FrameLenght_z/2;
1120 for($j=0;$j<$nCRM_z/2;$j++){#nCRM will give the collumn number (1 collumn per frame)
1121  place_OpDetsLateral($FrameCenter_y, $FrameCenter_z, $j);
1122  $FrameCenter_z+=$FrameLenght_z;
1123 }
1124 }
1125 
1126 } else { #membrane only PDS converage
1127 
1128 if($pdsconfig == 1){
1129 if ($nCRM_x==8) {
1130  $FrameCenter_y=$posZplane[0]; #anode position
1131  $FrameCenter_z=-19*$FrameLenght_z/2+(40-$nCRM_z)/2*$FrameLenght_z/2;
1132 for($j=0;$j<$nCRM_z/2;$j++){#nCRM will give the collumn number (1 collumn per frame)
1133  place_OpDetsMembOnly($FrameCenter_y, $FrameCenter_z, $j);
1134  $FrameCenter_z+=$FrameLenght_z;
1135 }
1136 }
1137 }
1138 
1139 }
1140 
1141  print CRYO <<EOF;
1142  </volume>
1143 </structure>
1144 </gdml>
1145 EOF
1146 
1147 close(CRYO);
1148 }
1149 
1150 
1151 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1152 #++++++++++++++++++++++++++++++++++++ place_OpDets +++++++++++++++++++++++++++++++++++++++
1153 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1154 
1155 sub place_OpDetsCathode()
1156 {
1157 
1158  $FrameCenter_x = $_[0];
1159  $FrameCenter_y = $_[1];
1160  $FrameCenter_z = $_[2];
1161  $Frame_x = $_[3];
1162  $Frame_z = $_[4];
1163 
1164 #Placing Arapucas over the Cathode
1165 #If there are both top and bottom volumes --> use double-sided:
1166 # <physvol>
1167 # <volumeref ref="volOpDetSensitive_ArapucaDouble_$Frame_x\-$Frame_z\-$ara"/>
1168 # <position name="posOpArapucaDouble$ara-Frame\-$Frame_x\-$Frame_z" unit="cm"
1169 # x="@{[$Ara_X]}"
1170 # y="@{[$Ara_Y]}"
1171 # z="@{[$Ara_Z]}"/>
1172 # </physvol>
1173 #else
1174 for ($ara = 0; $ara<4; $ara++)
1175 {
1176  # All Arapuca centers will have the same Y coordinate
1177  # X and Z coordinates are defined with respect to the center of the current Frame
1178 
1179  $Ara_Y = $FrameCenter_x+$list_posx_bot[$ara]; #GEOMETRY IS ROTATED: X--> Y AND Y--> X
1180  $Ara_X = $FrameCenter_y + $FrameToArapucaSpace;
1181  $Ara_Z = $FrameCenter_z+$list_posz_bot[$ara];
1182 
1183  print CRYO <<EOF;
1184  <physvol>
1185  <volumeref ref="volArapucaDouble_$Frame_x\-$Frame_z\-$ara"/>
1186  <position name="posArapucaDouble$ara-Frame\-$Frame_x\-$Frame_z" unit="cm"
1187  x="@{[$Ara_X]}"
1188  y="@{[$Ara_Y]}"
1189  z="@{[$Ara_Z]}"/>
1190  <rotation name="rPlus90AboutXPlus90AboutZ" unit="deg" x="90" y="0" z="90"/>
1191  </physvol>
1192  <physvol>
1193  <volumeref ref="volOpDetSensitive_ArapucaDouble_$Frame_x\-$Frame_z\-$ara"/>
1194  <position name="posOpArapucaDouble$ara-Frame\-$Frame_x\-$Frame_z" unit="cm"
1195  x="@{[$Ara_X+0.5*$ArapucaOut_y-0.5*$ArapucaAcceptanceWindow_y-0.01]}"
1196  y="@{[$Ara_Y]}"
1197  z="@{[$Ara_Z]}"/>
1198  </physvol>
1199 EOF
1200 
1201 }#end Ara for-loop
1202 
1203 }
1204 
1205 
1206 sub place_OpDetsLateral()
1207 {
1208 
1209  $FrameCenter_y = $_[0];
1210  $FrameCenter_z = $_[1];
1211  $Lat_z = $_[2];
1212 
1213 #Placing Arapucas on the laterals if nCRM_x=8 -- Single Sided
1214 for ($ara = 0; $ara<8; $ara++)
1215 {
1216  # Arapucas on laterals
1217  # All Arapuca centers on a given collumn will have the same Z coordinate
1218  # X coordinates are on the left and right laterals
1219  # Y coordinates are defined with respect to the cathode position
1220  # There are two collumns per frame on each side.
1221 
1222  if ($ara<4) {$Ara_Y = -0.5*$Argon_y + $FrameToArapucaSpaceLat;
1223  $Ara_YSens = ($Ara_Y+0.5*$ArapucaOut_y-0.5*$ArapucaAcceptanceWindow_y-0.01);
1224  $rot= "rIdentity"; }
1225  else {$Ara_Y = 0.5*$Argon_y - $FrameToArapucaSpaceLat;
1226  $Ara_YSens = ($Ara_Y-0.5*$ArapucaOut_y+0.5*$ArapucaAcceptanceWindow_y+0.01);
1227  $rot = "rPlus180AboutX";} #GEOMETRY IS ROTATED: X--> Y AND Y--> X
1228  if ($ara==0||$ara==4) {$Ara_X = $FrameCenter_y-40.0;} #first tile's center 40 cm bellow anode
1229  else {$Ara_X-=$VerticalPDdist;} #other tiles separated by VerticalPDdist
1230  $Ara_Z = $FrameCenter_z;
1231 
1232 # print "lateral arapucas: $Ara_X, $Ara_Y, $Ara_Z \n";
1233 
1234  print CRYO <<EOF;
1235  <physvol>
1236  <volumeref ref="volArapucaLat_$Lat_z\-$ara"/>
1237  <position name="posArapuca$ara-Lat\-$Lat_z" unit="cm"
1238  x="@{[$Ara_X]}"
1239  y="@{[$Ara_Y]}"
1240  z="@{[$Ara_Z]}"/>
1241  <rotationref ref="$rot"/>
1242  </physvol>
1243  <physvol>
1244  <volumeref ref="volOpDetSensitive_ArapucaLat_$Lat_z\-$ara"/>
1245  <position name="posOpArapuca$ara-Lat\-$Lat_z" unit="cm"
1246  x="@{[$Ara_X]}"
1247  y="@{[$Ara_YSens]}"
1248  z="@{[$Ara_Z]}"/>
1249  </physvol>
1250 EOF
1251 
1252 }#end Ara for-loop
1253 
1254 }
1255 
1256 
1257 sub place_OpDetsMembOnly()
1258 {
1259 
1260  $FrameCenter_y = $_[0];
1261  $FrameCenter_z = $_[1];
1262  $Lat_z = $_[2];
1263 
1264 #Placing Arapucas on the laterals if nCRM_x=8 -- Single Sided
1265 for ($ara = 0; $ara<18; $ara++)
1266 {
1267  # Arapucas on laterals
1268  # All Arapuca centers on a given collumn will have the same Z coordinate
1269  # X coordinates are on the left and right laterals
1270  # Y coordinates are defined with respect to the cathode position
1271  # There are two collumns per frame on each side.
1272 
1273  if($ara<9) {$Ara_Y = -0.5*$Argon_y + $FrameToArapucaSpaceLat;
1274  $Ara_YSens = ($Ara_Y+0.5*$ArapucaOut_y-0.5*$ArapucaAcceptanceWindow_y-0.01);
1275  $rot= "rIdentity"; }
1276  else {$Ara_Y = 0.5*$Argon_y - $FrameToArapucaSpaceLat;
1277  $Ara_YSens = ($Ara_Y-0.5*$ArapucaOut_y+0.5*$ArapucaAcceptanceWindow_y+0.01);
1278  $rot = "rPlus180AboutX";} #GEOMETRY IS ROTATED: X--> Y AND Y--> X
1279  if($ara==0||$ara==9) {$Ara_X = $FrameCenter_y-$ArapucaOut_x/2;} #first tile's center right below anode
1280  else {$Ara_X-=$ArapucaOut_x - $FrameToArapucaSpace;} #other tiles separated by minimal distance + buffer
1281  $Ara_Z = $FrameCenter_z;
1282 
1283 # print "lateral arapucas: $Ara_X, $Ara_Y, $Ara_Z \n";
1284 
1285  print CRYO <<EOF;
1286  <physvol>
1287  <volumeref ref="volArapucaLat_$Lat_z\-$ara"/>
1288  <position name="posArapuca$ara-Lat\-$Lat_z" unit="cm"
1289  x="@{[$Ara_X]}"
1290  y="@{[$Ara_Y]}"
1291  z="@{[$Ara_Z]}"/>
1292  <rotationref ref="$rot"/>
1293  </physvol>
1294  <physvol>
1295  <volumeref ref="volOpDetSensitive_ArapucaLat_$Lat_z\-$ara"/>
1296  <position name="posOpArapuca$ara-Lat\-$Lat_z" unit="cm"
1297  x="@{[$Ara_X]}"
1298  y="@{[$Ara_YSens]}"
1299  z="@{[$Ara_Z]}"/>
1300  </physvol>
1301 EOF
1302 
1303 }#end Ara for-loop
1304 
1305 
1306 
1307 }
1308 
1309 
1310 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1311 #+++++++++++++++++++++++++++++++++++++ gen_Enclosure +++++++++++++++++++++++++++++++++++++
1312 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1313 
1314 sub gen_Enclosure()
1315 {
1316 
1317 # Create the detector enclosure fragment file name,
1318 # add file to list of output GDML fragments,
1319 # and open it
1320  $ENCL = $basename."_DetEnclosure" . $suffix . ".gdml";
1321  push (@gdmlFiles, $ENCL);
1322  $ENCL = ">" . $ENCL;
1323  open(ENCL) or die("Could not open file $ENCL for writing");
1324 
1325 
1326 # The standard XML prefix and starting the gdml
1327  print ENCL <<EOF;
1328 <?xml version='1.0'?>
1329 <gdml>
1330 EOF
1331 
1332 
1333 # All the detector enclosure solids.
1334 print ENCL <<EOF;
1335 <solids>
1336 
1337  <box name="FoamPadBlock" lunit="cm"
1338  x="@{[$Cryostat_x + 2*$FoamPadding]}"
1339  y="@{[$Cryostat_y + 2*$FoamPadding]}"
1340  z="@{[$Cryostat_z + 2*$FoamPadding]}" />
1341 
1342  <subtraction name="FoamPadding">
1343  <first ref="FoamPadBlock"/>
1344  <second ref="Cryostat"/>
1345  <positionref ref="posCenter"/>
1346  </subtraction>
1347 
1348  <box name="SteelSupportBlock" lunit="cm"
1349  x="@{[$Cryostat_x + 2*$FoamPadding + 2*$SteelSupport_x]}"
1350  y="@{[$Cryostat_y + 2*$FoamPadding + 2*$SteelSupport_y]}"
1351  z="@{[$Cryostat_z + 2*$FoamPadding + 2*$SteelSupport_z]}" />
1352 
1353  <subtraction name="SteelSupport">
1354  <first ref="SteelSupportBlock"/>
1355  <second ref="FoamPadBlock"/>
1356  <positionref ref="posCenter"/>
1357  </subtraction>
1358 
1359  <box name="DetEnclosure" lunit="cm"
1360  x="$DetEncX"
1361  y="$DetEncY"
1362  z="$DetEncZ"/>
1363 
1364 </solids>
1365 EOF
1366 
1367 
1368 # Detector enclosure structure
1369  print ENCL <<EOF;
1370 <structure>
1371  <volume name="volFoamPadding">
1372  <materialref ref="fibrous_glass"/>
1373  <solidref ref="FoamPadding"/>
1374  </volume>
1375 
1376  <volume name="volSteelSupport">
1377  <materialref ref="AirSteelMixture"/>
1378  <solidref ref="SteelSupport"/>
1379  </volume>
1380 
1381  <volume name="volDetEnclosure">
1382  <materialref ref="Air"/>
1383  <solidref ref="DetEnclosure"/>
1384 
1385  <physvol>
1386  <volumeref ref="volFoamPadding"/>
1387  <positionref ref="posCryoInDetEnc"/>
1388  </physvol>
1389  <physvol>
1390  <volumeref ref="volSteelSupport"/>
1391  <positionref ref="posCryoInDetEnc"/>
1392  </physvol>
1393  <physvol>
1394  <volumeref ref="volCryostat"/>
1395  <positionref ref="posCryoInDetEnc"/>
1396  </physvol>
1397 EOF
1398 
1399 
1400 print ENCL <<EOF;
1401  </volume>
1402 EOF
1403 
1404 print ENCL <<EOF;
1405 </structure>
1406 </gdml>
1407 EOF
1408 
1409 close(ENCL);
1410 }
1411 
1412 
1413 
1414 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1415 #+++++++++++++++++++++++++++++++++++++++ gen_World +++++++++++++++++++++++++++++++++++++++
1416 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1417 
1418 sub gen_World()
1419 {
1420 
1421 # Create the WORLD fragment file name,
1422 # add file to list of output GDML fragments,
1423 # and open it
1424  $WORLD = $basename."_World" . $suffix . ".gdml";
1425  push (@gdmlFiles, $WORLD);
1426  $WORLD = ">" . $WORLD;
1427  open(WORLD) or die("Could not open file $WORLD for writing");
1428 
1429 
1430 # The standard XML prefix and starting the gdml
1431  print WORLD <<EOF;
1432 <?xml version='1.0'?>
1433 <gdml>
1434 EOF
1435 
1436 
1437 # All the World solids.
1438 print WORLD <<EOF;
1439 <solids>
1440  <box name="World" lunit="cm"
1441  x="@{[$DetEncX+2*$RockThickness]}"
1442  y="@{[$DetEncY+2*$RockThickness]}"
1443  z="@{[$DetEncZ+2*$RockThickness]}"/>
1444 </solids>
1445 EOF
1446 
1447 # World structure
1448 print WORLD <<EOF;
1449 <structure>
1450  <volume name="volWorld" >
1451  <materialref ref="DUSEL_Rock"/>
1452  <solidref ref="World"/>
1453 
1454  <physvol>
1455  <volumeref ref="volDetEnclosure"/>
1456  <position name="posDetEnclosure" unit="cm" x="$OriginXSet" y="$OriginYSet" z="$OriginZSet"/>
1457  </physvol>
1458 
1459  </volume>
1460 </structure>
1461 </gdml>
1462 EOF
1463 
1464 # make_gdml.pl will take care of <setup/>
1465 
1466 close(WORLD);
1467 }
1468 
1469 
1470 
1471 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1472 #++++++++++++++++++++++++++++++++++++ write_fragments ++++++++++++++++++++++++++++++++++++
1473 #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1474 
1475 sub write_fragments()
1476 {
1477  # This subroutine creates an XML file that summarizes the the subfiles output
1478  # by the other sub routines - it is the input file for make_gdml.pl which will
1479  # give the final desired GDML file. Specify its name with the output option.
1480  # (you can change the name when running make_gdml)
1481 
1482  # This code is taken straigh from the similar MicroBooNE generate script, Thank you.
1483 
1484  if ( ! defined $output )
1485  {
1486  $output = "-"; # write to STDOUT
1487  }
1488 
1489  # Set up the output file.
1490  $OUTPUT = ">" . $output;
1491  open(OUTPUT) or die("Could not open file $OUTPUT");
1492 
1493  print OUTPUT <<EOF;
1494 <?xml version='1.0'?>
1495 
1496 <!-- Input to Geometry/gdml/make_gdml.pl; define the GDML fragments
1497  that will be zipped together to create a detector description.
1498  -->
1499 
1500 <config>
1501 
1502  <constantfiles>
1503 
1504  <!-- These files contain GDML <constant></constant>
1505  blocks. They are read in separately, so they can be
1506  interpreted into the remaining GDML. See make_gdml.pl for
1507  more information.
1508  -->
1509 
1510 EOF
1511 
1512  foreach $filename (@defFiles)
1513  {
1514  print OUTPUT <<EOF;
1515  <filename> $filename </filename>
1516 EOF
1517  }
1518 
1519  print OUTPUT <<EOF;
1520 
1521  </constantfiles>
1522 
1523  <gdmlfiles>
1524 
1525  <!-- The GDML file fragments to be zipped together. -->
1526 
1527 EOF
1528 
1529  foreach $filename (@gdmlFiles)
1530  {
1531  print OUTPUT <<EOF;
1532  <filename> $filename </filename>
1533 EOF
1534  }
1535 
1536  print OUTPUT <<EOF;
1537 
1538  </gdmlfiles>
1539 
1540 </config>
1541 EOF
1542 
1543  close(OUTPUT);
1544 }
1545 
1546 
1547 print "Some of the principal parameters for this TPC geometry (unit cm unless noted otherwise)\n";
1548 print " CRM active area : $widthCRM_active x $lengthCRM_active\n";
1549 print " CRM total area : $widthCRM x $lengthCRM\n";
1550 print " Wire pitch in Y, Z : $wirePitchY, $wirePitchZ\n";
1551 print " TPC active volume : $driftTPCActive x $widthTPCActive x $lengthTPCActive\n";
1552 print " Argon volume : ($Argon_x, $Argon_y, $Argon_z) \n";
1553 print " Argon buffer : ($xLArBuffer, $yLArBuffer, $zLArBuffer) \n";
1554 print " Detector enclosure : $DetEncX x $DetEncY x $DetEncZ\n";
1555 print " TPC Origin : ($OriginXSet, $OriginYSet, $OriginZSet) \n";
1556 print " Field Cage : $FieldCage_switch \n";
1557 print " Cathode : $Cathode_switch \n";
1558 print " Workspace : $workspace \n";
1559 print " Wires : $wires_on \n";
1560 
1561 # run the sub routines that generate the fragments
1562 if ( $FieldCage_switch eq "on" ) { gen_FieldCage(); }
1563 #if ( $Cathode_switch eq "on" ) { gen_Cathode(); } #Cathode for now has the same geometry as the Ground Grid
1564 
1565 gen_Extend(); # generates the GDML color extension for the refactored geometry
1566 gen_Define(); # generates definitions at beginning of GDML
1567 gen_Materials(); # generates materials to be used
1568 gen_TPC(); # generate TPC for a given unit CRM
1569 gen_Cryostat(); #
1570 gen_Enclosure(); #
1571 gen_World(); # places the enclosure among DUSEL Rock
1572 write_fragments(); # writes the XML input for make_gdml.pl
1573  # which zips together the final GDML
1574 print "--- done\n\n\n";
1575 exit;