| @@ -0,0 +1,88 | |||
|
|
1 | subroutine guess(acf,tau,npts,zero,amin,te,tr) | |
|
|
2 | c | |
|
|
3 | c find zero crossing (zero), depth of minimum (amin), height of maximum | |
|
|
4 | c | |
|
|
5 | ||
|
|
6 | real acf(npts),tau(npts) | |
|
|
7 | c write(*,*) 'acf: ',acf | |
|
|
8 | c write(*,*) 'tau: ',tau | |
|
|
9 | c read(*,*) xx | |
|
|
10 | ||
|
|
11 | zero=0.0 | |
|
|
12 | amin=1.0 | |
|
|
13 | tmin=0.0 | |
|
|
14 | jmin=0 | |
|
|
15 | ||
|
|
16 | do i=npts,2,-1 | |
|
|
17 | if(acf(i)*acf(i-1).lt.0.0) then | |
|
|
18 | zero=(tau(i-1)*acf(i)-tau(i)*acf(i-1))/(acf(i)-acf(i-1)) | |
|
|
19 | end if | |
|
|
20 | if(acf(i).lt.amin) then | |
|
|
21 | amin=acf(i) | |
|
|
22 | jmin=i | |
|
|
23 | end if | |
|
|
24 | end do | |
|
|
25 | ||
|
|
26 | if(jmin.gt.0) then | |
|
|
27 | call parab1(tau(jmin-1),acf(jmin-1),a,b,c) | |
|
|
28 | tmin=-b/(2.0*a) | |
|
|
29 | amin=c+tmin*(b+tmin*a) | |
|
|
30 | end if | |
|
|
31 | ||
|
|
32 | tr=cdtr1(-amin) | |
|
|
33 | te=czte1(zero*1000.0,tr) | |
|
|
34 | c write(*,*) 'zero: ',zero | |
|
|
35 | c write(*,*) 'amin: ',amin | |
|
|
36 | c write(*,*) 'te: ',te | |
|
|
37 | c write(*,*) 'tr: ',tr | |
|
|
38 | c read(*,*) xx | |
|
|
39 | ||
|
|
40 | return | |
|
|
41 | end | |
|
|
42 | ||
|
|
43 | subroutine parab1(x,y,a,b,c) | |
|
|
44 | C----- | |
|
|
45 | dimension x(3),y(3) | |
|
|
46 | delta=x(1)-x(2) | |
|
|
47 | a=(y(1)-2.*y(2)+y(3))/(2.*delta*delta) | |
|
|
48 | b=(y(1)-y(2))/delta - a*(x(1)+x(2)) | |
|
|
49 | c=y(1)-a*x(1)*x(1)-b*x(1) | |
|
|
50 | return | |
|
|
51 | end | |
|
|
52 | ||
|
|
53 | real function cdtr1(depth) | |
|
|
54 | C-----convert depth to te/ti ratio | |
|
|
55 | dimension tr(4) | |
|
|
56 | C according to the curve published in farley et al 1967 | |
|
|
57 | c modified for 2004 conditions on axis | |
|
|
58 | c data tr/7.31081,3.53286,5.92271,.174/ | |
|
|
59 | data tr/9.5,4.0,8.5,.3/ | |
|
|
60 | data nt/4/ | |
|
|
61 | cdtr1=tr(1) | |
|
|
62 | do i=2,nt | |
|
|
63 | cdtr1=cdtr1*depth + tr(i) | |
|
|
64 | end do | |
|
|
65 | return | |
|
|
66 | end | |
|
|
67 | ||
|
|
68 | real function czte1(zlag,tr) | |
|
|
69 | C-----convert zero crossing point to te | |
|
|
70 | C according to the curve published in farley et al 1967 | |
|
|
71 | c modified for 2004 conditions on axis | |
|
|
72 | dimension dt(4) | |
|
|
73 | c data dt/0.00945025,-0.0774338,.203626,.812397/,nd/4/ | |
|
|
74 | data dt/0.00945025,-0.0774338,.2,0.9/,nd/4/ | |
|
|
75 | data t0/1000./ | |
|
|
76 | tr1=min(abs(tr),5.) | |
|
|
77 | if(zlag .eq. 0)then | |
|
|
78 | czte1=1000000. | |
|
|
79 | else | |
|
|
80 | dt0=dt(1) | |
|
|
81 | do i=2,nd | |
|
|
82 | dt0=dt0*tr1 + dt(i) | |
|
|
83 | end do | |
|
|
84 | czte1=t0*(dt0/zlag)**2 | |
|
|
85 | end if | |
|
|
86 | return | |
|
|
87 | end | |
|
|
88 | ||
| @@ -114,7 +114,8 setup(name='schainpy', | |||
|
|
114 | 114 | extra_f77_compile_args=["-fallow-argument-mismatch"]), |
|
|
115 | 115 | Extension("schainpy.model.proc.fitacf_guess", |
|
|
116 | 116 | sources=[ |
|
|
117 |
"schainf/Ffiles/fitacf_guess.pyf" |
|
|
|
117 | "schainf/Ffiles/fitacf_guess.pyf", | |
|
|
118 | "schainf/Ffiles/fitacf_guess.f",], | |
|
|
118 | 119 | extra_f77_compile_args=["-fallow-argument-mismatch"]), |
|
|
119 | 120 | Extension("schainpy.model.proc.fitacf_acf2", |
|
|
120 | 121 | sources = [ |
| @@ -227,6 +228,7 setup(name='schainpy', | |||
|
|
227 | 228 | "schainf/Ffiles/lmdif1.f", |
|
|
228 | 229 | "schainf/Ffiles/reader.c", |
|
|
229 | 230 | "schainf/Ffiles/cbesi.f", |
|
|
231 | "schainf/Ffiles/lagp.f", | |
|
|
230 | 232 | "schainf/Ffiles/i1mach.f", |
|
|
231 | 233 | "schainf/Ffiles/zeta.f", |
|
|
232 | 234 | "schainf/Ffiles/qc25f.f", |
General Comments 0
You need to be logged in to leave comments.
Login now