field_routines_pnc.F90
References to this file elsewhere.
1 !*------------------------------------------------------------------------------
2 !* Standard Disclaimer
3 !*
4 !* Forecast Systems Laboratory
5 !* NOAA/OAR/ERL/FSL
6 !* 325 Broadway
7 !* Boulder, CO 80303
8 !*
9 !* AVIATION DIVISION
10 !* ADVANCED COMPUTING BRANCH
11 !* SMS/NNT Version: 2.0.0
12 !*
13 !* This software and its documentation are in the public domain and
14 !* are furnished "as is". The United States government, its
15 !* instrumentalities, officers, employees, and agents make no
16 !* warranty, express or implied, as to the usefulness of the software
17 !* and documentation for any purpose. They assume no
18 !* responsibility (1) for the use of the software and documentation;
19 !* or (2) to provide technical support to users.
20 !*
21 !* Permission to use, copy, modify, and distribute this software is
22 !* hereby granted, provided that this disclaimer notice appears in
23 !* all copies. All modifications to this software must be clearly
24 !* documented, and are solely the responsibility of the agent making
25 !* the modification. If significant modifications or enhancements
26 !* are made to this software, the SMS Development team
27 !* (sms-info@fsl.noaa.gov) should be notified.
28 !*
29 !*----------------------------------------------------------------------------
30 !*
31 !* WRF NetCDF I/O
32 ! Author: Jacques Middlecoff jacquesm@fsl.noaa.gov
33 !* Date: October 6, 2000
34 !*
35 !*----------------------------------------------------------------------------
36 subroutine ext_pnc_RealFieldIO(IO,NCID,VarID,VStart,VCount,Data,Status)
37 use wrf_data
38 use ext_pnc_support_routines
39 implicit none
40 include 'wrf_status_codes.h'
41 # include "pnetcdf.inc"
42 character (*) ,intent(in) :: IO
43 integer ,intent(in) :: NCID
44 integer ,intent(in) :: VarID
45 integer ,dimension(NVarDims),intent(in) :: VStart
46 integer ,dimension(NVarDims),intent(in) :: VCount
47 real, dimension(*) ,intent(inout) :: Data
48 integer ,intent(out) :: Status
49 integer :: stat
50 !local
51 integer(KIND=MPI_OFFSET_KIND), dimension(NVarDims) :: VStart_mpi, VCount_mpi
52 VStart_mpi = VStart
53 VCount_mpi = VCount
54
55 if(IO == 'write') then
56 stat = NFMPI_PUT_VARA_REAL_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Data)
57 else
58 stat = NFMPI_GET_VARA_REAL_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Data)
59 endif
60 call netcdf_err(stat,Status)
61 if(Status /= WRF_NO_ERR) then
62 write(msg,*) 'NetCDF error in ',__FILE__,', line', __LINE__
63 call wrf_debug ( WARN , msg)
64 endif
65 return
66 end subroutine ext_pnc_RealFieldIO
67
68 subroutine ext_pnc_DoubleFieldIO(IO,NCID,VarID,VStart,VCount,Data,Status)
69 use wrf_data
70 use ext_pnc_support_routines
71 implicit none
72 include 'wrf_status_codes.h'
73 # include "pnetcdf.inc"
74 character (*) ,intent(in) :: IO
75 integer ,intent(in) :: NCID
76 integer ,intent(in) :: VarID
77 integer ,dimension(NVarDims),intent(in) :: VStart
78 integer ,dimension(NVarDims),intent(in) :: VCount
79 real*8 ,intent(inout) :: Data
80 integer ,intent(out) :: Status
81 integer :: stat
82 !local
83 integer(KIND=MPI_OFFSET_KIND), dimension(NVarDims) :: VStart_mpi, VCount_mpi
84 VStart_mpi = VStart
85 VCount_mpi = VCount
86
87 if(IO == 'write') then
88 stat = NFMPI_PUT_VARA_DOUBLE_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Data)
89 else
90 stat = NFMPI_GET_VARA_DOUBLE_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Data)
91 endif
92 call netcdf_err(stat,Status)
93 if(Status /= WRF_NO_ERR) then
94 write(msg,*) 'NetCDF error in ',__FILE__,', line', __LINE__
95 call wrf_debug ( WARN , msg)
96 endif
97 return
98 end subroutine ext_pnc_DoubleFieldIO
99
100 subroutine ext_pnc_IntFieldIO(IO,NCID,VarID,VStart,VCount,Data,Status)
101 use wrf_data
102 use ext_pnc_support_routines
103 implicit none
104 include 'wrf_status_codes.h'
105 # include "pnetcdf.inc"
106 character (*) ,intent(in) :: IO
107 integer ,intent(in) :: NCID
108 integer ,intent(in) :: VarID
109 integer ,dimension(NVarDims),intent(in) :: VStart
110 integer ,dimension(NVarDims),intent(in) :: VCount
111 integer ,intent(inout) :: Data
112 integer ,intent(out) :: Status
113 integer :: stat
114 !local
115 integer(KIND=MPI_OFFSET_KIND), dimension(NVarDims) :: VStart_mpi, VCount_mpi
116 VStart_mpi = VStart
117 VCount_mpi = VCount
118
119 if(IO == 'write') then
120 stat = NFMPI_PUT_VARA_INT_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Data)
121 else
122 stat = NFMPI_GET_VARA_INT_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Data)
123 endif
124 call netcdf_err(stat,Status)
125 if(Status /= WRF_NO_ERR) then
126 write(msg,*) 'NetCDF error in ',__FILE__,', line', __LINE__
127 call wrf_debug ( WARN , msg)
128 endif
129 return
130 end subroutine ext_pnc_IntFieldIO
131
132 subroutine ext_pnc_LogicalFieldIO(IO,NCID,VarID,VStart,VCount,Data,Status)
133 use wrf_data
134 use ext_pnc_support_routines
135 implicit none
136 include 'wrf_status_codes.h'
137 # include "pnetcdf.inc"
138 character (*) ,intent(in) :: IO
139 integer ,intent(in) :: NCID
140 integer ,intent(in) :: VarID
141 integer,dimension(NVarDims) ,intent(in) :: VStart
142 integer,dimension(NVarDims) ,intent(in) :: VCount
143 logical,dimension(VCount(1),VCount(2),VCount(3)),intent(inout) :: Data
144 integer ,intent(out) :: Status
145 integer,dimension(:,:,:),allocatable :: Buffer
146 integer :: stat
147 integer :: i,j,k
148 !local
149 integer(KIND=MPI_OFFSET_KIND), dimension(NVarDims) :: VStart_mpi, VCount_mpi
150 VStart_mpi = VStart
151 VCount_mpi = VCount
152
153 allocate(Buffer(VCount(1),VCount(2),VCount(3)), STAT=stat)
154 if(stat/= 0) then
155 Status = WRF_ERR_FATAL_ALLOCATION_ERROR
156 write(msg,*) 'Fatal ALLOCATION ERROR in ',__FILE__,', line', __LINE__
157 call wrf_debug ( FATAL , msg)
158 return
159 endif
160 if(IO == 'write') then
161 do k=1,VCount(3)
162 do j=1,VCount(2)
163 do i=1,VCount(1)
164 if(data(i,j,k)) then
165 Buffer(i,j,k)=1
166 else
167 Buffer(i,j,k)=0
168 endif
169 enddo
170 enddo
171 enddo
172 stat = NFMPI_PUT_VARA_INT_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Buffer)
173 else
174 stat = NFMPI_GET_VARA_INT_ALL(NCID,VarID,VStart_mpi,VCount_mpi,Buffer)
175 Data = Buffer == 1
176 endif
177 call netcdf_err(stat,Status)
178 if(Status /= WRF_NO_ERR) then
179 write(msg,*) 'NetCDF error in ',__FILE__,', line', __LINE__
180 call wrf_debug ( WARN , msg)
181 return
182 endif
183 deallocate(Buffer, STAT=stat)
184 if(stat/= 0) then
185 Status = WRF_ERR_FATAL_DEALLOCATION_ERR
186 write(msg,*) 'Fatal DEALLOCATION ERROR in ',__FILE__,', line', __LINE__
187 call wrf_debug ( FATAL , msg)
188 return
189 endif
190 return
191 end subroutine ext_pnc_LogicalFieldIO