source: PSPA/madxPSPA/src/trrun_coll.f90 @ 476

Last change on this file since 476 was 430, checked in by touze, 11 years ago

import madx-5.01.00

File size: 7.6 KB
Line 
1subroutine trcoll(flag, apx, apy, turn, sum, part_id, last_turn,  &
2     last_pos, last_orbit, z, ntrk,al_errors,offx,offy)
3
4  use twiss0fi
5  use name_lenfi
6  implicit none
7
8  !----------------------------------------------------------------------*
9  ! Purpose:                                                             *
10  !   test for collimator aperture limits.                               *
11  ! input:                                                               *
12  !   flag      (integer)   aperture type flag:                          *
13  !                         1: elliptic, 2: rectangular                  *
14  !   apx       (double)    x aperture or half axis                      *
15  !   apy       (double)    y aperture or half axis                      *
16  !   turn      (integer)   current turn number.                         *
17  !   sum       (double)    accumulated length.                          *
18  ! input/output:                                                        *
19  !   part_id   (int array) particle identification list                 *
20  !   last_turn (int array) storage for number of last turn              *
21  !   last_pos  (dp. array) storage for last position (= sum)            *
22  !   last_orbit(dp. array) storage for last orbit                       *
23  !   z(6,*)    (double)    track coordinates: (x, px, y, py, t, pt).    *
24  !   ntrk      (integer) number of surviving tracks.                    *
25  !----------------------------------------------------------------------*
26  integer flag,turn,part_id(*),last_turn(*),ntrk,i,n,nn
27  double precision apx,apy,sum,last_pos(*),last_orbit(6,*),z(6,*),  &
28       one,al_errors(align_max),offx,offy
29  parameter(one=1d0)
30  character(name_len) aptype
31
32  n = 1
3310 continue
34  do i = n, ntrk
35
36     !---- Is particle outside aperture?
37     if (flag .eq. 1.and.((z(1,i)-al_errors(11)- offx)/apx)**2             &
38          +((z(3,i)-al_errors(12)- offy) / apy)**2 .gt. one) then
39        go to 99
40     else if(flag .eq. 2                                             &
41          .and. (abs(z(1,i)-al_errors(11)- offx) .gt. apx                         &
42          .or. abs(z(3,i)-al_errors(12)- offy) .gt. apy)) then
43        go to 99
44        !***  Introduction of marguerite : two ellipses
45     else if(flag .eq. 3.and. ((z(1,i)-al_errors(11)- offx) / apx)**2      &
46          + ((z(3,i)-al_errors(12)- offy) / apy)**2 .gt. one .and.                &
47          ((z(1,i)-al_errors(11)- offx) / apy)**2 +                               &
48          ((z(3,i)-al_errors(12)- offy) / apx)**2 .gt. one) then
49        go to 99
50     endif
51     go to 98
5299   n = i
53     nn=name_len
54     call node_string('apertype ',aptype,nn)
55     call trkill(n, turn, sum, ntrk, part_id,                        &
56          last_turn, last_pos, last_orbit, z, aptype)
57     goto 10
5898   continue
59  enddo
60end subroutine trcoll
61
62subroutine trcoll1(flag, apx, apy, turn, sum, part_id, last_turn,  &
63     last_pos, last_orbit, z, ntrk,al_errors, apr,offx,offy)
64
65  use twiss0fi
66  use name_lenfi
67  implicit none
68
69  !----------------------------------------------------------------------*
70  ! Similar with trcoll, for racetrack type aperture
71  !-------------Racetrack type , added by Yipeng SUN 21-10-2008---
72  !----------------------------------------------------------------------*
73  integer flag,turn,part_id(*),last_turn(*),ntrk,i,n,nn
74  double precision apx,apy,sum,last_pos(*),last_orbit(6,*),z(6,*),  &
75       one,al_errors(align_max),apr,offx,offy
76  parameter(one=1d0)
77  character(name_len) aptype
78
79  n = 1
8010 continue
81  do i = n, ntrk
82
83     !---- Is particle outside aperture?
84     if (flag .eq. 4                                                 &
85          .and. (abs(z(1,i)-al_errors(11)- offx)) .gt. (apr+apx)                  &
86          .or. abs(z(3,i)-al_errors(12)- offy) .gt. (apy+apr) .or.                &
87          ((((abs(z(1,i)-al_errors(11)- offx)-apx)**2+                            &
88          (abs(z(3,i)-al_errors(12)- offy)-apy)**2) .gt. apr**2)                  &
89          .and. (abs(z(1,i)-al_errors(11)- offx)) .gt. apx                        &
90          .and. abs(z(3,i)-al_errors(12)- offy) .gt. apy)) then
91        go to 99
92     endif
93     go to 98
9499   n = i
95     nn=name_len
96     call node_string('apertype ',aptype,nn)
97     call trkill(n, turn, sum, ntrk, part_id,                        &
98          last_turn, last_pos, last_orbit, z, aptype)
99     goto 10
10098   continue
101  enddo
102end subroutine trcoll1
103subroutine trkill(n, turn, sum, jmax, part_id,                    &
104     last_turn, last_pos, last_orbit, z, aptype)
105
106  use name_lenfi
107  implicit none
108
109  !hbu--- kill particle:  print, modify part_id list
110  logical recloss
111  integer i,j,n,turn,part_id(*),jmax,last_turn(*),get_option
112  double precision sum, z(6,*), last_pos(*), last_orbit(6,*),       &
113       torb(6) !, theta, node_value, st, ct, tmp
114  character(name_len) aptype
115  !hbu
116  character(name_len) el_name
117
118  recloss = get_option('recloss ') .ne. 0
119
120  !!--- As elements might have a tilt we have to transform back
121  !!--- into the original coordinate system!
122  !      theta = node_value('tilt ')
123  !      if (theta .ne. 0.d0)  then
124  !          st = sin(theta)
125  !          ct = cos(theta)
126  !!--- rotate trajectory (at exit)
127  !            tmp = z(1,n)
128  !            z(1,n) = ct * tmp - st * z(3,n)
129  !            z(3,n) = ct * z(3,n) + st * tmp
130  !            tmp = z(2,n)
131  !            z(2,n) = ct * tmp - st * z(4,n)
132  !            z(4,n) = ct * z(4,n) + st * tmp
133  !      endif
134
135  last_turn(part_id(n)) = turn
136  last_pos(part_id(n)) = sum
137  do j = 1, 6
138     torb(j) = z(j,n)
139     last_orbit(j,part_id(n)) = z(j,n)
140  enddo
141
142  !hbu
143  call element_name(el_name,len(el_name))
144  !hbu
145  write(6,'(''particle #'',i6,'' lost turn '',i6,''  at pos. s ='', f10.2,'' element='',a,'' aperture ='',a)') &
146       part_id(n),turn,sum,el_name,aptype
147  print *,"   X=",z(1,n),"  Y=",z(3,n),"  T=",z(5,n)
148
149  if(recloss) then
150     call tt_ploss(part_id(n),turn,sum,torb,el_name)
151  endif
152
153  do i = n+1, jmax
154     part_id(i-1) = part_id(i)
155     do j = 1, 6
156        z(j,i-1) = z(j,i)
157     enddo
158  enddo
159  jmax = jmax - 1
160
161end subroutine trkill
162subroutine tt_ploss(npart,turn,spos,orbit,el_name)
163
164  use name_lenfi
165  implicit none
166
167  !hbu added spos
168  !----------------------------------------------------------------------*
169  !--- purpose: enter lost particle coordinates in table                 *
170  !    input:                                                            *
171  !    npart  (int)           particle number                            *
172  !    turn   (int)           turn number                                *
173  !    spos    (double)       s-coordinate when loss happens             *
174  !    orbit  (double array)  particle orbit                             *
175  !    orbit0 (double array)  reference orbit                            *
176  !----------------------------------------------------------------------*
177  integer npart,turn,j
178  double precision orbit(6),tmp,tt,tn
179  double precision energy,get_value
180  character(36) table
181  character(name_len) el_name
182  !hbu
183  double precision spos
184  !hbu
185  character(4) vec_names(7)
186  !hbu
187  data vec_names / 'x', 'px', 'y', 'py', 't', 'pt', 's' /
188  data table / 'trackloss' /
189
190  tn = npart
191  tt = turn
192
193  energy = get_value('probe ','energy ')
194
195
196  ! the number of the current particle
197  call double_to_table_curr(table, 'number ', tn)
198  ! the number of the current turn
199  call double_to_table_curr(table, 'turn ', tt)
200  !hbu spos
201
202  call double_to_table_curr(table,vec_names(7),spos)
203
204  do j = 1, 6
205     tmp = orbit(j)
206     call double_to_table_curr(table, vec_names(j), tmp)
207  enddo
208  call double_to_table_curr(table, 'e ', energy)
209  call string_to_table_curr(table, 'element ', el_name)
210
211  call augment_count(table)
212end subroutine tt_ploss
Note: See TracBrowser for help on using the repository browser.