Attachment '20220112_SQB1_third_version.f90'
Download 1 program KAGRA_Injection_2021
2 !!!This program simulates an optical system for SQB1
3 use optocad
4 use rsplot
5 !!!Variable declaration, allocation, and initialization
6 character(len=400),dimension(90,18) :: K = ''
7 real(8) :: pst1(7) = (/2.,1.,1.e-1,1.e-2,1.e-3,1.e-4,1.e-5/),&
8 pst2(7) = (/2.e-1,1.e-1,1.e-2,1.e-3,1.e-4,1.e-5,1.e-6/),&
9 shift2, shift3, angle
10 integer i, j
11 !!!Character arrays describe optical components (all components in this program, NO call external files)
12
13 !!!OMC chamber
14 ! ac, x, y, rd, ag, c, m r t # label
15 K(01,1)='o -143.289, 436.736' ! origin at Input Faraday Isolator crystals in OMC
16 K(02,1)='b 0.+shift3, 94.649, w=0.9047, ag=-90., p=1.0,z=0., rc=1.2809e3 # Input beam1 >>'
17 K(03,1)='d t, 0.+shift3, 94.649, rd=12.7, ag=-90.-56., m=1, t=1. # OMC PBS1 >> @lc 20,0' !Index_1.4496
18 K(04,1)='+ t, dm=6.35 '
19 K(05,1)='d t, 0., 0., 6., ag=-90., m=5, r=0., t=1. # OMC FI >> @lc20,0' !TGG crystal Index_1.954
20 K(06,1)='+ t, dm= 10.,'
21 K(07,1)='d t, 0, -48.503, rd=10. , ag=-90., m=1, r=0., t=1. # OMC HWP >> @lc20,0'
22 K(08,1)='+ t, dm=7.13'
23 K(09,1)='d srt, 0, -127.518, rd=12.7, ag=-45., m=1, r=0.99, t=0.01 # OMC PBS2 >> @lc 20,0' !Index_1.4496
24 K(10,1)='+ t, dm=6.35'
25 K(11,1)='c d, 2.789-.8079, -900. , 50.8, ag=90. #OMC Damper > @lc-100,-30' !Beam damper
26 K(12,1)='d r, -93.335-6.6366, -128.524+1.0058, 12.7, ag=-135., m=1, r=1., t=0. #OMC RM1 > @lc-65,-40' !Plane Mirror
27 K(13,1)='+ d, dm=5'
28 !!!OMMT chamber
29 K(14,1)='d r, -193.335-6.6366+100., 2358.649+59.459, 12.7*2, ag=135., m=1, r=1., t=0. #OMMT RM1 > @lc-100,40' !Plane Mirror
30 K(15,1)='+ d, dm=5'
31 K(16,1)='d r, 550.0284, 2358.649+59.459, 12.7*2, ag=-45., m=1, r=1., t=0. #OMMT RM2 > @lc-170,40' !Plane Mirror
32 K(17,1)='+ d, dm=5'
33 K(18,1)='d r, 550.0284, 2358.649+259.459, 12.7*2, ag=135., m=1, r=1., t=0. #OMMT RM3 > @lc-170,40' !Plane Mirror
34 K(19,1)='+ d, dm=5'
35 !!!SQB1a and SQB1b chamber
36 K(20,1)='d r, 3795.0444, 2618.108, 25.4*2, ag=-3., m=1,rc=-7.1e3, r=1.,t=0. #SQB1a RSM1 > @cc-80,70'
37 K(21,1)='+ d, dm=5'
38 K(22,1)='d r, 2809.553, 2721.6862, 25.4*2,ag=-3+180., m=1, rc=-7.1e3,r=1.,t=0. #SQB1a RSM2 > @cc10,70'
39 K(23,1)='+ d, dm=5'
40 K(24,1)='d r, 2809.553+940., 2721.6862, 12.7*2, ag=-45., m=1, r=1., t=0. #SQB1a RM1 > @lc-200,30' !Plane Mirror
41 K(25,1)='+ d, dm=5'
42 K(26,1)='d str, 2809.553+940., 2821.6862+102.481, rd=12.7, ag=45., m=1, r=.01,t=.99 # SQB1a PBS1 >> @lc 20,0' !Index_1.4496
43 K(27,1)='+ t, dm=6.35 '
44 K(28,1)='d t, 2809.553+941.982, 2821.6862+174.367, rd=10. , ag=90., m=1, r=0., t=1. # SQB1a HWP1 >> @lc20,0'
45 K(29,1)='+ t, dm=7.13'
46 K(30,1)='d t, 2809.553+941.982, 2821.6862+220., 6., ag=90., m=5, r=0., t=1. # SQB1a FI1 >> @lc20,0' !TGG crystal Index_1.954
47 K(31,1)='+ t, dm=10'
48 K(32,1)='d tr, 2809.553+941.982, 2821.6862+323.051, rd=12.7, ag=135., m=1, r=.01,t=.99 # SQB1a PBS2 >> @lc20,0' !Index_1.4496
49 K(33,1)='+ t, dm=6.35'
50 K(34,1)='d t, 2809.553+940, 2821.6862+2233.051, rd=10., ag=90., m=1, r=0., t=1. # SQB1b HWP1 >> @lc-220,20'
51 K(35,1)='+ t, dm=7.13'
52 K(36,1)='d rr, 2809.553+938.019, 2821.686+2665.126, rd=12.7, ag=45., m=1, t=0.,r=1. # SQB1b PBS2 >> @lc 20,0' !Index_1.4496
53 K(37,1)='+ t, dm=6.35'
54 K(38,1)='d t, 2809.553+938.019, 5389.705, 6., ag=90., m=5, r=0., t=1. # SQB1b FI1 >> @lc20,0' !TGG crystal Index_1.954
55 K(39,1)='+ t, dm=10'
56 K(40,1)='d t, 2809.553+938.019, 2821.686+2615.633, rd=10. , ag=90., m=1, r=0., t=1. # SQB1b HWP2 >> @lc20,0'
57 K(41,1)='+ t, dm=7.13'
58 K(42,1)='d tstr, 2809.553+940.001, 2821.686+2475.352, rd=12.7, ag=135., m=1, r=.01,t=.99 # SQB1b PBS1 >> @lc20,20' !Index_1.4496
59 K(43,1)='+ t, dm=6.35'
60 K(44,1)='d rrtt, 2809.553+638.019, 2821.686+2665.126, 12.7*4, ag=-135., m=1, r=0.99, t=0.01 # SQB1b DM @cb100,10'
61 K(45,1)='+ t, dm=5'
62 K(46,1)='c t, 2809.553+638.019, 6252, rd=25.4, ag=90., m=1, r=0., t=1. # SQB1b Test1 >> @lc20,50' !Index_1.4496
63 K(47,1)='d r, 2809.553+715.001, 2821.686+2477.332, rd=12.7*2, ag=135., m=1, r=1.,t=0. # SQB1b RM1 >> @lc-210,0' !Index_1.4496
64 K(48,1)='+ d, dm=5 '
65 K(49,1)='c t, 2809.553+715.001, 2103, rd=25.4, ag=90., m=1, r=0., t=1. # SQB1a Test1 >> @lc20,-30' !Index_1.4496
66 K(50,1)='d r, 2809.553+615.001, 2821.6862+325.031, rd=12.7*2, ag=135., m=1, r=1.,t=0. # SQB1a RM2 >> @lc-100,40' !Index_1.4496
67 K(51,1)='+ d, dm=5 '
68 K(52,1)='c t, 2609.553+815.001362, 2088, rd=25.4, ag=90., m=1, r=0., t=1. # SQB1a Test2 >> @lc-100,-30' !Index_1.4496
69 K(53,1)='d str, 2809.553+515.001, 2821.6862+102.481, rd=12.7, ag=135., m=1, r=.5,t=.5 # SQB1a BS1 >> @lc-110,40' !Index_1.4496
70 K(54,1)='+ t, dm=6.35'
71 K(55,1)='c d, 2809.553+515.001, 2821.6862+62.481, 41., ag=90. #SQB1a QPD1 > @lc-260,-10'
72 K(56,1)='c d, 2809.553+15.001, 2821.6862+102.481+1.981, 41., ag=0. #SQB1a QPD2 > @lc-100,60'
73 K(57,1)='c, t,870.0284, 2358.649+259.459, rd=25.4, ag=0., m=1, r=0., t=1. # OMMT Test >> @lc10,30' !Index_1.4496 '
74 !!!Simulate FDS signal path
75 K(01,2)='o -143.289+2809.553+638.019, 436.736+6252' ! origin at SQB1b Test1's output for FDS
76 K(02,2)='b 0, 0. , w=4.03498, ag=-90., p=1.0 ,z=0., rc=2322 # Input beam2 t plane >>'
77 !K(01,5)='o -143.289+2809.553+375.944, 436.736+3654.1672' ! origin at SQB1 Test2's output for FDS
78 !K(02,5)='b 0, 0. , w=5.29208 , ag=-90., p=1.0 ,z=0., rc=2050 # Input beam3 s plane >>'
79
80 !!! Simulate Green light from ESQB1
81 K(01,3)='o 2900, 2450' ! origin at SQB1 Test2's output for FDS
82 K(02,3)='b 0, 0. , w=1.0, ag=90., p=1.0 ,z=0., rc=1.0e6 # Input beam3 >>'
83 K(03,3)='c t, 0, 126., rd=25.4, ag=90., m=1, r=0., t=1. # SQB1a Test3 >> @lc-100,-50' !Index_1.4496
84 K(04,3)='d srt, 0, 3350., rd=12.7*2, ag=135., m=1, r=0.99,t=0.01 # SQB1b GM2 >> @lc -200,-30' !Index_1.4496
85 K(05,3)='+ t, dm=5 '
86 K(06,3)='d r, 400+2.723, 3350., rd=12.7*2, ag=-45., m=1, r=1.0,t=0. # SQB1b GM3 >> @lc -200,30' !Index_1.4496
87 K(07,3)='+ d, dm=5 '
88 K(08,3)='d str, -1.560, 3550-1.560, rd=12.7, ag=45., m=1, r=.5,t=.5 # SQB1b BS1 >> @lc0,40' !Index_1.4496
89 K(09,3)='+ t, dm=6.35'
90 K(10,3)='c d, 0, 3950, 41., ag=90. #SQB1b QPD1 > @lc0,-50'
91 K(11,3)='c d, -100, 3550-1.560, 41., ag=0. #SQB1b QPD2 > @lc-260,-10'
92 !K(08,3)='d str, 440.-1.981-1.981, 450.-1.981, rd=12.7, ag=-45., m=1, r=.99,t=.01 # SQB1a M14 >> @lc -100,60' !Index_1.4496
93 !K(09,3)='+ t, dm=6.35 '
94 !K(10,3)='c d, 530.-3.962+, 450.-1.981-1.981, 41., ag=0. # SQB1a Near2 >> @lc-100,-60' !Index_1.4496
95 !K(11,3)='d r, 440.-1.981-1.981, 550.-1.981, rd=12.7, ag=45., m=1, r=.99,t=.01 # SQB1a M15 >> @lc -100,30' !Index_1.4496
96 !K(12,3)='+ t, dm=6.35 '
97 !K(13,3)='c d, -110.-3.962, 550.-1.981, 41., ag=0. # SQB1a Far2 >> @lc-150,60' !Index_1.4496
98 !K(14,3)='c d, 130.+1.981+1.981, 450.-80., 41., ag=90. # SQB1a Near3 >> @lc-100,-80' !Index_1.4496
99
100 !!!Simulate the Green light from Filter cavity
101 !K(01,4)='o -143.289+2809.553+375.944, 436.736+3654.1672' ! origin at SQB1 Test2's output for FDS
102 !K(02,4)='b 0, 0. , w=1.17316, ag=-90., p=1.0 ,z=0., rc=6.624e3 # Input beam4 t plane >>'
103
104 !!! Simulate GW signal path
105 K(01,5)='o -143.289, 436.736' ! origin at Output Faraday Isolator
106 K(02,5)='b 0.+shift3, 94.649, w=0.9047, ag=90., p=1.0 ,z=0., rc=1.2809e3 # Input beam5 >>'
107 K(03,5)='d r, -25.0992+10*shift3, 2573.424 ,rd=50, ag=90.+3.07, m=1, r=1. #OMMT Mirror > @lc-100,100'
108 K(04,5)='+ d, dm=57.817'
109 K(05,5)='d r, 310.212-10*shift3, -25.088,rd=50, ag=-90.+3.07, m=1, r=1. #OMC Mirror > @lc60,-25'
110 K(06,5)='+ d, dm=57.817'
111
112 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
113 !!!!OptoCad initialization and setups
114 call oc_init(unit='mm',pform='A0_P')
115 call oc_frame(xbeg=-800., ybeg=-800., yend=10000., xend=1700.,xbl=100.,ybl=70., scale=.14,ax='XXYY',glp=0, unit='x')
116 call oc_frame(xbeg=1700., ybeg=2400., yend=10000., xend=4000.,xbl=500.,ybl=518., scale=.14,ax='XXYY',glp=0, unit='x')
117 !call oc_frame(xbeg=1700., ybeg=2400., yend=5000., xend=4000.,xbl=500.,ybl=518., scale=.14,ax='XXYY',glp=0, unit='x')
118 !call oc_frame(xbeg=-700., ybeg=2400., yend=3700., xend=1000.,xbl=100.,ybl=650., scale=.35,ax='XXYY',glp=0)
119 call oc_set(n=1, fslb=4., &
120 rix=(/1.4496_8, 1.4285_8, 1.5066_8, 2.2321, 1.954_8/), &
121 ! Fused Silica, CaF2, N-BK7, LiNb, TGG? (https://refractiveindex.info/)
122 print='rs s2 act w2t_6 w2s_6 z2t z2s R2t R2s an1 an2 pw ipp tpp lb dx dy ang med')
123
124
125 !!!!PS PLOT drawing section
126 call ps_quad(frm= 1,xu=-780.,yu=1000.,wu=200.,hu=70.,fill=14) !OMC optical bench
127 call ps_text(frm=1,str='OMC',xu=-780.,yu=1000.,fs=10., fn='HB')
128 call ps_arc(frm=1,xu= 0.,yu= 0.,rux= 649., lw=2., ci= 4) !OMC optical bench size
129 call ps_arc(frm=1,xu= 0.,yu= 0.,rux= 793.5, lw=2., ci= 1) !OMC chamber size
130 call ps_quad(frm=1,xu=-680.,yu=3800.,wu=250.,hu=70.,fill=14) !OMMT optical bench
131 call ps_text(frm=1,str='OMMT',xu=-680.,yu=3800.,fs=10., fn='HB')
132 call ps_arc(frm=1,xu=96.1,yu= 3014.1,rux= 499., lw=2., ci= 4) !OMMT optical bench size
133 call ps_arc(frm=1,xu=96.1,yu= 3014.1,rux= 632.8, lw=2., ci= 1) !OMMT chamber size
134 call ps_quad(frm=1,xu=-143.289-173.,yu=436.736-5-185.5,wu=346.,hu=371.,fill=14) !OMC_OFI
135 call ps_quad(frm=2,xu=2809.553+940-125-143.289,yu=2821.6862+436.736,wu=250.,hu=450.,fill=14) !SQB1a_OFI1
136 call ps_quad(frm=2,xu=2809.553+940-125-143.289,yu=3821.6862+1786.736,wu=250.,hu=450.,fill=14) !SQB1a_OFI2
137 call ps_line(frm=1,xu1=2585.,yu1=3000.,xu2=730.,yu2=3000.,lw=2., ci= 1)
138 call ps_line(frm=1,xu1=2400.,yu1=3400.,xu2=600.,yu2=3400.,lw=2., ci= 1)
139 call ps_arrow(1,10.,1000.,600.,pi/2,wd=.5,sag=.1,ls=100.,lw=.5,fill=2)
140 call ps_arrow(1,10.,1600.,0.,0.,wd=.5,sag=.1,ls=100.,lw=.5,fill=2)
141 call ps_text(frm=1,str='Y',xu=1000.,yu=600.,fs=10.,fn='HB')
142 call ps_text(frm=1,str='X',xu=1600.,yu=0.,fs=10.,fn='HB')
143 call ps_text(frm=1,str='Squeezing Injection Part',xl=0.,yt=25.,fs=15.,fn='HB')
144 call ps_line(frm=1,xu1=30.-398.4,yu1=2590.,xu2=30-398.4,yu2=850.-150.,lw=2., ci= 1)
145 call ps_line(frm=1,xu1=30.+398.4,yu1=2480.,xu2=30+398.4,yu2=850.-180.,lw=2., ci= 1)
146 call ps_line(frm=2,xu1=2470.,yu1=3000.,xu2=730.,yu2=3000.,lw=2., ci= 1)
147 call ps_line(frm=2,xu1=2400.,yu1=3400.,xu2=600.,yu2=3400.,lw=2., ci= 1)
148 call ps_arc(frm=2,xu=96.1+3100.,yu= 3014.1+299.,rux= 649., lw=2., ci= 4) !SQB1a optical bench size
149 call ps_arc(frm=2,xu=96.1+3100.,yu= 3014.1+299.,rux= 793.5, lw=2., ci= 1) !SQB1a chamber sizecall
150 call ps_arc(frm=2,xu=96.1+3100.,yu= 3014.1+299.,rux= 793.5, lw=2., ci= 1) !SQB1a chamber size
151 call ps_quad(frm=2,xu=1750.,yu=4000.,wu=280.,hu=70.,fill=14) !SQB1a optical bench
152 call ps_text(frm=2,str='SQB1a',xu=1750.,yu=4000.,fs=10., fn='HB')
153 call ps_quad(frm=2,xu=2809.553+515.001-143.289-41,yu=2821.6862+62.481+436.736-67,wu=82.,hu=67.,fill=14) !SQB1a QPD1 Near
154 call ps_quad(frm=2,xu=2809.553+15.001-143.289-67,yu=2821.6862+102.481+1.981+436.736-41,wu=67.,hu=82.,fill=14) !SQB1a QPD1 Far
155 call ps_arc(frm=2,xu=96.1+3100.,yu= 3014.1+99.+1200.+793.5+793.5,rux= 649., lw=2., ci= 4) !SQB1b optical bench size
156 call ps_arc(frm=2,xu=96.1+3100.,yu= 3014.1+99.+1200.+793.5+793.5,rux= 793.5, lw=2., ci= 1) !SQB1b chamber size
157 call ps_quad(frm=2,xu=1750.,yu=4000.+200.+793.5+793.5,wu=280.,hu=70.,fill=14) !SQB1b optical bench
158 call ps_text(frm=2,str='SQB1b',xu=1750.,yu=4000.+200.+793.5+793.5,fs=10., fn='HB')
159
160 !!!!"Shift" is a value used to compensate beam path shift caused by any optics.
161 call oc_bind('shift',shift(45._8,1.4496_8,0.5_8))
162 call oc_bind('shift2',shift(45._8,1.4496_8,0.5_8)*tan(45._8*3.14159265358979_8/180._8))
163 call oc_bind('shift3',shift(56._8,1.4496_8,6.35_8)) ! PBS1&PBS2 of OFI
164 call oc_bind('angle',atan(1.1_8/28.1_8)*180._8/3.14159265358979_8)
165 !!!!Beam color for 1.IR and 2.green.
166 ico1=ps_color((/1.,.8,.8/)) ! Plot outer part color ...
167 icm1=ps_color((/1.,.6,.6/)) ! ... middle part color ...
168 ici1=ps_color((/1.,.0,.0/)) ! .... inner part color of the beam
169 ico2=ps_color((/.8,1.,.8/))
170 icm2=ps_color((/.6,1.,.6/))
171 ici2=ps_color((/.0,1.,.0/))
172 !!!!Reading optical components data from character array K and plotting surfaces.
173 do i = 1, 18
174 call oc_input(K(:,i))
175 end do
176 call oc_trace(1, of='KAGRA.log') ! Save > or >> value in .log file
177 call oc_beam(3.0,fill=ico1,pst=pst1) ! Plot 3* waist range
178 call oc_beam(2.0,fill=icm1,pst=pst1) ! ... 2* waist ...
179 call oc_beam(1.0,fill=ici1,pst=pst1) ! ... 1* waist ...
180 call oc_beam(.0,1,5,.1,rnc=4,rns=.1) ! Plot the beam axes in dash--black
181 call oc_surf(lw=1.) ! Plot all surface with linewidth 0.2mm
182
183
184 call oc_trace(2, of='KAGRA.log')
185 call oc_beam(3.0,fill=ico1,pst=pst1) ! Plot 3* waist range
186 call oc_beam(2.0,fill=icm1,pst=pst1) ! ... 2* waist ...
187 call oc_beam(1.0,fill=ici1,pst=pst1) ! ... 1* waist ...
188 call oc_beam(.0,1,5,.1,rnc=4,rns=.1) ! Plot the beam axes in dash--black
189 call oc_surf(lw=1.) ! Plot all surface with linewidth 0.2mm
190
191 call oc_trace(3, of='KAGRA.log')
192 call oc_beam(3.0,fill=ico2,pst=pst2) ! Plot 3* waist range
193 call oc_beam(2.0,fill=icm2,pst=pst2) ! ... 2* waist ...
194 call oc_beam(1.0,fill=ici2,pst=pst2) ! ... 1* waist ...
195 call oc_beam(.0,1,5,.1,rnc=4,rns=.1) ! Plot the beam axes in dash--black
196 call oc_surf(lw=1.) ! Plot all surface with linewidth 0.2mm
197
198 call oc_trace(4, of='KAGRA.log')
199 call oc_beam(3.0,fill=ico1,pst=pst1) ! Plot 3* waist range
200 call oc_beam(2.0,fill=icm1,pst=pst1) ! ... 2* waist ...
201 call oc_beam(1.0,fill=ici1,pst=pst1) ! ... 1* waist ...
202 call oc_beam(.0,1,5,.1,rnc=4,rns=.1) ! Plot the beam axes in dash--black
203 call oc_surf(lw=1.) ! Plot all surface with linewidth 0.2mm
204
205 !call oc_trace(4, of='KAGRA.log')
206 !call oc_beam(3.0,fill=ico2,pst=pst2) ! Plot 3* waist range
207 !call oc_beam(2.0,fill=icm2,pst=pst2) ! ... 2* waist ...
208 !call oc_beam(1.0,fill=ici2,pst=pst2) ! ... 1* waist ...
209 !call oc_beam(.0,1,5,.1,rnc=4,rns=.1) ! Plot the beam axes in dash--black
210 !call oc_surf(lw=1.) ! Plot all surface with linewidth 0.2mm
211
212
213 !!!!Exit OptoCAD (Good job!)
214 call oc_exit
215 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
216 contains
217 real(8) function shift(theta1, n, a)
218 !!!This function calculates beam path shift
219 !!!caused by optical components (mainly beam splitter) of refractive index n.
220 !!!theta1[deg]: angle of incidence (normally 45 deg.).
221 !!!theta2[rad]: angle of refraction (normally less than theta1).
222 !!!n: refractive index of a component.
223 !!!a[cm]: thickness of a component.
224 real(8), intent(in) :: theta1, n, a
225 real(8) pi, theta2
226 pi = 3.14159265358979_8
227 theta2 = asin(sin(theta1*pi/180._8)/n)
228 shift = a*sin(theta1*pi/180._8) - a*cos(theta1*pi/180._8)*tan(theta2)
229 end function shift
230 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
231 end program KAGRA_Injection_2021
Attached Files
To refer to attachments on a page, use attachment:filename, as shown below in the list of files. Do NOT use the URL of the [get] link, since this is subject to change and can break easily.You are not allowed to attach a file to this page.