Commit 958b442e authored by ronesy's avatar ronesy
Browse files

Added global variable maxobs for default array size, implemented option zref...

Added global variable maxobs for default array size, implemented option zref and continuation to next file if current file not found
parent e0ad8a27
......@@ -30,6 +30,7 @@ module mod_var
integer, parameter :: max_path_len=200
integer, parameter :: max_name_len=100
integer, parameter :: maxobs=5000
real(kind=8) :: jdatei, jdatef
real(kind=8) :: jreldatei, jreldatef
integer :: nrec
......
......@@ -57,12 +57,12 @@ subroutine process_obs(settings, nobs, freq, jdobs, latobs, lonobs, altobs, conc
implicit none
type (settings_t), intent(in) :: settings
integer, intent(in out) :: nobs
real, intent(in) :: freq
real(kind=8), dimension(1000), intent(in) :: jdobs
real, dimension(1000), intent(in) :: latobs, lonobs, altobs, concobs, errobs
type (obs_t), intent(in out) :: obs
type (settings_t), intent(in) :: settings
integer, intent(in out) :: nobs
real, intent(in) :: freq
real(kind=8), dimension(maxobs), intent(in) :: jdobs
real, dimension(maxobs), intent(in) :: latobs, lonobs, altobs, concobs, errobs
type (obs_t), intent(in out) :: obs
integer :: cnt
integer :: c, i, m, n
......
......@@ -76,8 +76,8 @@ subroutine read_basic(settings, jd, nr, nobs, obs)
real :: dt, tmp, freq
real :: conc, err
real(kind=8) :: jdate
real(kind=8), dimension(1000) :: jdobs
real, dimension(1000) :: latobs, lonobs, altobs, concobs, errobs
real(kind=8), dimension(maxobs) :: jdobs
real, dimension(maxobs) :: latobs, lonobs, altobs, concobs, errobs
! freq determines length of the release period in the flexpart file RELEASES
! only used if not averaging
......@@ -87,12 +87,17 @@ subroutine read_basic(settings, jd, nr, nobs, obs)
print*, filelist
! find file matching receptor name
do n = 1, nfiles
do i = 1, len_trim(filelist(n))-2
if ( filelist(n)(i:(i+2)).eq.to_lower(trim(recname(nr))) ) go to 10
end do
end do
10 continue
if ( n.gt.nfiles ) then
print*, 'WARNING: file not found for receptor '//trim(recname(nr))
go to 20
endif
file_obs = filelist(n)
! open obs file
......
......@@ -64,7 +64,7 @@ subroutine read_noaa(settings, jd, nr, nobs, obs)
character(len=200) :: line
character(len=200), dimension(50) :: args
character(len=3) :: flag
real, dimension(12) :: temp
real, dimension(13) :: temp
logical :: lexist
integer :: ierr
integer :: cnt
......@@ -76,9 +76,9 @@ subroutine read_noaa(settings, jd, nr, nobs, obs)
real :: dt, tmp, freq
real :: conc, err
real(kind=8) :: jdate
real(kind=8), dimension(1000) :: jdobs
real, dimension(1000) :: latobs, lonobs, altobs, concobs, errobs
real(kind=8), dimension(1000,6) :: dataout
real(kind=8), dimension(maxobs) :: jdobs
real, dimension(maxobs) :: latobs, lonobs, altobs, concobs, errobs
real(kind=8), dimension(maxobs,6) :: dataout
! flask measurements assume cover 15 mins
freq = 0.25
......@@ -87,12 +87,17 @@ subroutine read_noaa(settings, jd, nr, nobs, obs)
print*, filelist
! find file matching receptor name
do n = 1, nfiles
do i = 1, len_trim(filelist(n))-2
if ( filelist(n)(i:(i+2)).eq.to_lower(trim(recname(nr))) ) go to 10
end do
end do
10 continue
if ( n.gt.nfiles ) then
print*, 'WARNING: file not found for receptor '//trim(recname(nr))
go to 20
endif
file_obs = filelist(n)
! open obs file
......@@ -132,6 +137,10 @@ subroutine read_noaa(settings, jd, nr, nobs, obs)
i = i + 1
read(args(n),*) temp(i)
end do
n = 27
i = i + 2
read(args(n),*,iostat=ierr) temp(i)
if ( ierr.ne.0 ) print*, 'WARNING: file does not appear to contain sample height -> using altitude'
yyyy = int(temp(1))
mm = int(temp(2))
dd = int(temp(3))
......@@ -145,7 +154,8 @@ subroutine read_noaa(settings, jd, nr, nobs, obs)
jjjjmmdd = yyyy*10000+mm*100+dd
hhmiss = hh*10000+mi*100+ss
jdate = juldate(jjjjmmdd, hhmiss)
if ( flag.ne.'...' ) cycle read_loop
! note last column of flag not used
if ( flag(1:2).ne.'..' ) cycle read_loop
if ( conc.le.-999. ) cycle read_loop
if ( (jdate.lt.jd).or.(jdate.ge.min(jreldatef,(jd+real(eomday,kind=8)))) ) cycle read_loop
cnt = cnt + 1
......@@ -154,7 +164,13 @@ subroutine read_noaa(settings, jd, nr, nobs, obs)
errobs(cnt) = err
latobs(cnt) = temp(9)
lonobs(cnt) = temp(10)
altobs(cnt) = temp(11)
if ( settings%zref.eq.1.and.ierr.eq.0 ) then
! height metres above ground
altobs(cnt) = temp(13)
else
! height metres above sea level
altobs(cnt) = temp(11)
endif
if ( errobs(cnt).le.-9.99 ) errobs(cnt) = 0.
end do read_loop
nobs = cnt
......
......@@ -74,18 +74,23 @@ subroutine read_obspack(settings, jd, nr, nobs, obs)
integer :: eomday, hloc
real :: conc, err, freq
real(kind=8) :: jdate
real(kind=8), dimension(1000) :: jdobs
real, dimension(1000) :: latobs, lonobs, altobs, concobs, errobs
real(kind=8), dimension(maxobs) :: jdobs
real, dimension(maxobs) :: latobs, lonobs, altobs, concobs, errobs
call caldate(jd, jjjjmmdd, hhmiss)
eomday = calceomday(jjjjmmdd/100)
! find file matching receptor name
do n = 1, nfiles
do i = 1, len_trim(filelist(n))-2
if ( filelist(n)(i:(i+2)).eq.to_lower(trim(recname(nr))) ) go to 10
end do
end do
10 continue
if ( n.gt.nfiles ) then
print*, 'WARNING: file not found for receptor '//trim(recname(nr))
go to 20
endif
file_obs = filelist(n)
! open obs file
......@@ -151,11 +156,23 @@ subroutine read_obspack(settings, jd, nr, nobs, obs)
if ( version.eq.'v2.1' ) then
latobs(cnt) = temp(10)
lonobs(cnt) = temp(11)
altobs(cnt) = temp(12)
if ( settings%zref.eq.1 ) then
! height metres above ground
altobs(cnt) = temp(14)
else if ( settings%zref.eq.2) then
! height metres above sea level
altobs(cnt) = temp(12)
endif
else if ( version.eq.'v3.2' ) then
latobs(cnt) = temp(12)
lonobs(cnt) = temp(13)
altobs(cnt) = temp(14)
if ( settings%zref.eq.1 ) then
! height metres above ground
altobs(cnt) = temp(16)
else if ( settings%zref.eq.2) then
! height metres above sea level
altobs(cnt) = temp(14)
endif
else
print*, 'ERROR: unknown ObsPack version'
stop
......@@ -167,6 +184,8 @@ subroutine read_obspack(settings, jd, nr, nobs, obs)
end do read_loop
nobs = cnt
print*, 'lat, lon, alt = ',latobs(1),lonobs(1),altobs(1)
! close input file
close(100)
......
......@@ -71,13 +71,13 @@ subroutine read_wdcgg(settings, jd, nr, nobs, obs)
integer :: cnt
integer :: narg, nflag
integer :: cs, rem
integer :: c, i, j, m, n, p, q
integer :: c, i, j, m, n, p, q, r
integer :: jjjjmmdd, hhmiss, yyyy, mm, dd, hh, mi, ss
integer :: eomday, hloc, num
real :: conc, sd, err, lon, lat, alt, freq
real(kind=8) :: jdate, jdateutc
real(kind=8), dimension(1000) :: jdobs
real, dimension(1000) :: latobs, lonobs, altobs, concobs, errobs
real(kind=8), dimension(maxobs) :: jdobs
real, dimension(maxobs) :: latobs, lonobs, altobs, concobs, errobs
call caldate(jd, jjjjmmdd, hhmiss)
eomday = calceomday(jjjjmmdd/100)
......@@ -85,6 +85,7 @@ subroutine read_wdcgg(settings, jd, nr, nobs, obs)
freq = 1.
! find file matching receptor name
do n = 1, nfiles
do i = 1, len_trim(filelist(n))-2
if ( filelist(n)(i:(i+2)).eq.to_lower(trim(recname(nr))) ) then
......@@ -95,6 +96,10 @@ subroutine read_wdcgg(settings, jd, nr, nobs, obs)
end do
end do
10 continue
if ( n.gt.nfiles ) then
print*, 'WARNING: file not found for receptor '//trim(recname(nr))
go to 20
endif
file_obs = filelist(n)
! open obs file
......@@ -108,7 +113,8 @@ subroutine read_wdcgg(settings, jd, nr, nobs, obs)
! read header
n = len_trim('LATITUDE:')
m = len_trim('LONGITUDE:')
p = len_trim('ALTITUDE:')
p = len_trim('SAMPLING HEIGHTS:')
r = len_trim('ALTITUDE:')
q = len_trim('DATE TIME')
do while (ierr.eq.0)
read (100, fmt='(A)', iostat=ierr) line
......@@ -121,9 +127,18 @@ subroutine read_wdcgg(settings, jd, nr, nobs, obs)
call parse_string(line, ":", args(:), narg)
read(args(2),*) lon
endif
if ( line(5:(p+4)) == 'ALTITUDE:' ) then
call parse_string(line, ":", args(:), narg)
read(args(2),*) alt
if ( settings%zref.eq.1 ) then
! height metres above ground
if ( line(5:(p+4)) == 'SAMPLING HEIGHTS:' ) then
call parse_string(line, ":", args(:), narg)
read(args(2),*) alt
endif
else if ( settings%zref.eq.2) then
! height metres above sea level
if ( line(5:(r+4)) == 'ALTITUDE:' ) then
call parse_string(line, ":", args(:), narg)
read(args(2),*) alt
endif
endif
end do
......@@ -140,7 +155,8 @@ subroutine read_wdcgg(settings, jd, nr, nobs, obs)
read(args(1:3),*) yyyy, mm, dd
call parse_string(atime, ":", args(:), narg)
read(args(1:2),*) hh, mi
read(flag,*) nflag
! read(flag,*) nflag
nflag = 0.
! calculate julian date
jjjjmmdd = yyyy*10000+mm*100+dd
hhmiss = hh*10000+mi*100
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment