-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathpastesoilcode.f90
115 lines (74 loc) · 2.81 KB
/
pastesoilcode.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
program pastesoilcode
use iso_fortran_env
use netcdf
implicit none
integer, parameter :: i1 = int8
character(100) :: infile
character(100) :: outfile
integer :: status
integer :: ncid
integer :: dimid
integer :: varid
integer :: xlen
integer :: ylen
integer(i1) :: imissing
integer(i1), dimension(2) :: actual_range
integer(i1), allocatable, dimension(:,:) :: var
character(40) :: varname
!------------------------------------------------------------------------------------------------------------
call getarg(1,infile)
call getarg(2,outfile)
call getarg(3,varname)
!-------------------------------------------------------
! read input file
! write(0,*)'reading'
status = nf90_open(infile,nf90_nowrite,ncid)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_inq_dimid(ncid,'lon',dimid)
if (status == nf90_ebaddim) status = nf90_inq_dimid(ncid,'x',dimid)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_inquire_dimension(ncid,dimid,len=xlen)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_inq_dimid(ncid,'lat',dimid)
if (status == nf90_ebaddim) status = nf90_inq_dimid(ncid,'y',dimid)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_inquire_dimension(ncid,dimid,len=ylen)
if (status /= nf90_noerr) call handle_err(status)
allocate(var(xlen,ylen))
status = nf90_inq_varid(ncid,'Band1',varid)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_get_var(ncid,varid,var)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_get_att(ncid,varid,'_FillValue',imissing)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_close(ncid)
if (status /= nf90_noerr) call handle_err(status)
!-------------------------------------------------------
! write to output file
! write(0,*)'writing'
status = nf90_open(outfile,nf90_write,ncid)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_inq_varid(ncid,varname,varid)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_put_var(ncid,varid,var)
if (status /= nf90_noerr) call handle_err(status)
actual_range(1) = minval(var,mask=var/=imissing)
actual_range(2) = maxval(var,mask=var/=imissing)
status = nf90_put_att(ncid,varid,'actual_range',actual_range)
if (status /= nf90_noerr) call handle_err(status)
status = nf90_close(ncid)
if (status /= nf90_noerr) call handle_err(status)
write(0,*)trim(varname),' range:',actual_range
!-------------------------------------------------------
contains
subroutine handle_err(status)
! Internal subroutine - checks error status after each netcdf call,
! prints out text message each time an error code is returned.
integer, intent (in) :: status
if(status /= nf90_noerr) then
print *, trim(nf90_strerror(status))
stop
end if
end subroutine handle_err
!-------------------------------------------------------
end program pastesoilcode