123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276 |
- module read_input_file
- !-----------------------------------------------------------------------
- !
- ! contains routines for reading a file to interpolate
- !
- !-----------------------------------------------------------------------
- use kinds_mod ! defines common data types
- use constants ! defines useful constants
- use grids ! includes all grid information
- use netcdf_mod ! module with netcdf vars and utilities
- use remap_vars ! module for all required remapping variables
- implicit none
- real (kind=dbl_kind), dimension(:,:), allocatable ::
- & grid1_array
- integer (kind=int_kind) ::
- & ntime, nlat, nlon, natts, nvaratts, nglobatts
- integer (kind=int_kind) ::
- & nc_infile_id, nc_time_id, nc_invar_id, dimtime
- real (kind=dbl_kind), dimension(:), allocatable ::
- & time
- logical :: invertlat
- contains
- !***********************************************************************
- subroutine read_input(infile, varname)
- !-----------------------------------------------------------------------
- !
- ! input variables
- !
- !-----------------------------------------------------------------------
- character(char_len), intent(in) ::
- & infile, varname ! input filename an varname
- !-----------------------------------------------------------------------
-
- character(char_len) ::
- & timename, name2
- real (kind=dbl_kind), dimension(:,:,:), allocatable ::
- & grid1_array_tmp
- integer, dimension(:,:,:), allocatable ::
- & grid1_array_tmp1
-
- real, dimension(:,:,:), allocatable ::
- & grid1_array_tmp2
- double precision, dimension(:,:,:), allocatable ::
- & grid1_array_tmp3
- integer, dimension(:), allocatable ::
- & time1
- real, dimension(:), allocatable ::
- & time2
- double precision, dimension(:), allocatable ::
- & time3
- integer (kind=int_kind) ::
- & ncstat, nc_var_type, ndims, ji
- integer (kind=int_kind), dimension(:), allocatable :: ! netCDF ids
- & nc_dims_ids
- !-----------------------------------------------------------------------
- !
- ! read dimension information
- !
- !-----------------------------------------------------------------------
- ncstat = nf_open(infile, NF_NOWRITE, nc_infile_id)
- call netcdf_error_handler(ncstat)
- ncstat = nf_inq_varid(nc_infile_id, varname, nc_invar_id)
- call netcdf_error_handler(ncstat)
-
- ncstat = nf_inq_varndims(nc_infile_id, nc_invar_id, ndims)
- call netcdf_error_handler(ncstat)
- allocate(nc_dims_ids(ndims))
- ncstat = nf_inq_vardimid(nc_infile_id, nc_invar_id, nc_dims_ids)
- call netcdf_error_handler(ncstat)
- if ( ndims < 1 .or. ndims > 3) then
- stop "Input files should have (lon) or (lon,lat) or (lon,time)
- & or (lon,lat,time) dimensions"
- endif
- ncstat = nf_inq_dimlen(nc_infile_id, nc_dims_ids(1), nlon)
- call netcdf_error_handler(ncstat)
- if (nlon .ne. grid1_dims(1)) then
- print*," Input file : ",nlon," longitudes"
- print*," Weight file : ",grid1_dims(1) ," longitudes"
- stop " Inconsistent number of longitudes. We stop "
- endif
- name2='none'
- nlat=1
- if ( ndims >= 2 ) then
- ncstat = nf_inq_dimname(nc_infile_id, nc_dims_ids(2), name2)
-
- if ( name2 /= 'time' .and. name2 /= 'time_counter' .and.
- & name2 /= 't' ) then
-
- ncstat = nf_inq_dimlen(nc_infile_id, nc_dims_ids(2), nlat)
- call netcdf_error_handler(ncstat)
- if (nlat .ne. grid1_dims(2)) then
- print*," Input file : ",nlat," latitudes"
- print*," Weight file : ",grid1_dims(2)," latitudes"
- stop " Inconsistent number of latitudes. We stop "
- endif
- dimtime=3
- else
- dimtime=2
- endif
- endif
- !-----------------------------------------------------------------------
- !
- ! read time dimension
- !
- !-----------------------------------------------------------------------
- nc_time_id=-1
- if ( ndims == 3 .or. dimtime == 2 ) then
- ncstat = nf_inq_dim(nc_infile_id, nc_dims_ids(dimtime),
- & timename, ntime)
- call netcdf_error_handler(ncstat)
-
- allocate(time(ntime))
- ncstat = nf_inq_varid(nc_infile_id, timename, nc_time_id)
- if (ncstat == NF_NOERR) then
- ncstat = nf_inq_vartype(nc_infile_id, nc_time_id, nc_var_type)
- call netcdf_error_handler(ncstat)
-
- allocate(time1(ntime))
- allocate(time2(ntime))
- allocate(time3(ntime))
- if ( nc_var_type == NF_INT ) then
- ncstat = nf_get_var_int(nc_infile_id, nc_time_id, time1)
- if ( ncstat == NF_NOERR ) then
- time = real(time1,kind=dbl_kind)
- else
- stop "Problem reading input data"
- endif
- elseif (nc_var_type == NF_REAL) then
- ncstat = nf_get_var_real(nc_infile_id, nc_time_id, time2)
- if ( ncstat == NF_NOERR ) then
- time = real(time2,kind=dbl_kind)
- else
- stop "Problem reading input data"
- endif
- elseif ( nc_var_type == NF_DOUBLE) then
- ncstat = nf_get_var_double(nc_infile_id,nc_time_id, time3)
- if ( ncstat == NF_NOERR ) then
- time = real(time3,kind=dbl_kind)
- else
- stop"Problem reading input data"
- endif
- else
- stop"Problem with input data type"
- endif
- deallocate(time1, time2, time3)
- ncstat = nf_inq_varnatts(nc_infile_id, nc_time_id, natts)
- call netcdf_error_handler(ncstat)
- else
- do ji=1,ntime
- time(ji)=ji
- enddo
- natts=0
- endif
- else
- ntime=1
- allocate(time(1))
- time=1.
- natts=0
- endif
- !-----------------------------------------------------------------------
- !
- ! allocate arrays
- !
- !-----------------------------------------------------------------------
- allocate( grid1_array_tmp (nlon,nlat,ntime))
- allocate( grid1_array_tmp1 (nlon,nlat,ntime))
- allocate( grid1_array_tmp2 (nlon,nlat,ntime))
- allocate( grid1_array_tmp3 (nlon,nlat,ntime))
- !-----------------------------------------------------------------------
- !
- ! read variable
- !
- !-----------------------------------------------------------------------
-
- ncstat = nf_inq_vartype(nc_infile_id, nc_invar_id, nc_var_type)
- call netcdf_error_handler(ncstat)
- if ( nc_var_type == NF_INT ) then
- ncstat = nf_get_var_int(nc_infile_id, nc_invar_id,
- & grid1_array_tmp1)
- if ( ncstat == NF_NOERR ) then
- grid1_array_tmp = real(grid1_array_tmp1,kind=dbl_kind)
- else
- stop "Problem reading input data"
- endif
- elseif (nc_var_type == NF_REAL) then
- ncstat = nf_get_var_real(nc_infile_id, nc_invar_id,
- & grid1_array_tmp2)
- if ( ncstat == NF_NOERR ) then
- grid1_array_tmp = real(grid1_array_tmp2,kind=dbl_kind)
- else
- stop "Problem reading input data"
- endif
- elseif ( nc_var_type == NF_DOUBLE) then
- ncstat = nf_get_var_double(nc_infile_id,nc_invar_id,
- & grid1_array_tmp3)
- if ( ncstat == NF_NOERR ) then
- grid1_array_tmp = real(grid1_array_tmp3,kind=dbl_kind)
- else
- stop"Problem reading input data"
- endif
- else
- stop"Problem with input data type"
- endif
- ncstat = nf_inq_varnatts(nc_infile_id, nc_invar_id, nvaratts)
- call netcdf_error_handler(ncstat)
-
- deallocate(grid1_array_tmp1, grid1_array_tmp2, grid1_array_tmp3)
- !-----------------------------------------------------------------------
- !
- ! reshape input file
- !
- !-----------------------------------------------------------------------
-
- allocate( grid1_array (nlon*nlat,ntime))
- if (invertlat) then
- grid1_array = reshape ( grid1_array_tmp (:,nlat:1:-1,:),
- & (/ nlat*nlon, ntime /) )
- else
- grid1_array = reshape (grid1_array_tmp, (/ nlat*nlon, ntime /) )
- endif
- deallocate(grid1_array_tmp)
- ncstat = nf_inq_varnatts(nc_infile_id, NF_GLOBAL, nglobatts)
- call netcdf_error_handler(ncstat)
- !-----------------------------------------------------------------------
- end subroutine read_input
- !***********************************************************************
- end module read_input_file
- !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|