-
Notifications
You must be signed in to change notification settings - Fork 0
/
M20.g4bl
687 lines (533 loc) · 33.6 KB
/
M20.g4bl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
########################################################################################################################################################
#
# This is a simulation of the M20 beamline
# written by Tunok Mondol, Winter 2021
# Additional files needed to run this g4bl file:
# Kickerfield.txt, 12dq12map.txt, 12q12map.txt and
# BOpera9ASep.txt
#
########################################################################################################################################################
# Simulation to date April, 2021.
# Things that need to be resolved:
# a) The Separator fields (the physical dimensions of plate gap and magnet gap are that of M20) are from that of M9A Separators.
# b) The Septum Magnet was simply implemented with 2 "idealsectorbend" one beside the other and length was optimized with GMINUIT for reducing "misalignment". A # correct fieldmap should be implemented.
#
# Simulation is tuned only for muon focusing and max luminosity. Non-spin rotation, particle decays, Pion production, Cloud muons. etc not considered.
# 12q12 Quadrupole and 12dq12 QuadSteering fields introduced, however the conversion of 12dq12 & 12q12 into current for steering is required.
# Intermediate beam pipes: Q2pipe,ChicagoPipe,Bend1_Sep1pipe,Sep2_Bend2pipe,Q8_9pipe,Q10_12pipe are designed and placed with guesses.
# Genericbends are placed in between two 'corner' commands with each having the half-angle of the bending angle, as per suggestions by Fred Jones.
#
########################################################################################################################################################
# M20 Geometry along the CENTERLINE (Overall Assembly Drawing: TB21633)
#
# Center (mm) Center (inch) (TB21633) Center (mm/inch) (Estimated)
# T2 0 0
# Q1 1310.0304 51.576
# Q2 1963.023 77.245
# Scraper 3200/125.9842
# B1 4622.038 181.97
# S1 7258.05 285.75
# Q3 8761.8832 344.958
# Q4 9333.8658 367.477
# Q5 9866.8848 388.462
# Q6 10397.8464 409.368
# Q7 10,970.8958 431.927
# XSLIT0 11300.9894/444.9208
# S2 12430.9894 489.411
# YSLIT1 13560.9894/533.8972
# B2 15069.058 593.27
# Q8 15938.9064 627.516
# Q9 16499.9162 649.603
# K1 18330.0624 721.656
# Driftbox 20487.6908/806.602
# Septum 21429.1418 843.667
# XSLIT1D/C 22429.1418/883.037
# Q10-D/C 23953.9018 943.067
# Q11-D/C 24515.9022 965.193
# Q12-D/C 25077.9026 987.319
#
########################################################################################################################################################
#
# Drawings used:
# Overall M20 Assembly - TB21633
# Q1 - 8QN16M/7
# Q1,Q2 Beampipe - TB 20867D
# Q2,3,7,8,9,10,11,12 - QPMA 019003(Revision -2)
# Q4,5,6 - Chicago 7Q7
# Scraper - TB21790, TB21835
# B1, B2 - Scanditronix Technical Documentation - Dipole Magnet
# B1, B2 Entrance, Exit - TB21759
# Separator/Wien Filter - 7103002458(E-plate and vacBox dimensions), RFProp-M20-SEP-v2.odt (Request for Proposal)
# Kicker - 7103002410, 7103002412, 7103002416
# DriftBox - TB21810
# Septum - Scanditronix Technical Documentation - Septum, SEMA 001010, SEMA 001021, SEMA 001022
# Slits, Gates, Bellows - VAT 12148-PA44-0001, VAT 16248-PA41-0001, TB21579, TB21607, TB21634
#
########################################################################################################################################################
physics QGSP_BIC spinTracking=1
trackcuts keep=mu+,e+ #,proton,pi+,e+,anti_nu_mu,nu_e
param histoFile=M20.root
param -unset PROFILE_D=./m20bl_profileD.txt PROFILE_C=./m20bl_profileC.txt
# parameter 's' is the Cumulative S-Layout as described in "Notes on Using G4Beamline" by Fred Jones. zFocus is the position of the final focus. targetNo can be chosen between # 0 and 1. Refer to "Target" section.
param -unset inch=25.4 s=0 zFocus=26278 nEvents=300000 targetNo=0
# Beamline Parameters
param -unset BL=1 # BL can be 1,2 or any other value where '1' represents M20 with only D leg, '2' does M20 with only C leg and otherwise it will produce the entire beamline.
# Values of tuned current/voltage in 'Amp'/'Volt' for the quads
if $BL==1 # Tunes for only D leg
param -unset Chisq=26391.2 P=28.64796 dP=3 dX=10 dY=5 dXp=0.05 dYp=0.05 Q1=-33.954834 Q2=72.803772 B1=91.715524 Q3=84.832049 Q4=-79.526047 Q5=69.238529 Q6=-89.643757 Q7=96.532851 B2=92.310794 Q8T=0.20620102 Q8B=41.603141 Q9=-21.425507 Q10DB=-0.01663766 Q10DT=0.68211732 Q10CB=-0.0087409701 Q10CT=2.7418251 Q11DL=-0.61942555 Q11DR=-75.573346 Q11CL=0.051537266 Q11CR=-77.658069 Q12D=111.17829 Q12C=104.79 SM1D=591.77515 SM1C=537.82754 K1=9317.2666 XSLIT0=170 YSLIT1=170 XSLIT1D=170 XSLIT1C=170
# NSR
#Chisq=20171.9 P=28.830843 dP=3 dX=10 dY=5 dXp=0.05 dYp=0.05 Q1=-35.092869 Q2=74.153194 B1=91.659896 Q3=87.03341 Q4=-78.701606 Q5=69.864461 Q6=-87.828481 Q7=98.737528 B2=92.417056 Q8T=0.029556707 Q8B=42.720258 Q9=-23.36234 Q10DB=-0.012752652 Q10DT=1.2561858 Q10CB=-0.0087409701 Q10CT=2.7418251 Q11DL=-0.61812826 Q11DR=-75.548041 Q11CL=0.051537266 Q11CR=-77.658069 Q12D=110.99466 Q12C=104.79 SM1D=595.77564 SM1C=537.82754 K1=9493.5595 XSLIT0=170 YSLIT1=170 XSLIT1D=170 XSLIT1C=170
elseif $BL==2 # Tunes for only C leg
param -unset Chisq=32774.9 P=28.088473 dP=3 dX=10 dY=5 dXp=0.05 dYp=0.05 Q1=-34.27704 Q2=73.4977 B1=89.964951 Q3=83.004018 Q4=-76.528024 Q5=70.167314 Q6=-90.179082 Q7=96.133368 B2=91.498756 Q8T=0.19307551 Q8B=41.243296 Q9=-21.643377 Q10DB=-0.01663766 Q10DT=0.68211732 Q10CB=-0.0072477314 Q10CT=4.1414194 Q11DL=-0.61942555 Q11DR=-75.573346 Q11CL=0.075910099 Q11CR=-76.435973 Q12D=111.17829 Q12C=103.61998 SM1D=591.77515 SM1C=564.769 K1=-8329.602 XSLIT0=170 YSLIT1=170 XSLIT1D=170 XSLIT1C=170
else
param -unset Chisq=26391.2 P=28.64796 dP=3 dX=10 dY=5 dXp=0.05 dYp=0.05 Q1=-33.954834 Q2=72.803772 B1=91.715524 Q3=84.832049 Q4=-79.526047 Q5=69.238529 Q6=-89.643757 Q7=96.532851 B2=92.310794 Q8T=0.20620102 Q8B=41.603141 Q9=-21.425507 Q10DB=-0.01663766 Q10DT=0.68211732 Q10CB=-0.0087409701 Q10CT=2.7418251 Q11DL=-0.61942555 Q11DR=-75.573346 Q11CL=0.051537266 Q11CR=-77.658069 Q12D=111.17829 Q12C=104.79 SM1D=591.77515 SM1C=537.82754 K1=9317.2666 XSLIT0=170 YSLIT1=170 XSLIT1D=170 XSLIT1C=170
#Chisq=34300699226 P=28.59065 dP=2 dX=10 dY=5 dXp=0.025 dYp=0.025 Q1=-33.954834 Q2=72.803772 B1=91.715524 Q3=81.159428 Q4=-54.086886 Q5=38.910438 Q6=-81.120617 Q7=93.801147 B2=92.842371 Q8T=0.061972332 Q8B=49.142003 Q9=-2.6113863 Q10DB=0.0021351042 Q10DT=-0.80486429 Q10CB=-0.011088366 Q10CT=3.8929161 Q11DL=-0.062298353 Q11DR=-72.977205 Q11CL=0.07716875 Q11CR=-79.319223 Q12D=103.84862 Q12C=106.57949 SM1D=573.41212 SM1C=586.60508 K1=677.15475 XSLIT0=170 YSLIT1=170 XSLIT1D=170 XSLIT1C=170
endif
# Elements' maximum current capacity and maximum field strength
param -unset RadHardMaxCurr=750 RadHardMaxStng=7.08 DipoleMaxCurr=120 DipoleMaxStng=0.125 12Q12MaxCurr=200 12Q12MaxStng=0.98 ChicagoMaxCurr=1000 ChicagoMaxStng=6.3 SeptumMaxCurr=775 SeptumMaxStng=0.114
trackcolor mu+=1,1,0 e+=1,0,0 proton=0.2,0.2,0.2 pi+=0,0,1 reference=1,1,1
##########################################################
# #
# TARGET #
# #
##########################################################
if $targetNo==0
box T2 width=100.0 height=100.0 length=10.0 material=Vacuum color=1,0,1,0.3
particlefilter Filter radius=500 require='(sqrt(Px^2 + Py^2 + Pz^2)) < 29.85' # This puts a physical object that looks just like a virtualdetector that instead acts \
as filter. This is used to replicate the initial momentum distribution of surface muons after production (0 MeV/c - Surface Muons Edge) with peak at 29.8 MeV/c.
### This part was copied from the g4bl script for T2 target with M20 and M9A first quads, acquired from one of Syd's folders.
### It could be called by "targetNo=1" in the command line, but it wasn't used to develop this simulation or for tuning.
elseif $targetNo==1 # This was not used to test or minimize anything for the simulation itself, just incorporated in case one is interested. Can be put in the simulation using "targetNo=1" from the command line, otherwise the simpler version of beamsource will be chosen to run the simulation.
param targetMaterial='Be'
param -unset xPBeamT2Surface=2.54 targetLen=50
param xTarget=-.1*25.4+$xPBeamT2Surface
#the Be target is .2" thick in the x-direction
param TargetNoDescription='H2OBe'
* the Standard SS encased Water Jacketed Be T2 target ... the Be is $targetLen long
#see drawing TBP0301D
#water jacket ss enclosure endRadius(outside)=.234" midSecLength=.762-.468 wallthickness=.010" length=$targetLen+.25"
#Be target has a rectangular cross section .58" vertical, .2" horizontal
tubs EndSec_SS \
outerRadius=.234*25.4 initialPhi=0 finalPhi=180 length=($targetLen+.145*25.4) material=Stainless304 color=.7,.7,.7,.5
trap MidSec_SS \
height=(.762-.468)*25.4 upperWidth=.468*25.4 lowerWidth=.468*25.4 length=($targetLen+.145*25.4) material=Stainless304 color=.7,.7,.7,.5
tubs EndSec_Water \
outerRadius=.224*25.4 initialPhi=0 finalPhi=180 length=($targetLen+.125*25.4) material=WATER color=0,0,1,.5
trap MidSec_Water \
height=(.742-.448)*25.4 upperWidth=.448*25.4 lowerWidth=.448*25.4 length=($targetLen+.125*25.4) material=WATER color=0,0,1,.5
trap Be_target \
height=.58*25.4 upperWidth=.2*25.4 lowerWidth=.2*25.4 length=$targetLen material=$targetMaterial color=0,1,1,.5
place Be_target height=.5*(.58+.448-.742)*25.4 y=.25*(.58+.448-.742)*25.4 parent=EndSec_Water rename=+_Be_target
place EndSec_Water parent=EndSec_SS rename=+_Water
place Be_target height=(.742-.448)*25.4 parent=MidSec_Water rename=+_Be_target
place MidSec_Water parent=MidSec_SS rename=+_Water
group T2target material=Vacuum length=($targetLen+.145*25.4)
place EndSec_SS rename=+_EndSec_SS_Up z=0 y=(.762-.468)*25.4*.5
place MidSec_SS z=0 rename=+_MidSec_SS
place EndSec_SS rename=+_EndSec_SS_Dn z=0 y=-(.762-.468)*25.4*.5 rotation=Z180
endgroup
particlefilter Filter_ends radius=70 length=1 color=1,1,1,0.4 keep=mu+ require='(sqrt(Px^2 + Py^2 + Pz^2)) < 29.85'
particlefilter Filter_wall innerRadius=70 radius=71 length=400 keep=mu+ color=1,1,1,0.4 require='(sqrt(Px^2 + Py^2 + Pz^2)) < 29.85'
box Cyc width=10.0 height=10 length=1 material=Vacuum color=1,0.7,0.8
group T2Assembly length=900
place T2target z=0
place Filter_ends z=-300
place Filter_ends z=100
place Filter_wall z=-100
place Cyc z=-200
endgroup
endif
##########################################################
# #
# BEAM #
# #
##########################################################
if $targetNo==0
beam gaussian z=0 sigmaX=$dX sigmaY=$dY P=$P \
particle=mu+ sigmaXp=$dXp sigmaYp=$dYp sigmaP=$dP \
polarization=0,0,-1 nEvents=$nEvents/2
beam gaussian z=0 sigmaX=$dX sigmaY=$dY P=$P \
particle=e+ sigmaXp=$dXp sigmaYp=$dYp sigmaP=$dP \
nEvents=$nEvents/2
elseif $targetNo==1
param -unset pKE=480
#KE=sqrt(p^2+m^2)-m, m_p=938.27; m_n=939.565, m_e=.511. m_pi=139.57
param Pp=sqrt($pKE^2+2*938.2*$pKE)
#proton rest mass=938.2MeV, TRIUMF bl1A proton KE=490MeV,P=sqrt(KE^2+2*KE*m) & KE=sqrt(P^2+m^2)-m [energy in MeV & momentum in MeV/c]
beam gaussian z=0 x=200 rotation=Y-90 P=$Pp sigmaP=10 beamHeight=0.5 beamWidth=0.5 nEvents=625000000 particle=proton
endif
reference particle=mu+ referenceMomentum=$P noEloss=1 noEfield=0 beamX=0 beamY=0 beamZ=0
# the 'tune' command can be used to tune the Dipole bending magnets' strengths given that we know the momentum of the beam. It has to be used in conjunction with the 'reference' particle command.
#tune B1y z0=0 z1=5000 initial=0.1 step=0.01 expr=x1 tolerance=0.000001
# definition of the materials to use in the elements.
material PiMuVacuum Vacuum,1.0 density=0 keep=pi+,mu+,e+,proton,tau+
material MuVacuum Vacuum,1.0 density=0 keep=mu+
material AlKill Al,1.0 density=2.70 kill=mu+,pi+,e+,proton,tau+
material FeKill Fe,1.0 density=7.874 kill=mu+,pi+,e+,proton,tau+
material WKill W,1.0 density=19.25 kill=mu+,pi+,e+,proton,tau+
##########################################################
# #
# Detectors, Scraper, Pipes #
# #
##########################################################
virtualdetector Det radius=500 length=1 material=Vacuum format=rootExtended
beamlossntuple Lost
tubs FocusOuter outerRadius=25 innerRadius=0 length=0.1 material=Vacuum color=1,0,0,0.5
tubs Focus3rd outerRadius=20 innerRadius=0 length=0.1 material=Vacuum color=1,0.5,0,0.5
tubs Focus2nd outerRadius=15 innerRadius=0 length=0.1 material=Vacuum color=1,1,0,0.5
tubs Focus1st outerRadius=10 innerRadius=0 length=0.1 material=Vacuum color=0,0,1,0.5
tubs FocusCenter outerRadius=5 innerRadius=0 length=0.1 material=Vacuum color=0,.8,0,0.5
group FocalArea
place FocusCenter parent=Focus1st kill=1
place Focus1st parent=Focus2nd kill=1
place Focus2nd parent=Focus3rd kill=1
place Focus3rd parent=FocusOuter kill=1
place FocusOuter kill=1
endgroup
# SCRAPER construction. The vertices are put in with a rough guess and adjusted from TB21790 and TB21835. The 12" entrance tube (ScraperEnt) length and the 8" exit tube (ScraperExit) length are estimated to extend to the next element right after and before.
extrusion ScraperIron vertices=1385,692.15;1450,600;500,-692.15;-720,-692.15;-1058,-50;-651,692.15 length=54.5*25.4 material=Fe color=1,1,1,.5
tubs ScraperEnt innerRadius=(12/2)*$inch outerRadius=(12/2)*$inch+50 length=18*25.4 material=Vacuum color=0.6,0.2,0.9,0.5 # actual length=7.28*25.4
tubs ScraperUps outerRadius=(10/2)*$inch-3 length=31.143*$inch material=Vacuum color=0.6,0.6,0.6
tubs ScraperDns outerRadius=(8/2)*25.4-3 length=(54.5-31.143)*$inch material=Vacuum color=0.6,0.6,0.6
tubs ScraperExit outerRadius=8*25.4/2-3 length=2*$inch material=Vacuum color=0.6,0.6,0.6 # the pipe in between B1 and Scraper.
place ScraperUps z=0 y=((-54.5/2)+(31.143/2))*$inch rotation=X-90 parent=ScraperIron rename=+Ups
place ScraperDns z=0 y=((54.5/2)-(23.357/2))*$inch rotation=X-90 parent=ScraperIron rename=+Dns
group Scraper #length=(54.5+3+18)*$inch
place ScraperEnt z=(18/2)*$inch kill=1
place ScraperIron z=(18+54.5/2)*$inch rotation=X90
place ScraperExit z=(18+54.5+2/2)*$inch
endgroup
# These cylinders are defined and placed to fill up intermediate spaces to ensure no beam spills, which would otherwise hinder optimization.
tubs Q2pipe outerRadius=(12/2)*$inch innerRadius=(11.712/2)*$inch length=20.34*$inch material=WKill color=1,1,1,0.4
tubs ChicagoPipe outerRadius=100 innerRadius=95 length=1700 material=Vacuum color=1,1,1,0.4
tubs Bend1_Sep1pipe outerRadius=5*$inch innerRadius=4*$inch length=34*$inch material=WKill color=1,1,1,0.4
tubs Sep2_Bend2pipe outerRadius=155 innerRadius=150 length=25*$inch material=WKill color=1,1,1,0.4
tubs Q8_9pipe outerRadius=155 innerRadius=150 length=60*$inch material=WKill color=1,1,1,0.4
tubs Q10_12pipe outerRadius=155 innerRadius=150 length=120*$inch material=WKill color=1,1,1,0.4 #length=2666.9492
##########################################################
# #
# QUADS and DIPOLES #
# #
##########################################################
# The Rad-Hard Quadrupole magnet box, to be used as a parent of genericquad 'Q1'.
box Quad1box width=635 height=889 length=698.5 material=FeKill color=0,.6,0,0.5
genericquad Quad1 fieldLength=490.5 apertureRadius=(8.25*25.4)/2 ironRadius=635/2 ironLength=698.5 ironColor=0,.6,0,0.5 kill=1
place Quad1 gradient=($Q1/750)*7.08 z=0 parent=Quad1box
tubs Q1pipe outerRadius=(8.25/2)*$inch length=698.5 material=Vacuum color=1,1,1,0.3
#actual pipe dimension outerRadius=(8.25/2)*$inch innerRadius=(7.76/2)*$inch length=21.74*$inch
place Q1pipe z=0 parent=Quad1box
genericquad Chicago fieldLength=322 apertureRadius=(7.935*$inch)/2 ironRadius=479.4 ironLength=7*$inch ironColor=0,.6,0,0.5 kill=1
genericquad 12Q12 fieldLength=322 apertureRadius=(12.250*$inch)/2 ironRadius=479.4 ironLength=12*$inch ironColor=0,.6,0,0.5 kill=1
# These definitions are for same Chicago and 12Q12 Quadrupoles but with 'rounded +' aperture.
#genericquad Chicago fieldLength=322 poleTipRadius=(7.935*25.4)/2 coilRadius=180 \
# coilHalfWidth=50 ironRadius=479.4 ironLength=501.65 ironColor=0,.6,0,0.5 kill=1
# coilRadius,coilHalfWidth ARE NOT CORRECT. irontRadius is the distance from the center of the magnet to its diagonal tip.
#genericquad 12Q12 fieldLength=322 poleTipRadius=(12.250*25.4)/2 coilRadius=245 \
# coilHalfWidth=50 ironRadius=479.4 ironLength=501.65 ironColor=0,.6,0,0.5 kill=1
# coilRadius,coilHalfWidth are not correct. irontRadius is the distance from the center of the magnet to its diagonal tip.
tubs BendGate innerRadius=4*$inch outerRadius=12.3*$inch length=10 material=FeKill
tubs BendGate2 innerRadius=4*$inch outerRadius=20.3*$inch length=10 material=FeKill # outerRadius=20" was chosen so that it is big enough to kill any possible outgoing muon.
genericbend Bend fieldWidth=790 fieldHeight=229 fieldLength=600 \
ironColor=1,0,0,0.5 ironWidth=930 ironHeight=485 ironLength=660.68 ironMaterial=FeKill
# Steering Quad Fieldmaps
fieldmap 12q12map file=12q12map.txt
fieldmap 12dq12map file=12dq12map.txt
##########################################################
# #
# WIEN FILTER #
# #
##########################################################
# Followed Syd's footsteps in making this not so complicated physical obstruction with bits and pieces of information from 7103002458(E-plate and vacBox dimensions) and
# RFProp-M20-SEP-v2.odt (Request for Proposal).
#Wien Filters have 10" windows
#E-Feild plates have a 80(113)mm gap, ~ 1500(1463)mm long (end-end), ~240(271.5)mm tall (recall rounded ends), ~44(64)mm thick (not counting the ends)
#B-Feild poles have a 560(530)mm gap, 121(105)mm thick, 560(570)mm width, 1500(1300) mm length
#Vacuum Chamber has height 560(530), width 510(636), length 1960(1790)
box VacBox height=560 width=510 length=1960 material=Vacuum color=0.9,0.2,1,.5
tubs Face outerRadius=510/2 innerRadius=254/2 length=10 material=FeKill color=0.9,0.2,0.5,.5 kill=1
box BPole height=121 width=560 length=1500 material=FeKill color=1,0,0,0.5 kill=1
box EPlate height=240 width=44 length=1500 material=FeKill color=0.2,.2,1,.5 kill=1
place EPlate rotation=Z90 z=0 y=-.5*(80+44) parent=VacBox rename=+Neg
place EPlate rotation=Z90 z=0 y=+.5*(80+44) parent=VacBox rename=+Pos
place Face z=-.5*(1960-10) parent=VacBox
place Face z=+.5*(1960-10) parent=VacBox
group WienRot length=1960 height=805 width=900
place BPole z=0 x=.5*(510+121) rotation=Z90 rename=+Neg
place VacBox z=0
place BPole z=0 x=-.5*(510+121) rotation=Z-90 rename=+Pos
endgroup
fieldmap Separatorfield file=BEoperaM9ASep.txt
#the amount of spin precession is governed by the rest frame values of the magnetic field and the time spent in that field.
#for a boost vector b_vec E'||=E|| B'||=B||, E'_|=gamma*(E_vec+b_vec*X*B_vec)_| B'_|= gamma*(B_vec-b_vec*X*E_vec/c^2)_| *X*=crossproduct; || or _| = parallel or perp component to b_vec
param -unset WFBy=.004 WFVxppt=-1.66
param m_mu=105.65837 C=299792458 #MeV/c^2, m/s
param Velmu=sqrt($P^2/($m_mu^2+$P^2)) #Wien Filter E/B ratio is Velmu*c
param WF_EB_ratio=$Velmu*$C/1E6 #80.619 for 29.5MeV/c muons, this is about right, for B=500G=.05T need E~400kV/10cm=4MeV/m (the M15 values)
param WFVx=$WF_EB_ratio*(1+.001*$WFVxppt)*$WFBy
##########################################################
# #
# KICKER #
# #
##########################################################
# Electrodes have width 50mm, height 296.56mm, and length 1913mm. Dimensions taken from Kicker Drawings.
# VacuumBox have length 2223mm, width 500mm, height 500mm
box KickerVac height=500 width=500 length=2223 material=Vacuum color=0.1,0.5,0,0.5
box Electrode height=296.56 width=50 length=1913 material=FeKill color=0.2,.2,1,.5 kill=1
tubs KickerGate outerRadius=540/2 innerRadius=8*25.4/2 length=10 material=FeKill color=0.2,0.2,1,0.5
place Electrode z=0 x=-(75+25) rename=+Left parent=KickerVac
place Electrode z=0 x=(75+25) rename=+Right parent=KickerVac
fieldmap Kickerfield file=Kickerfield.txt
group Kicker length=2273 width=540 height=540
place KickerGate z=-0.5*(2273-10)
place Kickerfield z=0 current=1 gradient=$K1/10000 # Half length of the kicker is 2273/2=1136.5 #
place KickerVac z=0
place KickerGate z=0.5*(2273-10)
endgroup
box Driftwall height=10*$inch width=16.09*$inch length=(58.13+11)*$inch material=FeKill color=1,0,0,0.5
box Driftvac height=9.58*$inch width=15.09*$inch length=(58.13+11)*$inch material=Vacuum color=1,1,1,0.3
place Driftvac z=0 parent=Driftwall
tubs Driftgate outerRadius=(12/2)*$inch+100 innerRadius=(8/2)*$inch length=2*$inch color=0,0,1,0.5 material=FeKill
group Driftbox length=(60.13+11)*$inch width=18.11*$inch+200 height=12.60*$inch+340
place Driftwall z=$inch
place Driftgate z=((-60.13-11+2)/2)*$inch
endgroup
##########################################################
# #
# SEPTUM #
# #
##########################################################
param -unset h=100 # MINUIT calculated the radius of curvature needs to be increased by ~100mm for a "nicer" beam #
box SMEnt height=320 width=460 length=25 material=FeKill color=1,0,1
box Gate height=226 width=168.5 length=25 material=Vacuum color=1,1,1,0.8
if $BL==1
place Gate z=0 x=-99.5 rename=+right parent=SMEnt
elseif $BL==2
place Gate z=0 x=99.5 rename=+left parent=SMEnt
else
place Gate z=0 x=-99.5 rename=+right parent=SMEnt
place Gate z=0 x=99.5 rename=+left parent=SMEnt
endif
idealsectorbend SM_D angle=-30 fieldCenterRadius=1181-$h fieldInnerRadius=1091.64+6-$h \
fieldOuterRadius=1267.64-6-$h fieldHeight=226 ironInnerRadius=1091.64-$h \
ironOuterRadius=1267.64-$h ironHeight=230 ironColor=1,0,0,.5
idealsectorbend SM_C angle=30 fieldCenterRadius=1181-$h fieldInnerRadius=1091.64+6-$h \
fieldOuterRadius=1267.64-6-$h fieldHeight=226 ironInnerRadius=1091.64-$h \
ironOuterRadius=1267.64-$h ironHeight=230 ironColor=1,0,0,.5
##########################################################
# #
# VALVES, SLITS #
# #
##########################################################
box XSLIT0 height=600 width=600 length=10 material=WKill color=1,0,1,.3
box YSLIT1 height=600 width=600 length=10 material=WKill color=1,0,1,.3
box XSLIT1D height=600 width=600 length=10 material=WKill color=1,0,1,.3
box XSLIT1C height=600 width=600 length=10 material=WKill color=1,0,1,.3
box Slit height=600 width=600 length=10 material=Vacuum color=1,1,1,.3
place Slit z=0 parent=XSLIT0 width=$XSLIT0 height=300 length=10
place Slit z=0 parent=YSLIT1 width=300 height=$YSLIT1 length=10
place Slit z=0 parent=XSLIT1D width=$XSLIT1D height=300 length=10
place Slit z=0 parent=XSLIT1C width=$XSLIT1C height=300 length=10
# Valves+Gates near Separator 1, Separator 2 and after Septum. Radius and length guessed from drawings and pictures. PDFs: VAT 12148-PA44-0001, VAT 16248-PA41-0001, TB21579, TB21607, TB21634. The lengths of the apparatuses are guessed to nearest round figure.
tubs IV innerRadius=125 outerRadius=250 length=100 material=AlKill # Length around 4" which is around 100mm
tubs WV innerRadius=101.5 outerRadius=250 length=100 material=AlKill # Valve body has ID of 250 but 203mm opening is used when in operation.
tubs Bellow innerRadius=(8.66*$inch)/2 outerRadius=((8.66*$inch)/2)+150 length=115 material=AlKill # length little more than 4"
group Vve_Sep1Ent #length=100+100+115=315
place Bellow z=115/2
place WV z=115+100/2
place IV z=115+100+100/2
endgroup
group Vve_Sep1Exit #length=100+100+115=315
place IV z=100/2
place WV z=100+100/2
place Bellow z=100+100+115/2
endgroup
group Vve_Sep2Ent #length=100+100+115=315
place Bellow z=115/2
place XSLIT0 z=115+100/2
place IV z=115+100+100/2
endgroup
group Vve_Sep2Exit #length=100+100+100+100+115=515
place IV z=100/2
place YSLIT1 z=100+100/2
place WV z=100+100+100+100/2 # there's a Slit-Pendulum Valve Valve which is around 4inch~100mm
place Bellow z=100+100+100+100+115/2
endgroup
group Vve_SeptumD #length=100+100+115=315
place XSLIT1D z=100/2
place IV z=100+100/2
place Bellow z=100+100+115/2
endgroup
group Vve_SeptumC #length=100+100+115=315
place XSLIT1C z=100/2
place IV z=100+100/2
place Bellow z=100+100+115/2
endgroup
##########################################################
# #
# Etc. #
# #
##########################################################
#fieldlines exit=0 center=0,0,4622.038 nLines=100 dl=20 radius=500
##########################################################
# #
# PLACEMENT of ELEMENTS #
# #
##########################################################
### T2 TARGET/BEAM SOURCE ###
if $targetNo==0
place T2 z=$s
place Filter z=50
elseif $targetNo==1
place T2Assembly z=0 rotation=Y-90
endif
#place Det z=800 rename=Det_pre_Q1
### QUADRUPOLE MAGNET 1 ###
param s=$s+51.576*$inch
place Quad1box z=$s rename=Quad1 # s=1310.0304
#place Det z=$s+375 rename=Det_post_Q1
### QUADRUPOLE MAGNET 2 ###
param s=$s+25.669*$inch # s=1962.023
place Q2pipe z=$s
place 12q12map z=$s current=($Q2/200)*0.98
place 12Q12 gradient=0 z=$s rename=Quad2
#place Det z=$s+(20.5*$inch)/2 rename=Det_post_Q2
### SCRAPER ###
place Scraper z=3200 ### PLACED WITH ESTIMATE ###
### BENDING DIPOLE MAGNET 1 ###
param s=181.97*$inch # s=4622.038, position of B1. #
place Bend By=($B1/120)*0.125 z=$s rotation=Y-17.5 rename=B1 #### By maximum be 0.125 ####
#place Det z=$s-500 rename=Det_post_Scraper
place BendGate z=$s-450
corner z=$s-300 rotation=Y-17.5 radiusCut=700
corner z=$s+300 rotation=Y-17.5 radiusCut=700
place BendGate2 z=$s+450 ### PLACED WITH ESTIMATE ###
place Bend1_Sep1pipe z=$s+900
#place Det z=$s+900+(34*$inch)/2 rename=Det_post_B1
### SEPARATOR/WIEN FILTER 1 ###
param s=$s+103.78*$inch
place Vve_Sep1Ent z=$s-1960/2-315/2
place Separatorfield z=$s rotation=Z-90 current=$WFBy gradient=$WFVx
place WienRot z=$s rename=Sep1 # s=7258.05, center of Separator S1
place Vve_Sep1Exit z=$s+1960/2+315/2
### QUADRUPOLE MAGNET 3-7 ###
param s=$s+59.208*$inch # s=8761.8832, position of QUAD#3
#place Det z=$s-7*$inch rename=Det_pre_Q3
place 12q12map z=$s current=($Q3/200)*0.98
place 12Q12 gradient=0 z=$s rename=Quad3
#place Det z=$s+10*$inch rename=Det_post_Q3
place Chicago gradient=($Q4/1000)*6.3 z=$s+(43.504-20.985)*$inch rename=Quad4
#place Det z=$s+(43.504-20.985)*$inch+7*$inch rename=Det_post_Q4
#place ChicagoPipe z=$s+43.504*$inch
place Chicago gradient=($Q5/1000)*6.3 z=$s+43.504*$inch rename=Quad5
#place Det z=$s+43.504*$inch+7*$inch rename=Det_post_Q5
place Chicago gradient=($Q6/1000)*6.3 z=$s+(43.504+20.904)*$inch rename=Quad6
#place Det z=$s+(43.504+20.904)*$inch+7*$inch rename=Det_post_Q6
param s=$s+86.969*$inch # s=10,970.8958, position of Quad7
place 12q12map z=$s current=($Q7/200)*0.98
place 12Q12 gradient=0 z=$s rename=Q7
#place Det z=$s+6.1*$inch rename=Det_post_Q7
### SEPARATOR/WIEN FILTER 2 ###
param s=$s+57.484*$inch # s=12430.9894, center of Separator S2
place Vve_Sep2Ent z=$s-1960/2-315/2
place Separatorfield z=$s rotation=Z-90 current=$WFBy gradient=$WFVx
place WienRot z=$s
place Vve_Sep2Exit z=$s+1960/2+515/2
### BENDING DIPOLE MAGNET 2 ###
param s=(181.97+411.30)*$inch ### s=15069.058, position of B2. Measured from z=0. ###
#place Det z=$s-800-(25.2*$inch)/2 rename=Det_post_WF2
place Sep2_Bend2pipe z=$s-800
place Bend By=-($B2/120)*0.125 z=$s rotation=Y17.5 rename=B2
place BendGate z=$s-450 ### NOT PLACED CORRECTLY ###
corner z=$s-300 rotation=Y17.5 radiusCut=99999
corner z=$s+300 rotation=Y17.5 radiusCut=99999
place BendGate2 z=$s+450 ### NOT PLACED CORRECTLY ###
#place Det z=$s+500 rename=Det_post_B2
### QUADRUPOLE MAGNET 8-9 ###
param s=$s+34.246*$inch # s=15938.9064, position of Quad8. Measured from B2. #
place 12q12map z=$s current=($Q8B/200)*0.98 rename=M20Q8F #the fieldmap
place 12dq12map z=$s current=($Q8T/200)*0.98 rotation=z-90 rename=M20VQ8F
place 12Q12 gradient=0 z=$s rename=Quad8
param s=$s+22.087*$inch # s=16499.9162, position of Quad9. Measured from Quad8.
place 12q12map z=$s current=($Q9/200)*0.98 #rotation=Z-90
place 12Q12 gradient=0 z=$s rename=Quad9
place Q8_9pipe z=$s-100
### KICKER ###
param s=15069.058+128.386*$inch # s=18330.0624, position of Kicker. Measured from B2. #
#place Det z=$s-1140 rename=Det_post_Q9
place Kicker z=$s
#place Det z=$s+1140 rename=Det_post_Kicker
### KICKER DRIFT BOX ###
place Driftbox z=21429.1418-37.065*$inch # The kicker is placed with a guess
### SEPTUM MAGNET ###
param s=15069+250.397*$inch # s=21429.1418, position of Septum. Measured from B2. #
place SMEnt z=$s-12.5
place SM_D z=$s x=-99.5 By=($SM1D/775)*0.114
place SM_C z=$s x=99.5 By=-($SM1C/775)*0.114
if $BL==1
### QUADRUPOLE MAGNETS D-Channel 10-12 ### # All position of the Quads10-12 are measured from Septum. #
corner z=$s rotation=Y-30 radiusCut=99999 ### its put at the center of the septum ###
#place Det z=$s+1000-160 rename=Det_post_SeptumD
place Vve_SeptumD z=$s+1000 # Valves with Slit XSLITD placed with a rough estimate from drawing #
#place Det z=$s+1000+160 rename=Det_post_XSLIT1D
place 12q12map z=$s+99.4*$inch current=($Q10DT/200)*0.98 rename=M20Q10DF #the fieldmap
place 12dq12map z=$s+99.4*$inch current=($Q10DB/200)*0.98 rotation=Z-90 rename=M20VQ10DF
place 12Q12 gradient=0 z=$s+99.4*$inch rename=Q10-D # z=23953.9018
place 12q12map z=$s+(99.4+22.126)*$inch current=($Q11DR/200)*0.98 rename=M20Q11DF #the fieldmap
place 12dq12map z=$s+(99.4+22.126)*$inch current=($Q11DL/200)*0.98 rename=M20VQ11DF
place 12Q12 gradient=0 z=$s+(99.4+22.126)*$inch rename=Q11-D # z=24515.9022
place Q10_12pipe z=$s+(99.4+10)*$inch
place 12q12map z=$s+(99.4+22.126*2)*$inch current=($Q12D/200)*0.98
place 12Q12 gradient=0 z=$s+(99.4+22.126*2)*$inch rename=Q12-D # z=25077.9026
#place Det z=$zFocus rename=Det_post_Q12D
place FocalArea z=$zFocus
# This next profile is for gminuit.
profile file=$PROFILE_D particle=mu+ zloop=0:26600:100 z=$zFocus
elseif $BL==2
### QUADRUPOLE MAGNETS C-Channel 10-12 ### # All position of the Quads10-12 are measured from Septum. #
corner z=$s rotation=Y30 radiusCut=99999
place Vve_SeptumC z=$s+1000 # Valves with Slit XSLITC placed with a rough estimate from drawing #
place 12q12map z=$s+99.4*$inch current=($Q10CT/200)*0.98 rename=M20Q10CF #the fieldmap
place 12dq12map z=$s+99.4*$inch current=($Q10CB/200)*0.98 rotation=Z-90 rename=M20VQ10CF
place 12Q12 gradient=0 z=$s+99.4*$inch rename=Q10-C
place 12q12map z=$s+(99.4+22.126)*$inch current=($Q11CR/200)*0.98 rename=M20Q11CF #the fieldmap
place 12dq12map z=$s+(99.4+22.126)*$inch current=($Q11CL/200)*0.98 rename=M20VQ11CF
place 12Q12 gradient=0 z=$s+(99.4+22.126)*$inch rename=Q11-C
place Q10_12pipe z=$s+(99.4+10)*$inch
place 12q12map z=$s+(99.4+22.126*2)*$inch current=($Q12C/200)*0.98
place 12Q12 gradient=0 z=$s+(99.4+22.126*2)*$inch rename=Q12-C
#place Det z=$zFocus rename=Det#
#place FocalArea z=$zFocus
# This next profile is for gminuit.
profile file=$PROFILE_C particle=mu+ zloop=0:26600:100 z=$zFocus
else
corner z=$s rotation=Y-30 radiusCut=99999 ### put at the position of the septum ###
#place Det z=$s+1000-160 rename=Det_post_SeptumD
place Vve_SeptumD z=$s+1000
#place Det z=$s+1000+160 rename=Det_post_XSLIT1D
place 12q12map z=$s+99.4*$inch current=($Q10DT/200)*0.98 rename=M20Q10DF #the fieldmap
place 12dq12map z=$s+99.4*$inch current=($Q10DB/200)*0.98 rotation=Z-90 rename=M20VQ10DF
place 12Q12 gradient=0 z=$s+99.4*$inch rename=Q10-D # z=23953.9018
place 12q12map z=$s+(99.4+22.126)*$inch current=($Q11DR/200)*0.98 rename=M20Q11DF #the fieldmap
place 12dq12map z=$s+(99.4+22.126)*$inch current=($Q11DL/200)*0.98 rename=M20VQ11DF
place 12Q12 gradient=0 z=$s+(99.4+22.126)*$inch rename=Q11-D # z=24515.9022
place Q10_12pipe z=$s+(99.4+10)*$inch
place 12q12map z=$s+(99.4+22.126*2)*$inch current=($Q12D/200)*0.98
place 12Q12 gradient=0 z=$s+(99.4+22.126*2)*$inch rename=Q12-D # z=25077.9026
#place Det z=$zFocus rename=Det_post_Q12D
place FocalArea z=$zFocus
corner z=$s rotation=Y60 radiusCut=99999
place Vve_SeptumC z=$s+1000
place 12q12map z=$s+99.4*$inch current=($Q10CT/200)*0.98 rename=M20Q10CF #the fieldmap
place 12dq12map z=$s+99.4*$inch current=($Q10CB/200)*0.98 rotation=Z-90 rename=M20VQ10CF
place 12Q12 gradient=0 z=$s+99.4*$inch rename=Q10-C
place 12q12map z=$s+(99.4+22.126)*$inch current=($Q11CR/200)*0.98 rename=M20Q11CF #the fieldmap
place 12dq12map z=$s+(99.4+22.126)*$inch current=($Q11CL/200)*0.98 rename=M20VQ11CF
place 12Q12 gradient=0 z=$s+(99.4+22.126)*$inch rotation=Z-90 rename=Q11-C
place Q10_12pipe z=$s+(99.4+10)*$inch
place 12q12map z=$s+(99.4+22.126*2)*$inch current=($Q12C/200)*0.98
place 12Q12 gradient=0 z=$s+(99.4+22.126*2)*$inch rename=Q12-C
#place Det z=$zFocus rename=Det#
place FocalArea z=$zFocus
endif
# The profile 'virtualdetector' is placed at the positions of the virtualdetectors all of which have been commented out. This is to check the beam status before and after each element. Here, the doublet Q8-Q9 and the triplet Q10-Q12 are treated as single element.
#profile file=virtualdetector.txt particle=mu+ z=800,1685.03,2222.37,4122.04,5915.74,8584.08,9015.88,9511.67,10044.68,10575.65,11125.84,13949.02,15569.06,17190.06,19470.06,22269.14,22589.14,26278