source: trunk/fsource/pydiffax.for @ 2236

Last change on this file since 2236 was 2236, checked in by vondreele, 7 years ago

allow delete last layer & fraction input for stacking probability
File size: 9.9 KB
Line 
1      SUBROUTINE PYLOADSCF(NATP,ATYPES,SFDAT,DEBG)
2       
3Cf2py intent(in) NATP
4Cf2py intent(in) ATYPES
5Cf2py intent(in) SFDAT
6cf2py depend(NATP) ATYPES,SFDAT
7cf2py intent(in) DEBG
8           
9      INCLUDE 'DIFFaXsubs/DIFFaX.par'
10      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
11               
12      INTEGER*4 NATP,I,J
13      CHARACTER*4 ATYPES(NATP)
14      REAL*4  SFDAT(9,NATP)
15      LOGICAL DEBG
16               
17C fill common x-ray scattering factors
18      debug = DEBG
19      DO J=1,NATP
20        WRITE(atom_l(J),'(A4)') ATYPES(J)
21        DO I=1,9
22          x_sf(I,J) = SFDAT(I,J)
23        END DO
24        if (debug) print '(1x,a4,9f10.6)',ATYPES(J),(x_sf(I,J),I=1,9)
25      END DO
26      intp_F = .TRUE.
27      n_atoms = NATP
28      RETURN
29      END
30       
31      SUBROUTINE PYGETINST(LAMB,TTHMIN,TTHMAX,DELTTH)
32       
33Cf2py intent(in) LAMB
34Cf2py intent(in) TTHMIN
35Cf2py intent(in) TTHMAX
36Cf2py intent(in) DELTTH
37   
38      INCLUDE 'DIFFaXsubs/DIFFaX.par'
39      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
40
41      REAL*8 LAMB,TTHMIN,TTHMAX,DELTTH
42       
43      lambda = lamb
44      th2_min = TTHMIN*DEG2RAD
45      th2_max = TTHMAX*DEG2RAD
46      d_theta = DELTTH*HALF*DEG2RAD
47       
48      RETURN
49      END
50       
51      SUBROUTINE PYGETCLAY(CNTRLS,LAUESYM,WDTH,NST,STSEQ)
52       
53Cf2py intent(in) CNTRLS
54Cf2py intent(in) LAUESYM
55Cf2py intent(in) WDTH
56Cf2py intent(in) NST
57Cf2py intent(in) STSEQ
58cf2py depend(NST) STSEQ
59     
60      INCLUDE 'DIFFaXsubs/DIFFaX.par'
61      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
62
63      CHARACTER*12 LAUESYM
64      INTEGER*4 CNTRLS(7),NST,STSEQ(NST),I
65      REAL*8 WDTH(2)                 
66      LOGICAL*4 ok,GETLAY
67      EXTERNAL GETLAY
68                                     
69      PI = FOUR * atan(ONE)
70      PI2 = TWO * PI
71      DEG2RAD = PI / ONE_EIGHTY
72      RAD2DEG = ONE_EIGHTY / PI
73      rad_type = X_RAY
74      lambda = 0.1
75      trim_origin = .TRUE.
76      blurring = NONE
77      loglin = 1
78      tolerance = 0.01
79      finite_width = .TRUE.
80      Wa = WDTH(1)*10000.
81      Wb = WDTH(2)*10000.
82      IF (Wa.GE.10000.) finite_width = .FALSE.
83      WRITE(pnt_grp,'(A12)') LAUESYM
84      SymGrpNo = CNTRLS(1)
85      check_sym = .TRUE.
86      full_shrp = 1
87      full_brd = 1
88C CNTRLS = [laueId,planeId,lmax,mult,StkType,StkParm,ranSeed]
89      bitdepth = 16
90      ok = .TRUE.
91      scaleint = FLOAT(CNTRLS(4))
92C fill in stacking seq stuff                 
93      IF (CNTRLS(5).NE.0) THEN
94        xplcit = .TRUE.
95        recrsv = .FALSE.
96        IF (CNTRLS(6).EQ.1) THEN
97            rndm = .TRUE.
98        ELSE
99            rndm = .FALSE.
100            l_cnt = NST
101            DO I=1,NST
102              l_seq(I) = STSEQ(I)
103            END DO       
104        END IF   
105      ELSE
106        recrsv = .TRUE.
107        xplcit = .FALSE.
108        IF (CNTRLS(6).NE.0) THEN
109            l_cnt = CNTRLS(6)
110            inf_thick = .FALSE.
111        ELSE
112            inf_thick = .TRUE.
113        END IF
114      END IF
115      IF (rndm) ok = GETLAY()
116      RETURN
117      END
118           
119      SUBROUTINE PYCELLAYER(CELL,NATM,ATMTP,ATMXOU,NU,LSYM,NL,LNUM)
120                   
121Cf2py intent(in) CELL
122Cf2py intent(in) NATM
123Cf2py intent(in) ATMTP
124Cf2py intent(in) ATMXOU
125cf2py depend(NATM) ATMTP,ATMXOU
126Cf2py intent(in) NU
127Cf2py intent(in) LSYM
128Cf2py depend(NU) LSYM
129Cf2py intent(in) NL
130Cf2py intent(in) LNUM
131Cf2py depend(NL) LNUM
132                       
133      INCLUDE 'DIFFaXsubs/DIFFaX.par'
134      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
135
136      INTEGER*4 NATM,NL,LNUM(NL),NU,LSYM(NU)
137      CHARACTER*4 ATMTP(NATM)
138      REAL*8  CELL(4),ATMXOU(8,NATM)
139      INTEGER*4 I,J,K,IL,IA
140
141C fill Common - cell stuff & finish symmetry stuff
142      cell_a = CELL(1)
143      cell_b = CELL(2)
144      cell_c = CELL(3)
145      cell_gamma = CELL(4)*DEG2RAD
146C fill common layer stuff - atoms & symm
147      DO I=1,NATM
148        IL = NINT(ATMXOU(1,I))
149        IA = NINT(ATMXOU(2,I))       
150        a_type(IA,IL) = NINT(ATMXOU(3,I))
151        a_number(IA,IL) = IA
152        WRITE(a_name(IA,IL),'(A4)') ATMTP(I)
153        DO K=1,3
154            a_pos(K,IA,IL) = ATMXOU(K+3,I)
155        END DO
156        high_atom(IL) = max(high_atom(IL),a_pos(3,IA,IL))
157        low_atom(IL) = min(low_atom(IL),a_pos(3,IA,IL))
158        IF (LSYM(IL).EQ.CENTRO) THEN
159            high_atom(IL) = MAX(high_atom(IL),-a_pos(3,IA,IL))
160            low_atom(IL) = MIN(low_atom(IL),-a_pos(3,IA,IL))
161        END IF
162        a_occup(IA,IL) = ATMXOU(7,I)
163        a_B(IA,IL) = ATMXOU(8,I)
164        l_n_atoms(IL) = IA
165        l_symmetry(IL) = LSYM(IL)
166      END DO
167      n_actual = IL
168      n_layers = NL
169      DO I=1,NL
170        l_actual(I) = LNUM(I)
171        DO J=1,NL
172            Bs_zero(J,I) = .TRUE.
173        END DO
174      END DO
175      all_Bs_zero = .TRUE.
176      RETURN
177      END
178
179      SUBROUTINE PYGETTRANS(NL,TRP,TRX)
180     
181Cf2py intent(in) NL
182Cf2py intent(in) TRP
183Cf2py intent(in) TRX
184Cf2py depend(NL) TRP,TRX
185     
186     
187      INCLUDE 'DIFFaXsubs/DIFFaX.par'
188      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
189       
190      INTEGER*4 I,J,K
191      INTEGER*4 NL
192      REAL*4  TRP(NL,NL),TRX(NL,NL,3)
193                               
194C fill common transitions stuff
195      DO J=1,NL
196        DO I=1,NL
197          l_alpha(J,I) = TRP(I,J)
198          DO K=1,3
199            l_r(K,J,I) = TRX(J,I,K)
200          END DO
201        END DO
202      END DO
203      RETURN
204      END
205       
206      SUBROUTINE PYGETSPC(CNTRLS,NSADP,SADP)
207       
208Cf2py intent(in) CNTRLS
209Cf2py intent(in) NSADP
210Cf2py intent(in/out) SADP
211Cf2py depend(NSADP) SADP
212           
213      INCLUDE 'DIFFaXsubs/DIFFaX.par'
214      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
215
216      INTEGER*4 CNTRLS(7),NSADP,I,j,k
217      REAL*8 SADP(NSADP),AGLQ16
218      LOGICAL GETSPC,ok
219       
220      EXTERNAL AGLQ16,GETSPC
221
222      DoSymDump = .FALSE.
223   
224      ok = .TRUE.
225      CALL SPHCST()
226      CALL DETUN()
227      CALL OPTIMZ('GSAS-II',ok)
228      If (debug) then
229        print *,cell_a,cell_b,cell_c,cell_gamma,pnt_grp,SymGrpNo
230        DoSymDump = .TRUE.
231        print *,n_actual,(l_n_atoms(i),i=1,n_actual)
232        do j=1,n_actual
233          do i=1,l_n_atoms(j)
234            print *,a_name(i,j),(a_pos(k,i,j),k=1,3)
235          end do
236        end do
237        do i=1,n_layers
238        print *,' layer',i
239           do j=1,n_layers
240              print *,'layer',j,l_alpha(i,j),(l_r(k,i,j),k=1,3)
241           end do
242        end do
243        print *, recrsv,inf_thick,xplcit,rndm,l_cnt,has_l_mirror
244       
245C      print *,lambda,max_angle,h_bnd,k_bnd,l_bnd,no_trials,
246C     1  rad_type,X_RAY,n_atoms
247C      print *,(l_g(j),j=1,n_layers)
248C      do j=1,n_layers
249C        print *,(hx_ky(i,j),i=1,l_n_atoms(j))
250C        print *,(mat(i,j),i=1,n_layers)
251C        print *,(mat1(i,j),i=1,n_layers)
252C        print *,(l_phi(i,j),i=1,n_layers)
253C      end do
254      end if
255       
256      ok = GETSPC(AGLQ16,'GSAS-II')             
257      DO I=1,NSADP
258        SADP(I) = spec(I)
259      END DO
260      RETURN
261      END
262       
263      SUBROUTINE PYPROFILE(U,V,W,HW,BLUR,NBRD,BRDSPC)
264       
265Cf2py intent(in) U
266Cf2py intent(in) V
267Cf2py intent(in) W
268Cf2py intent(in) HW
269Cf2py intent(in) NBRD
270Cf2py intent(in/out) BRDSPC
271Cf2py depend(NBRD) BRDSPC
272               
273      INCLUDE 'DIFFaXsubs/DIFFaX.par'
274      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
275       
276      INTEGER*4 BLUR,i,NBRD   
277      REAL*8 U,V,W,HW,BRDSPC(NBRD),tth_min
278       
279      tth_min = ZERO
280           
281      if (blur.eq.GAUSS) then
282        FWHM = HW
283        call GAUSSN(tth_min)
284      else if (blur.eq.PV_GSS) then
285        pv_u = U
286        pv_v = V
287        pv_w = W
288        pv_gamma = ZERO
289        call PV(tth_min)
290      end if
291      do i=1,NBRD
292        BRDSPC(i) = brd_spc(i)
293      end do
294       
295      RETURN
296      END
297       
298      SUBROUTINE PYGETSADP(CNTRLS,NSADP,SADP,HKLIM,INCR,NBLK)
299       
300Cf2py intent(in) CNTRLS
301Cf2py intent(in) NSADP
302Cf2py intent(in/out) SADP
303Cf2py depend(NSADP) SADP
304Cf2py intent(out) HKLIM
305Cf2py intent(out) INCR
306Cf2py intent(out) NBLK
307   
308      INCLUDE 'DIFFaXsubs/DIFFaX.par'
309      INCLUDE 'DIFFaXsubs/DIFFaX.inc'
310
311      INTEGER*4 CNTRLS(7),NSADP,i_plane,hk_lim,i,j,k
312      INTEGER*4 HKLIM,NBLK
313      REAL*8 SADP(NSADP),AGLQ16,l_upper,INCR
314      LOGICAL ok
315       
316      EXTERNAL AGLQ16                 
317                   
318      i_plane = CNTRLS(2)
319      l_upper = CNTRLS(3)
320      DoSymDump = .FALSE.
321      if (debug) then
322          print *,cell_a,cell_b,cell_c,cell_gamma
323          print *,pnt_grp,SymGrpNo
324          DoSymDump = .TRUE.
325          print *,n_actual,(l_n_atoms(i),i=1,n_actual)
326          do j=1,n_actual
327            do i=1,l_n_atoms(j)
328              print *,a_name(i,j),(a_pos(k,i,j),k=1,3)
329            end do
330          end do
331          do i=1,n_layers
332          print *,' layer',i
333             do j=1,n_layers
334                print *,'layer',j,l_alpha(i,j),(l_r(k,i,j),k=1,3)
335             end do
336          end do
337          print *, recrsv,inf_thick,xplcit,rndm,l_cnt,has_l_mirror
338      end if
339      ok = .TRUE.
340       
341      CALL SPHCST()
342      CALL DETUN()
343      CALL OPTIMZ('GSAS-II',ok)
344C      print *,lambda,max_angle,h_bnd,k_bnd,l_bnd,no_trials,
345C     1  rad_type,X_RAY,n_atoms
346C      print *,(l_g(j),j=1,n_layers)
347C      do j=1,n_layers
348C        print *,(hx_ky(i,j),i=1,l_n_atoms(j))
349C        print *,(mat(i,j),i=1,n_layers)
350C        print *,(mat1(i,j),i=1,n_layers)
351C        print *,(l_phi(i,j),i=1,n_layers)
352C      end do
353      CALL GETSAD(AGLQ16,i_plane,l_upper,hk_lim,'GSAS-II',ok)
354      NBLK = sadblock
355      HKLIM = hk_lim+1
356      INCR = dble(SADSIZE/2)/l_upper
357      if (i_plane.eq.1) then
358        INCR = INCR*sqrt(a0/c0)
359      else if (i_plane.eq.2) then
360        INCR = INCR*sqrt(b0/c0)
361      else if (i_plane.eq.3) then
362        INCR = INCR*sqrt((a0+b0+d0)/c0)
363      else if (i_plane.eq.4) then
364        INCR = INCR*sqrt((a0+b0-d0)/c0)
365      end if
366      do I=1,NSADP
367        SADP(i) = spec(i)
368      end do
369      RETURN
370      END
371
372           
Note: See TracBrowser for help on using the repository browser.