da_be4_scale_length.f90
References to this file elsewhere.
1 program da_scale_length
2
3 ! use da_control
4
5 implicit none
6
7 integer, parameter :: nt=1000 ! maximum time series
8 character*10 :: start_date, end_date ! Starting and ending dates.
9 character*10 :: date, new_date ! Current date (ccyymmddhh).
10 character*10 :: variable ! Variable name
11 character*3 :: be_method ! Be method ('NMC', or 'ENS')
12 character*80 :: filename ! Input filename.
13 character*1 :: k_label ! = 'k' if model level, 'm' if mode output.
14 character*2 :: ck ! Level index -> character.
15 character*3 :: ce ! Member index -> character.
16 integer :: ni, nj, nk ! Dimensions read in.
17 integer :: ni0, nj0 ! Dimensions read in.
18 integer :: vertical_level
19 integer :: sdate, cdate, edate ! Starting, current ending dates.
20 integer :: interval ! Period between dates (hours).
21 integer :: ne ! Number of ensemble members.
22 integer :: i, j, k, k1, k2, b, member ! Loop counters.
23
24 real, allocatable :: field(:,:,:) ! Field projected into EOF space.
25
26 real :: cut_dist_km,resolution_km
27
28 namelist / scale_length_nl / start_date, end_date, interval, variable, &
29 be_method, ne, &
30 ni, nj, nk ,vertical_level, &
31 cut_dist_km, resolution_km
32 logical :: use_global_eofs, data_on_levels
33 integer :: num
34
35 write(6,'(a)')' [1] Initialize namelist variables and other scalars.'
36
37
38 ! vertical_ip = 0
39
40 start_date = '2004030312'
41 end_date = '2004033112'
42 interval = 24
43 variable = 'psi'
44 be_method = 'NMC'
45 ne = 1
46
47 use_global_eofs = .false.
48 data_on_levels = .false.
49
50 open(unit=namelist_unit, file='scale_length_nl.nl', &
51 form='formatted', status='old', action='read')
52 read(namelist_unit, scale_length_nl)
53 close(namelist_unit)
54
55 read(start_date(1:10), fmt='(i10)')sdate
56 read(end_date(1:10), fmt='(i10)')edate
57 write(6,'(4a)')' Computing vertical error statistics for dates ', start_date, ' to ', end_date
58 write(6,'(a,i8,a)')' Interval between dates = ', interval, 'hours.'
59 write(6,'(a,i8)')' Number of ensemble members at each time = ', ne
60
61 date = start_date
62 cdate = sdate
63
64 allocate(field(1:ni, 1:nj, 1:nt))
65
66 num=0
67 do while ( cdate <= edate )
68 do member = 1, ne
69 num=num+1
70 write(6,'(5a,i4)')' Date = ', date, ', variable ', trim(variable), &
71 ' and member ', member
72
73 write(ce,'(i3)')member
74 if ( member < 10 ) ce = '00'//ce(3:3)
75 if ( member >= 10 .and. member < 100 ) ce = '0'//ce(2:3)
76
77 ! Output fields (split into 2D files to allow parallel horizontal treatment):
78 ! read(ck,'(i2)') nk
79 write(ck,'(i2)') nk
80
81 if ( nk < 10 ) ck = '0'//ck(2:2)
82 k_label = 'm'
83
84 ! Assumes variable directory has been created by script:
85 filename = trim(variable)//'/'//date(1:10)//'.'//trim(variable)
86 ! ols.101204
87 ! if (trim(variable).eq.'ps_u') then
88 ! filename = trim(filename)//'.'//trim(be_method)//'.e'//ce
89 ! else
90 ! filename = trim(filename)//'.'//trim(be_method)//'.e'//ce//'.'//k_label//ck
91 ! end if
92 filename = trim(filename)//'.'//trim(be_method)//'.e'//ce//'.'//ck
93
94 open (iunit, file = filename, form='unformatted')
95 read(iunit) ni0, nj0, k
96 if (ni.ne.ni0.or.nj.ne.nj0) then
97 write(*,'(a)') 'dimension problem',ni,ni0,nj,nj0
98 stop
99 end if
100 read(iunit) data_on_levels, use_global_eofs
101 read(iunit)field(1:ni,1:nj,num)
102 end do ! End loop over members.
103
104 ! Calculate next date:
105 call da_advance_cymdh( date, interval, new_date )
106 date = new_date
107 read(date(1:10), fmt='(i10)')cdate
108 end do
109
110 num=num-1
111
112 call da_process_single_variable(field,variable,cut_dist_km,&
113 resolution_km,num,ni,nj,ck,nt)
114
115 if (nk.eq.vertical_level) then
116 call da_merge_scale_length(variable, nk)
117 end if
118
119 end program da_scale_length