-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathNORTRIP_read_t2m500yr_netcdf4.f90
415 lines (337 loc) · 21.6 KB
/
NORTRIP_read_t2m500yr_netcdf4.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
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
subroutine NORTRIP_read_analysismeteo_netcdf4
!Reads yr temperature 66 hour forecast data
use NORTRIP_multiroad_index_definitions
!Update to netcdf 4 and 64 bit in this version 2 of NORTRIP_read_t2m500yr_netcdf
use netcdf
implicit none
!include 'netcdf.inc'
!Local variables
integer status_nc2 !Error message
integer id_nc2
integer dim_id_nc2(num_dims_nc2)
integer xtype_nc2(num_var_nc2)
integer natts_nc2(num_var_nc2)
integer var_id_nc2(num_var_nc2)
integer dim_length_metcoop_nc2(num_dims_nc2+1)
integer dim_start_metcoop_nc2(num_dims_nc2+1)
character(256) dimname_temp
integer i
integer i_grid_mid,j_grid_mid
real dlat_nc2
integer exists
integer ii,jj,tt,t
integer var_id_nc2_projection
character(256) pathname_nc2_in,filename_nc2_in
integer new_start_date_input(num_date_index)
double precision temp_date
double precision date_to_number
double precision, allocatable :: var2d_nc2_dp(:,:)
double precision, allocatable :: var3d_nc2_dp(:,:)
logical dim_read_flag
write(unit_logfile,'(A)') '================================================================'
write(unit_logfile,'(A)') 'Reading additional meteorological data (NORTRIP_read_analysismeteo_netcdf4)'
write(unit_logfile,'(A)') '================================================================'
pathfilename_nc2=trim(pathname_nc2)//trim(filename_nc2)
pathname_nc2_in=pathname_nc2
filename_nc2_in=filename_nc2_template
new_start_date_input=start_date_input
if (.not.allocated(meteo_nc2_available)) allocate (meteo_nc2_available(n_hours_input))
if (.not.allocated(meteo_var_nc2_available)) allocate (meteo_var_nc2_available(n_hours_input,num_var_nc2))
meteo_var_nc2_available=.true.
dim_read_flag=.false.
!Loop through the number of time steps nad read in data when available
do t=1,n_hours_input
temp_date=date_to_number(start_date_input,ref_year)
call number_to_date(temp_date+(t-1)/dble(24.),new_start_date_input,ref_year)
write(unit_logfile,'(a,7i)') 'Date array: ',t,new_start_date_input(1:6)
call date_to_datestr_bracket(new_start_date_input,filename_nc2_in,filename_nc2)
call date_to_datestr_bracket(new_start_date_input,pathname_nc2_in,pathname_nc2)
pathfilename_nc2=trim(pathname_nc2)//trim(filename_nc2)
!Test existence of the filename. If does not exist then skip the time index
inquire(file=trim(pathfilename_nc2),exist=exists)
if (.not.exists) then
write(unit_logfile,'(A,A)') ' WARNING: Meteo netcdf2 file does not exist: ', trim(pathfilename_nc2)
meteo_nc2_available(t)=.false.
else
meteo_nc2_available(t)=.true.
endif
!Open the netcdf file for reading
if (meteo_nc2_available(t)) then
write(unit_logfile,'(2A)') ' Opening netcdf meteo file: ',trim(pathfilename_nc2)
status_nc2 = NF90_OPEN (pathfilename_nc2, NF90_NOWRITE, id_nc2)
if (status_nc2 .NE. NF90_NOERR) write(unit_logfile,'(A,I)') 'ERROR opening netcdf file: ',status_nc2
!Find the projection. If no projection then in lat lon coordinates
status_nc2 = NF90_INQ_VARID (id_nc2,trim(projection_name_nc2),var_id_nc2_projection)
if (status_nc2.eq.NF90_NOERR) then
!If there is a projection then read in the attributes. All these are doubles
!status_nc = nf90_inquire_variable(id_nc, var_id_nc_projection, natts = numAtts_projection)
status_nc2 = nf90_get_att(id_nc2, var_id_nc2_projection, 'standard_parallel', meteo_nc2_projection_attributes(1:2))
status_nc2 = nf90_get_att(id_nc2, var_id_nc2_projection, 'longitude_of_central_meridian', meteo_nc2_projection_attributes(3))
status_nc2 = nf90_get_att(id_nc2, var_id_nc2_projection, 'latitude_of_projection_origin', meteo_nc2_projection_attributes(4))
status_nc2 = nf90_get_att(id_nc2, var_id_nc2_projection, 'earth_radius', meteo_nc2_projection_attributes(5))
meteo_nc2_projection_type=LCC_projection_index
write(unit_logfile,'(A,5f12.2)') 'Reading lambert_conformal_conic projection. ',meteo_nc2_projection_attributes(1:5)
else
meteo_nc2_projection_type=LL_projection_index
endif
if (.not.dim_read_flag) then
!Find out the x,y and time dimmensions of the file
status_nc2 = NF90_INQ_DIMID (id_nc2,dim_name_nc2(x_index2),dim_id_nc2(x_index2))
status_nc2 = NF90_INQUIRE_DIMENSION (id_nc2,dim_id_nc2(x_index2),dimname_temp,dim_length_nc2(x_index2))
status_nc2 = NF90_INQ_DIMID (id_nc2,dim_name_nc2(y_index2),dim_id_nc2(y_index2))
status_nc2 = NF90_INQUIRE_DIMENSION (id_nc2,dim_id_nc2(y_index2),dimname_temp,dim_length_nc2(y_index2))
call NORTRIP_reduce_meteo_region2(id_nc2)
dim_read_flag=.true.
endif
status_nc2 = NF90_INQ_DIMID (id_nc2,dim_name_nc2(time_index2),dim_id_nc2(time_index2))
status_nc2 = NF90_INQUIRE_DIMENSION (id_nc2,dim_id_nc2(time_index2),dimname_temp,dim_length_nc2(time_index2))
write(unit_logfile,'(A,3I)') ' Size of dimensions (x,y,t): ',dim_length_nc2
dim_length_nc2(time_index)=number_of_time_steps
write(unit_logfile,'(A,3I)') ' WARNING: Reducing dimensions of (t) to save space: ',dim_length_nc2(time_index)
!Allocate the nc arrays for reading
!write(*,*) dim_length_nc2(x_index2),dim_length_nc2(y_index2),dim_length_nc2(time_index2)
if (.not.allocated(var1d_nc2)) allocate (var1d_nc2(num_dims_nc2,maxval(dim_length_nc2))) !x and y and time maximum dimmensions
if (.not.allocated(var3d_nc2)) allocate (var3d_nc2(num_var_nc2,dim_length_nc2(x_index2),dim_length_nc2(y_index2),n_hours_input))
if (.not.allocated(var2d_nc2)) allocate (var2d_nc2(num_var_nc2,dim_length_nc2(x_index2),dim_length_nc2(y_index2))) !Lat and lon and elevation
if (.not.allocated(var3d_nc2_dp)) allocate (var3d_nc2_dp(dim_length_nc2(x_index2),dim_length_nc2(y_index2)))
if (.not.allocated(var2d_nc2_dp)) allocate (var2d_nc2_dp(dim_length_nc2(x_index2),dim_length_nc2(y_index2))) !Lat and lon
!Set the number of hours to be read
!Read the x, y and time values
do i=1,num_dims_nc2
status_nc2 = NF90_INQ_VARID (id_nc2, trim(dim_name_nc2(i)), var_id_nc2(i))
!status_nc2 = NF_GET_VARA_REAL (id_nc2, var_id_nc2(i), dim_start_nc2(i), dim_length_nc2(i), var1d_nc2(i,:))
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var1d_nc2(i,1:dim_length_nc2(i)), start=(/dim_start_nc2(i)/), count=(/dim_length_nc2(i)/))
if (i.eq.time_index2) then
write(unit_logfile,'(3A,2i12)') ' ',trim(dim_name_nc2(i)),' (min, max in hours): ' &
,minval(int((var1d_nc2(i,1:dim_length_nc2(i))-var1d_nc2(i,dim_start_nc2(i)))/3600.+.5)+1) &
,maxval(int((var1d_nc2(i,1:dim_length_nc2(i))-var1d_nc2(i,dim_start_nc2(i)))/3600.+.5)+1)
else
write(unit_logfile,'(3A,2f12.2)') ' ',trim(dim_name_nc2(i)),' (min, max in km): ' &
,minval(var1d_nc2(i,1:dim_length_nc2(i))),maxval(var1d_nc2(i,1:dim_length_nc2(i)))
endif
enddo
!MetCoOp data has 3 dimensions. z is 1 so must adjust this.
!dim_length_metcoop_nc(1:2)=dim_length_nc(1:2)
!dim_length_metcoop_nc(3)=1
!dim_length_metcoop_nc(4)=dim_length_nc(time_index)
!dim_start_metcoop_nc(1:2)=dim_start_nc(1:2)
!dim_start_metcoop_nc(3)=1
!dim_start_metcoop_nc(4)=dim_start_nc(time_index)
!Read through the variables in a loop
do i=1,num_var_nc2
!write(*,*) i,trim(var_name_nc(i))
status_nc2 = NF90_INQ_VARID (id_nc2, trim(var_name_nc2(i)), var_id_nc2(i))
!write(*,*) 'Status1: ',status_nc2,id_nc2,var_id_nc2(i),trim(var_name_nc2(i)),NF_NOERR
!write(*,*) 'Status1: ',dim_start_metcoop_nc
!write(*,*) 'Status1: ',dim_length_metcoop_nc
if (status_nc2.eq.NF90_NOERR) then
if (i.eq.lat_index2.or.i.eq.lon_index2) then
!write(*,*) id_nc2, var_id_nc2(i), dim_start_nc2(1:2), dim_length_nc2(1:2)
!status_nc = NF_GET_VARA_REAL (id_nc, var_id_nc(i), dim_start_metcoop_nc(1:2), dim_length_metcoop_nc(1:2), var2d_nc(i,:,:))
!status_nc2 = NF_GET_VARA_DOUBLE (id_nc2, var_id_nc2(i), dim_start_nc2(1:2), dim_length_nc2(1:2), var2d_nc2_dp);var2d_nc2(i,:,:)=real(var2d_nc2_dp)
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var2d_nc2_dp(:,:), start=(/dim_start_nc2/), count=(/dim_length_nc2/));var2d_nc2(i,:,:)=real(var2d_nc2_dp)
write(unit_logfile,'(A,i3,2A,2f16.4)') ' ',status_nc2,trim(var_name_nc2(i)),' (min, max): ' &
,minval(var2d_nc2(i,:,:)),maxval(var2d_nc2(i,:,:))
elseif (i.eq.elevation_index2) then
!status_nc2 = NF_GET_VARA_DOUBLE (id_nc2, var_id_nc2(i), dim_start_nc2(1:2), dim_length_nc2(1:2), var2d_nc2_dp);var2d_nc2(i,:,:)=real(var2d_nc2_dp)
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var2d_nc2_dp(:,:), start=(/dim_start_nc2/), count=(/dim_length_nc2/));var2d_nc2(i,:,:)=real(var2d_nc2_dp)
write(unit_logfile,'(A,i3,2A,2f16.4)') ' ',status_nc2,trim(var_name_nc2(i)),' (min, max): ' &
,minval(var2d_nc2(i,:,:)),maxval(var2d_nc2(i,:,:))
!write(*,*) maxval(var2d_nc2(i,:,:))
else
!write(*,*) id_nc2, var_id_nc2(i)
!write(*,*) dim_start_nc2
!write(*,*) dim_length_nc2
!Due to memory problems must loop through time on this variable
!do t=1,dim_length_nc2(time_index2)
dim_start_nc2(time_index2)=1
dim_length_nc2(time_index2)=1
!status_nc2 = NF_GET_VARA_DOUBLE (id_nc2, var_id_nc2(i), dim_start_nc2, dim_length_nc2, var3d_nc2_dp);var3d_nc2(i,:,:,t)=real(var3d_nc2_dp)
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var3d_nc2_dp(:,:), start=(/dim_start_nc2/), count=(/dim_length_nc2/));var3d_nc2(i,:,:,t)=real(var3d_nc2_dp)
!enddo
dim_start_nc2(time_index2)=1
dim_length_nc2(time_index2)=dim_length_nc2(time_index2)
write(unit_logfile,'(A,i3,2A,2f16.2)') ' ',status_nc2,trim(var_name_nc2(i)),' (min, max): ' &
,minval(var3d_nc2(i,:,:,t)),maxval(var3d_nc2(i,:,:,t))
endif
else
write(unit_logfile,'(8A,8A)') ' Cannot read ',trim(var_name_nc2(i))
meteo_var_nc2_available(t,i)=.false.
endif
if (i.eq.precip_index2) then
!Don't allow precip below the cutoff value
where (var3d_nc2(i,:,:,t).lt.precip_cutoff) var3d_nc2(i,:,:,t)=0.
endif
enddo
status_nc2 = NF90_CLOSE (id_nc2)
!Put in some basic data checks to see if file is corrupt
if (abs(maxval(var3d_nc2(temperature_index2,:,:,:))).gt.500) then
!write(unit_logfile,'(A,e12.2)') ' ERROR: out of bounds temperature: ', maxval(var3d_nc2(temperature_index,:,:,:))
!write(unit_logfile,'(A)') ' STOPPING'
write(unit_logfile,'(A,e12.2)') ' WARNING: out of bounds temperature. Will not use these data but will continue calculations: ', maxval(var3d_nc2(temperature_index2,:,:,:))
meteo_var_nc2_available(temperature_index2,t)=.false.
!stop
endif
if (abs(maxval(var3d_nc2(precip_index2,:,:,:))).gt.1000) then
!write(unit_logfile,'(A,e12.2)') ' ERROR: out of bounds temperature: ', maxval(var3d_nc2(temperature_index,:,:,:))
!write(unit_logfile,'(A)') ' STOPPING'
write(unit_logfile,'(A,e12.2)') ' WARNING: out of bounds precipitation. Will not use these data but will continue calculations: ', maxval(var3d_nc2(precip_index2,:,:,:))
meteo_var_nc2_available(precip_index2,t)=.false.
!stop
endif
if (abs(maxval(var3d_nc2(relhumidity_index2,:,:,:))).gt.1.10) then
write(unit_logfile,'(A,e12.2)') ' WARNING: out of bounds humidity. Will not use these data but will continue calculations: ', maxval(var3d_nc2(relhumidity_index2,:,:,:))
meteo_var_nc2_available(relhumidity_index2,t)=.false.
!stop
endif
i_grid_mid=int(dim_length_nc2(x_index)/2)
j_grid_mid=int(dim_length_nc2(y_index)/2)
dgrid_nc2(x_index2)=var1d_nc2(x_index2,i_grid_mid)-var1d_nc2(x_index2,i_grid_mid-1)
dgrid_nc2(y_index2)=var1d_nc2(y_index2,j_grid_mid)-var1d_nc2(y_index2,j_grid_mid-1)
!If the coordinates are in km instead of metres then change to metres (assuming the difference is not going to be > 100 km
if (dgrid_nc2(x_index).lt.100.and.meteo_nc2_projection_type.ne.LL_projection_index) then
dgrid_nc2=dgrid_nc2*1000.
var1d_nc2(x_index2,:)=var1d_nc2(x_index2,:)*1000.
var1d_nc2(y_index2,:)=var1d_nc2(y_index2,:)*1000.
endif
write(unit_logfile,'(A,2f12.1)') ' Grid spacing X and Y (m): ', dgrid_nc2(x_index2),dgrid_nc2(y_index2)
!Set the array dimensions to the available ones. Can be changed later based on input information, particularly for time
end_dim_nc2=dim_length_nc2
start_dim_nc2=dim_start_nc2
endif !End if exists
enddo !end t loop
if (allocated(var3d_nc2_dp)) deallocate (var3d_nc2_dp)
if (allocated(var2d_nc2_dp)) deallocate (var2d_nc2_dp)
end subroutine NORTRIP_read_analysismeteo_netcdf4
subroutine NORTRIP_read_t2m500yr_netcdf4
!Reads yr temperature 66 hour forecast data
!This is an older version and is no longer produced
use NORTRIP_multiroad_index_definitions
!Update to netcdf 4 and 64 bit in this version 2 of NORTRIP_read_t2m500yr_netcdf
use netcdf
implicit none
!include 'netcdf.inc'
!Local variables
integer status_nc2 !Error message
integer id_nc2
integer dim_id_nc2(num_dims_nc2)
integer xtype_nc2(num_var_nc2)
integer natts_nc2(num_var_nc2)
integer var_id_nc2(num_var_nc2)
integer dim_length_metcoop_nc2(num_dims_nc2+1)
integer dim_start_metcoop_nc2(num_dims_nc2+1)
character(256) dimname_temp
integer i
integer i_grid_mid,j_grid_mid
real dlat_nc2
integer exists
integer ii,jj,tt,t
double precision, allocatable :: var2d_nc2_dp(:,:)
double precision, allocatable :: var3d_nc2_dp(:,:)
write(unit_logfile,'(A)') '================================================================'
write(unit_logfile,'(A)') 'Reading meteorological data (NORTRIP_read_t2m500yr_netcdf v2)'
write(unit_logfile,'(A)') '================================================================'
!pathname_nc='C:\BEDRE BYLUFT\NORTRIP implementation\test\';
!filename_nc='AROME_1KM_OSLO_20141028_EPI.nc'
pathfilename_nc2=trim(pathname_nc2)//trim(filename_nc2)
!Test existence of the filename. If does not exist then use default
inquire(file=trim(pathfilename_nc2),exist=exists)
if (.not.exists) then
write(unit_logfile,'(A,A)') ' ERROR: Meteo netcdf file does not exist: ', trim(pathfilename_nc2)
write(unit_logfile,'(A)') ' STOPPING'
!write(*,'(A,A)') ' ERROR: Meteo netcdf file does not exist. Stopping: ', trim(pathfilename_nc)
stop 30
endif
!Open the netcdf file for reading
write(unit_logfile,'(2A)') ' Opening netcdf meteo file: ',trim(pathfilename_nc2)
status_nc2 = NF90_OPEN (pathfilename_nc2, NF90_NOWRITE, id_nc2)
if (status_nc2 .NE. NF90_NOERR) write(unit_logfile,'(A,I)') 'ERROR opening netcdf file: ',status_nc2
!Find out the x,y and time dimmensions of the file
status_nc2 = NF90_INQ_DIMID (id_nc2,dim_name_nc2(x_index2),dim_id_nc2(x_index2))
status_nc2 = NF90_INQUIRE_DIMENSION (id_nc2,dim_id_nc2(x_index2),dimname_temp,dim_length_nc2(x_index2))
status_nc2 = NF90_INQ_DIMID (id_nc2,dim_name_nc2(y_index2),dim_id_nc2(y_index2))
status_nc2 = NF90_INQUIRE_DIMENSION (id_nc2,dim_id_nc2(y_index2),dimname_temp,dim_length_nc2(y_index2))
status_nc2 = NF90_INQ_DIMID (id_nc2,dim_name_nc2(time_index2),dim_id_nc2(time_index2))
status_nc2 = NF90_INQUIRE_DIMENSION (id_nc2,dim_id_nc2(time_index2),dimname_temp,dim_length_nc2(time_index2))
write(unit_logfile,'(A,3I)') ' Size of dimensions (x,y,t): ',dim_length_nc2
if (number_of_time_steps.ne.0) then
dim_length_nc2(time_index)=number_of_time_steps
write(unit_logfile,'(A,3I)') ' WARNING: Reducing dimensions of (t) to save space: ',dim_length_nc2(time_index)
endif
!Allocate the nc arrays for reading
!write(*,*) dim_length_nc2(x_index2),dim_length_nc2(y_index2),dim_length_nc2(time_index2)
allocate (var1d_nc2(num_dims_nc2,maxval(dim_length_nc2))) !x and y and time maximum dimmensions
allocate (var3d_nc2(temperature_index2,dim_length_nc2(x_index2),dim_length_nc2(y_index2),dim_length_nc2(time_index2)))
allocate (var2d_nc2(elevation_index2,dim_length_nc2(x_index2),dim_length_nc2(y_index2))) !Lat and lon and elevation
allocate (var3d_nc2_dp(dim_length_nc2(x_index2),dim_length_nc2(y_index2)))
allocate (var2d_nc2_dp(dim_length_nc2(x_index2),dim_length_nc2(y_index2))) !Lat and lon
!Set the number of hours to be read
!Read the x, y and time values
do i=1,num_dims_nc2
status_nc2 = NF90_INQ_VARID (id_nc2, trim(dim_name_nc2(i)), var_id_nc2(i))
!status_nc2 = NF_GET_VARA_REAL (id_nc2, var_id_nc2(i), dim_start_nc2(i), dim_length_nc2(i), var1d_nc2(i,:))
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var1d_nc2(i,1:dim_length_nc2(i)), start=(/dim_start_nc2(i)/), count=(/dim_length_nc2(i)/))
if (i.eq.time_index2) then
write(unit_logfile,'(3A,2i12)') ' ',trim(dim_name_nc2(i)),' (min, max in hours): ' &
,minval(int((var1d_nc2(i,1:dim_length_nc2(i))-var1d_nc2(i,dim_start_nc2(i)))/3600.+.5)+1) &
,maxval(int((var1d_nc2(i,1:dim_length_nc2(i))-var1d_nc2(i,dim_start_nc2(i)))/3600.+.5)+1)
else
write(unit_logfile,'(3A,2f12.2)') ' ',trim(dim_name_nc2(i)),' (min, max in km): ' &
,minval(var1d_nc2(i,1:dim_length_nc2(i))),maxval(var1d_nc2(i,1:dim_length_nc2(i)))
endif
enddo
!MetCoOp data has 3 dimensions. z is 1 so must adjust this.
!dim_length_metcoop_nc(1:2)=dim_length_nc(1:2)
!dim_length_metcoop_nc(3)=1
!dim_length_metcoop_nc(4)=dim_length_nc(time_index)
!dim_start_metcoop_nc(1:2)=dim_start_nc(1:2)
!dim_start_metcoop_nc(3)=1
!dim_start_metcoop_nc(4)=dim_start_nc(time_index)
!Read through the variables in a loop
do i=1,num_var_nc2
!write(*,*) i,trim(var_name_nc(i))
status_nc2 = NF90_INQ_VARID (id_nc2, trim(var_name_nc2(i)), var_id_nc2(i))
!write(*,*) 'Status1: ',status_nc2,id_nc2,var_id_nc2(i),trim(var_name_nc2(i)),NF_NOERR
!write(*,*) 'Status1: ',dim_start_metcoop_nc
!write(*,*) 'Status1: ',dim_length_metcoop_nc
if (status_nc2.eq.NF90_NOERR) then
if (i.eq.lat_index2.or.i.eq.lon_index2) then
!write(*,*) id_nc2, var_id_nc2(i), dim_start_nc2(1:2), dim_length_nc2(1:2)
!status_nc = NF_GET_VARA_REAL (id_nc, var_id_nc(i), dim_start_metcoop_nc(1:2), dim_length_metcoop_nc(1:2), var2d_nc(i,:,:))
!status_nc2 = NF_GET_VARA_DOUBLE (id_nc2, var_id_nc2(i), dim_start_nc2(1:2), dim_length_nc2(1:2), var2d_nc2_dp);var2d_nc2(i,:,:)=real(var2d_nc2_dp)
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var2d_nc2_dp(:,:), start=(/dim_start_nc2/), count=(/dim_length_nc2/));var2d_nc2(i,:,:)=real(var2d_nc2_dp)
write(unit_logfile,'(A,i3,2A,2f16.4)') ' ',status_nc2,trim(var_name_nc2(i)),' (min, max): ' &
,minval(var2d_nc2(i,:,:)),maxval(var2d_nc2(i,:,:))
elseif (i.eq.elevation_index2) then
!status_nc2 = NF_GET_VARA_DOUBLE (id_nc2, var_id_nc2(i), dim_start_nc2(1:2), dim_length_nc2(1:2), var2d_nc2_dp);var2d_nc2(i,:,:)=real(var2d_nc2_dp)
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var2d_nc2_dp(:,:), start=(/dim_start_nc2/), count=(/dim_length_nc2/));var2d_nc2(i,:,:)=real(var2d_nc2_dp)
write(unit_logfile,'(A,i3,2A,2f16.4)') ' ',status_nc2,trim(var_name_nc2(i)),' (min, max): ' &
,minval(var2d_nc2(i,:,:)),maxval(var2d_nc2(i,:,:))
!write(*,*) maxval(var2d_nc2(i,:,:))
else
!write(*,*) id_nc2, var_id_nc2(i)
!write(*,*) dim_start_nc2
!write(*,*) dim_length_nc2
!Due to memory problems must loop through time on this variable
do t=1,dim_length_nc2(time_index2)
dim_start_nc2(time_index2)=t
dim_length_nc2(time_index2)=1
!status_nc2 = NF_GET_VARA_DOUBLE (id_nc2, var_id_nc2(i), dim_start_nc2, dim_length_nc2, var3d_nc2_dp);var3d_nc2(i,:,:,t)=real(var3d_nc2_dp)
status_nc2 = NF90_GET_VAR (id_nc2, var_id_nc2(i), var3d_nc2_dp(:,:), start=(/dim_start_nc2/), count=(/dim_length_nc2/));var3d_nc2(i,:,:,t)=real(var3d_nc2_dp)
enddo
dim_start_nc2(time_index2)=1
dim_length_nc2(time_index2)=dim_length_nc2(time_index2)
write(unit_logfile,'(A,i3,2A,2f16.2)') ' ',status_nc2,trim(var_name_nc2(i)),' (min, max): ' &
,minval(var3d_nc2(i,:,:,:)),maxval(var3d_nc2(i,:,:,:))
endif
else
write(unit_logfile,'(8A,8A)') ' Cannot read ',trim(var_name_nc2(i))
endif
enddo
status_nc2 = NF90_CLOSE (id_nc2)
!Set the array dimensions to the available ones. Can be changed later based on input information, particularly for time
end_dim_nc2=dim_length_nc2
start_dim_nc2=dim_start_nc2
if (allocated(var3d_nc2_dp)) deallocate (var3d_nc2_dp)
if (allocated(var2d_nc2_dp)) deallocate (var2d_nc2_dp)
end subroutine NORTRIP_read_t2m500yr_netcdf4