Home

Context Navigation

calc_obstruction.f90

Last change on this file was 4870, checked in by suehring, 4 years ago
preprocessor to calculate UV-obstruction added
Property svn:keywords set to `Id`
File size: 21.1 KB

Line
1	!> @file calc_obstruction.f90
2	!------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the
6	! terms of the GNU General Public License as published by the Free Software
7	! Foundation, either version 3 of the License, or (at your option) any later
8	! version.
9	!
10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
13	!
14	! You should have received a copy of the GNU General Public License along with
15	! PALM. If not, see <http://www.gnu.org/licenses/>.
16	!
17	! Copyright 1997-2021 Leibniz Universitaet Hannover
18	!------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! -----------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: calc_obstruction.f90 4870 2021-02-08 18:09:41Z banzhafs $
27	! Initial revision
28	!
29	!
30	! Authors:
31	! --------
32	! @author Michael Schrempf (main author),
33	! Matthias SÃŒhring (first major revision of the code)
34	!
35	!> @todo - Code formatting
36	!------------------------------------------------------------------------------!
37	! Description:
38	! ------------
39	!> Calculation of sky-obstruction.
40	!------------------------------------------------------------------------------!
41
42	PROGRAM obstruction
43	#if defined ( __netcdf )
44	USE netcdf
45	#endif
46
47	IMPLICIT NONE
48
49	CHARACTER(LEN=10) :: b(3)
50	CHARACTER(LEN=300) :: input_file
51	CHARACTER(LEN=300) :: path_to_expo_files
52
53	INTEGER:: spx, spy, counter, cc, sza, rza, azi, str_i, str_j, ii, &
54	vector_length_x, vector_length_y, id, id_var, nc_stat, i, j, k
55	INTEGER(KIND=1) :: bitvar
56	INTEGER, DIMENSION(1:8):: date_time
57	INTEGER, DIMENSION(0:99):: valy, valx
58
59	INTEGER :: obs_ny = 9 !< number of zenith values
60	INTEGER :: obs_nx = 35 !< number of azimuth values
61	INTEGER :: obs_nz = 44 !< ?
62
63	INTEGER(KIND=2), DIMENSION(:,:,:), ALLOCATABLE :: obs
64	INTEGER(KIND=2), DIMENSION(:,:,:), ALLOCATABLE :: obsi
65
66	! INTEGER(KIND=2), DIMENSION(0:1159,0:1159,0:44):: obs
67	! INTEGER(KIND=2), DIMENSION(0:959,0:959,0:44):: obsi
68
69	INTEGER :: nx !< number of grid points in x-direction
70	INTEGER :: ny !< number of grid points in y-direction
71	INTEGER :: nsza
72	REAL :: dx !< grid spacing in x
73	REAL :: dy !< grid spacing in y
74	REAL :: dz !< grid spacing in z
75
76	REAL :: delta = 10.0, fill !< delta of zenith and azimuth angles (namelist param)
77
78	REAL:: dtor, pi = 3.141592654, t1 = 0.0, t2 = 0.0
79	REAL, DIMENSION(0:99):: value_y, value_x, value_z, value_zn, distance_x, distance_y
80	! REAL, DIMENSION(0:1159,0:1159):: topo
81	! REAL, DIMENSION(0:959,0:959):: topo_0
82	REAL, DIMENSION(:,:), ALLOCATABLE :: topo
83	REAL, DIMENSION(:,:), ALLOCATABLE :: topo_0
84
85	REAL(kind=8), DIMENSION(0:35,0:9) :: integration
86	REAL, DIMENSION(0:2,0:90) :: irrad
87	REAL, DIMENSION(0:35,0:9) :: tmp
88	REAL, DIMENSION(0:35,0:9,0:90) :: cube_rad
89	REAL, DIMENSION(0:35,0:9,0:2) :: pfs
90
91	CHARACTER(LEN=15) :: output ='Data_input_all'
92
93	CALL obstruct_parin
94
95	nsza = 90
96
97	dtor = pi / 180.
98	vector_length_x = INT( 99.0 / dx )
99	vector_length_y = INT( 99.0 / dy )
100
101
102	ALLOCATE( obs(0:nx-1+(2vector_length_x),0:ny-1+(2vector_length_y),0:44) )
103	ALLOCATE( obsi(0:nx-1,0:ny-1,0:44) )
104	ALLOCATE( topo(0:nx-1+(2vector_length_x),0:ny-1+(2vector_length_y)) )
105	ALLOCATE( topo_0(0:nx-1,0:ny-1) )
106
107	obs(:,:,:) = 9
108	obsi(:,:,:) = 0
109	topo(:,:) = 0
110
111	call date_and_time(b(1), b(2), b(3), date_time)
112	t1 = date_time(5)+ date_time(6)/60. + (date_time(7)+(date_time(8)/1000.))/3600.
113
114	distance_x = 0.0
115	distance_y = 0.0
116
117	DO ii = 0, 99
118	distance_x(ii) = ii * dx
119	distance_y(ii) = ii * dy
120	ENDDO
121	distance_x = distance_x + dx
122	distance_y = distance_y + dy
123
124	#if defined ( __netcdf )
125	nc_stat = NF90_OPEN( input_file, NF90_NOWRITE, id )
126	nc_stat = NF90_INQ_VARID( id, 'buildings_2d', id_var )
127	nc_stat = NF90_GET_VAR( id, id_var, topo_0, start = (/ 1, 1 /), count = (/ nx, ny /) )
128	nc_stat = NF90_GET_ATT( id, id_var, '_FillValue', fill )
129	nc_stat = NF90_CLOSE( id )
130	#endif
131
132	! OPEN(11, FILE='/home/matthias/UNI/MOSAIK/UV/test/Topo_input/topo_Reuter.txt', FORM='FORMATTED')
133	! DO rza = 0, nx-1
134	! READ(11,'(960(F5.1,1X))') topo_0(:,rza) !(nx-1)-rza
135	! ENDDO
136	! CLOSE(11)
137
138	topo(vector_length_x:nx-1+vector_length_x,vector_length_y:ny-1+vector_length_y) = &
139	MERGE( topo_0, 0.0, topo_0 /= fill )
140
141
142	Do spy = vector_length_y, ny-1+vector_length_y
143	DO spx = vector_length_x, nx-1+vector_length_x
144	print*, spy, spx, topo(spx,spy)
145	IF (topo(spx,spy) == 0.0 ) THEN
146	! print*, "in if"
147	obs(spx,spy,:)=0
148	cc = 0
149	bitvar=2
150	DO sza = 0,90, 10
151	! print*, sza
152	DO azi = 0,359, 10
153	! print*, azi
154	value_zn = SIN((sza+0.0000001)*dtor)
155
156	value_x = ANINT((distance_x(:)SIN(azidtor)) *value_zn(:))
157	value_y = ANINT((distance_y(:)COS(azidtor)) *value_zn(:))
158	! print*, value_x
159	! print*, value_y
160	valx = INT(value_x)
161	valy = INT(value_y)
162	! print*, "valx", valx
163	!value_z = (distance(:)COS(szadtor))
164	value_z = COS(szadtor) SQRT( value_x2 + value_y2) / value_zn
165	! print, szadtor, sza
166	!print*, value_z
167	! stop
168	counter=0
169	DO str_i = 0, vector_length_x-1
170	DO str_j = 0, vector_length_y-1
171	IF ( topo( spx+valx(str_i), spy+valy(str_j) ) >= (value_z(str_i)*dz) ) THEN
172	IF ( topo( spx+valx(str_i), spy+valy(str_j) ) == 0) THEN
173	counter=counter + 0
174	ELSE
175	counter=counter + 1
176	ENDIF
177	ENDIF
178	ENDDO
179	ENDDO
180	IF (sza == 90 ) Counter = 1
181	IF (counter == 0) bitvar = IBSET(bitvar,MOD(cc,8))
182	IF (counter > 0) bitvar = IBCLR(bitvar,MOD(cc,8))
183	!IF (counter .GT. 0) bitvar = IBSET(bitvar,MOD(cc,8))
184
185
186	IF ( MOD(cc,8) == 7 ) THEN
187	obs(spx,spy,cc/8) = bitvar
188	bitvar = 2
189	ENDIF
190	cc = cc+1
191	ENDDO
192	ENDDO
193	ENDIF
194	ENDDO
195	ENDDO
196
197	call date_and_time(b(1), b(2), b(3), date_time)
198	t2 = date_time(5)+ date_time(6)/60. + (date_time(7)+(date_time(8)/1000.))/3600.
199
200	obsi(:,:,:)=obs(vector_length_x:nx-1+vector_length_x,vector_length_y:ny-1+vector_length_y,:)
201
202
203
204	OPEN(11, FILE= TRIM( path_to_expo_files ) // 'Integration_solid_angles.txt', FORM='FORMATTED')
205	DO rza = 0,9
206	READ(11,'(36(F10.8,2X))') integration(:,rza)
207	ENDDO
208	CLOSE(11)
209
210	OPEN(11, FILE= TRIM( path_to_expo_files ) // 'Irradiance_ROM.txt', FORM='FORMATTED')
211	Do sza = 0,90
212	READ(11,'(3(F10.6,2X))') irrad(:,sza)
213	ENDDO
214	CLOSE(11)
215
216
217
218	OPEN(11, FILE= TRIM( path_to_expo_files ) // 'Projection_areas_CL0.txt', FORM='FORMATTED')
219	DO rza = 0,9
220	READ(11,'(36(F10.6,2X))') tmp(:,rza)
221	ENDDO
222	CLOSE(11)
223	pfs(:,:,0) = tmp(:,:)
224
225	OPEN(11, FILE=TRIM( path_to_expo_files ) // 'Projection_areas_CL1.txt', FORM='FORMATTED')
226	DO rza = 0,9
227	READ(11,'(36(F10.6,2X))') tmp(:,rza)
228	ENDDO
229	CLOSE(11)
230	pfs(:,:,1) = tmp(:,:)
231
232	OPEN(11, FILE= TRIM( path_to_expo_files ) // 'Projection_areas_CL3.txt', FORM='FORMATTED')
233	DO rza = 0,9
234	READ(11,'(36(F10.6,2X))') tmp(:,rza)
235	ENDDO
236	CLOSE(11)
237	pfs(:,:,2) = tmp(:,:)
238
239
240	OPEN(11, FILE= TRIM( path_to_expo_files ) // 'Radiance_ROM.txt', FORM='FORMATTED')
241	Do sza = 0,90
242	DO rza = 0,9
243	READ(11,'(36(F10.6,2X))') tmp(:,rza)
244	ENDDO
245	cube_rad(:,:,sza) = tmp(:,:)
246	ENDDO
247	CLOSE(11)
248
249	!- Finalize output
250	CALL data_output('open', output)
251	CALL data_output('write', output)
252	CALL data_output('close', output)
253
254
255	!- End
256	WRITE(*, '(A)') ' --- creation finished! :-) ---'
257
258
259	! open(unit=12, status='replace',file=&
260	! '/home/matthias/UNI/MOSAIK/UV/test/result_shading/obstruction_Berlin_reuter.bin', FORM='UNFORMATTED')
261	! write(12) obsi
262	! close(12)
263
264	PRINT *,'Time_Start = ',t1
265	PRINT *,'Time_ENDE = ',t2
266	PRINT ,'elapsed time used (dateT) = ',(t2-t1)3600.
267
268
269	CONTAINS
270
271	SUBROUTINE obstruct_parin
272
273	IMPLICIT NONE
274
275	INTEGER :: file_id_parin = 90
276
277	NAMELIST /uv_parin/ dx, dy, dz, nx, ny, input_file, path_to_expo_files
278	!
279	!-- Open namelist file.
280	OPEN( file_id_parin, FILE='uv_parin', STATUS='OLD', FORM='FORMATTED')
281	!
282	!-- Read namelist.
283	READ ( file_id_parin, uv_parin )
284	!
285	!-- Close namelist file.
286	CLOSE( file_id_parin )
287
288	END SUBROUTINE obstruct_parin
289
290
291	SUBROUTINE data_output(action,output)
292	! Data output to NetCDF file
293
294	CHARACTER(LEN=*), INTENT(IN) :: action !< flag to steer routine (open/close/write)
295	CHARACTER(LEN=*), INTENT(IN) :: output !< output file name prefix
296	! CHARACTER(LEN=10) :: output ='Data_input'
297
298	INTEGER, SAVE :: do_count=0 !< counting output
299	INTEGER, SAVE :: id_dim_time !< ID of dimension time
300	INTEGER, SAVE :: id_dim_x !< ID of dimension x
301	INTEGER, SAVE :: id_dim_xirr !< ID of dimension x
302	INTEGER, SAVE :: id_dim_y !< ID of dimension y
303	INTEGER, SAVE :: id_dim_z !< ID of dimension z
304	INTEGER, SAVE :: id_dim_obsx !< ID of dimension x
305	INTEGER, SAVE :: id_dim_obsy !< ID of dimension y
306	INTEGER, SAVE :: id_dim_obsz !< ID of dimension z
307	INTEGER, SAVE :: id_file !< ID of NetCDF file
308	INTEGER, SAVE :: id_var_int !< ID of geopotential
309	INTEGER, SAVE :: id_var_irr !< ID of geopotential
310	INTEGER, SAVE :: id_var_obs !< ID of geopotential
311	INTEGER, SAVE :: id_var_obsx !< ID of x
312	INTEGER, SAVE :: id_var_obsy !< ID of y
313	INTEGER, SAVE :: id_var_obsz !< ID of y
314
315	INTEGER, SAVE :: id_var_pfs !< ID of geopotential
316	INTEGER, SAVE :: id_var_rad !< ID of geopotential
317	INTEGER, SAVE :: id_var_x !< ID of x
318	INTEGER, SAVE :: id_var_xirr !< ID of x
319	INTEGER, SAVE :: id_var_y !< ID of y
320	INTEGER, SAVE :: id_var_z !< ID of y
321	INTEGER, SAVE :: nc_stat !< status flag of NetCDF routines
322
323	REAL, DIMENSION(:), ALLOCATABLE :: netcdf_data_1d !< 1D output array
324
325	REAL, DIMENSION(:,:), ALLOCATABLE :: netcdf_data_2d !< 2D output array
326
327	REAL, DIMENSION(:,:,:), ALLOCATABLE :: netcdf_data_3d !< 3D output array
328
329
330	SELECT CASE (TRIM(action))
331
332	!- Initialize output file
333	CASE ('open')
334
335	!- Delete any pre-existing output file
336	OPEN(20, FILE=TRIM(output)//'.nc')
337	CLOSE(20, STATUS='DELETE')
338
339	!- Open file
340	#if defined ( __netcdf )
341	nc_stat = NF90_CREATE(TRIM(output)//'.nc', NF90_NOCLOBBER, id_file)
342	IF (nc_stat /= NF90_NOERR) PRINT*, '+++ netcdf error'
343
344	!- Write global attributes
345	nc_stat = NF90_PUT_ATT(id_file, NF90_GLOBAL, &
346	'Conventions', 'COARDS')
347	nc_stat = NF90_PUT_ATT(id_file, NF90_GLOBAL, &
348	'title', 'UVEM-model')
349
350	!- Define time coordinate
351	nc_stat = NF90_DEF_DIM(id_file, 'TIME', NF90_UNLIMITED, &
352	id_dim_time)
353
354
355	!- Define spatial coordinates
356	nc_stat = NF90_DEF_DIM(id_file, 'AZI', obs_nx+1, id_dim_x)
357	nc_stat = NF90_DEF_VAR(id_file, 'AZI', NF90_DOUBLE, id_dim_x, id_var_x)
358	nc_stat = NF90_PUT_ATT(id_file, id_var_x, 'units', 'Grad')
359
360	nc_stat = NF90_DEF_DIM(id_file, 'ZEN', obs_ny+1, id_dim_y)
361	nc_stat = NF90_DEF_VAR(id_file, 'ZEN', NF90_DOUBLE, &
362	id_dim_y, id_var_y)
363	nc_stat = NF90_PUT_ATT(id_file, id_var_y, 'units', 'Grad')
364
365	nc_stat = NF90_DEF_DIM(id_file, 'sza', nsza+1, id_dim_z)
366	nc_stat = NF90_DEF_VAR(id_file, 'sza', NF90_DOUBLE, &
367	id_dim_z, id_var_z)
368	nc_stat = NF90_PUT_ATT(id_file, id_var_z, 'units', 'Grad')
369
370	!- Define spatial coordinates
371	nc_stat = NF90_DEF_DIM(id_file, 'COM', 3, id_dim_xirr)
372	nc_stat = NF90_DEF_VAR(id_file, 'COM', NF90_DOUBLE, id_dim_xirr, id_var_xirr)
373	nc_stat = NF90_PUT_ATT(id_file, id_var_xirr, 'units', 'Komponente')
374
375
376	!- Define spatial coordinates for obstruction file
377	nc_stat = NF90_DEF_DIM(id_file, 'X', nx-1+1, id_dim_obsx)
378	nc_stat = NF90_DEF_VAR(id_file, 'X', NF90_DOUBLE, id_dim_obsx, id_var_obsx)
379	nc_stat = NF90_PUT_ATT(id_file, id_var_obsx, 'units', 'meters')
380
381	nc_stat = NF90_DEF_DIM(id_file, 'Y', ny-1+1, id_dim_obsy)
382	nc_stat = NF90_DEF_VAR(id_file, 'Y', NF90_DOUBLE, id_dim_obsy, id_var_obsy)
383	nc_stat = NF90_PUT_ATT(id_file, id_var_obsy, 'units', 'meters')
384
385	nc_stat = NF90_DEF_DIM(id_file, 'directions', nsza+1, id_dim_obsz)
386	nc_stat = NF90_DEF_VAR(id_file, 'directions', NF90_DOUBLE, id_dim_obsz, id_var_obsz)
387	nc_stat = NF90_PUT_ATT(id_file, id_var_obsz, 'units', 'ibits')
388
389
390
391	!- Define output arrays
392	nc_stat = NF90_DEF_VAR(id_file, 'radiance', NF90_DOUBLE, &
393	(/id_dim_x, id_dim_y, id_dim_z/), id_var_rad)
394	nc_stat = NF90_PUT_ATT(id_file, id_var_rad, &
395	'long_name', 'vitD weighetd radiance at 300DU')
396	nc_stat = NF90_PUT_ATT(id_file, id_var_rad, &
397	'short_name', 'vitD radiance')
398	nc_stat = NF90_PUT_ATT(id_file, id_var_rad, 'units', 'W/m2 sr')
399
400
401
402	nc_stat = NF90_DEF_VAR(id_file, 'int_factors', NF90_DOUBLE, &
403	(/id_dim_x, id_dim_y/), id_var_int)
404	nc_stat = NF90_PUT_ATT(id_file, id_var_int, &
405	'long_name', 'integration factors')
406	nc_stat = NF90_PUT_ATT(id_file, id_var_int, &
407	'short_name', 'integration factors')
408	! nc_stat = NF90_PUT_ATT(id_file, id_var_int, 'units', '-')
409
410
411	!- Define output arrays
412	nc_stat = NF90_DEF_VAR(id_file, 'projarea', NF90_DOUBLE, &
413	(/id_dim_x, id_dim_y, id_dim_time/), id_var_pfs)
414	nc_stat = NF90_PUT_ATT(id_file, id_var_pfs, &
415	'long_name', 'projection areas of a human')
416	nc_stat = NF90_PUT_ATT(id_file, id_var_pfs, &
417	'short_name', 'projection areas')
418	nc_stat = NF90_PUT_ATT(id_file, id_var_pfs, 'units', 'm2')
419
420
421	!- Define output arrays
422	nc_stat = NF90_DEF_VAR(id_file, 'irradiance', NF90_DOUBLE, &
423	(/id_dim_xirr, id_dim_z/), id_var_irr)
424	nc_stat = NF90_PUT_ATT(id_file, id_var_irr, &
425	'long_name', 'global, diffuse and direct solar irradiance')
426	nc_stat = NF90_PUT_ATT(id_file, id_var_irr, &
427	'short_name', 'irradiance values')
428	nc_stat = NF90_PUT_ATT(id_file, id_var_irr, 'units', 'W/m2')
429
430
431	!- Define output arrays for obstruction file
432	nc_stat = NF90_DEF_VAR(id_file, 'obstruction', NF90_DOUBLE, &
433	(/id_dim_obsx, id_dim_obsy, id_dim_obsz/), id_var_obs)
434	nc_stat = NF90_PUT_ATT(id_file, id_var_obs, &
435	'long_name', 'calculated obstructions in different directions')
436	nc_stat = NF90_PUT_ATT(id_file, id_var_obs, &
437	'short_name', 'obstructions')
438	nc_stat = NF90_PUT_ATT(id_file, id_var_obs, 'units', '-')
439
440
441
442
443	!- Leave define mode
444	nc_stat = NF90_ENDDEF(id_file)
445
446	!- Write axis to file
447	ALLOCATE(netcdf_data_1d(0:nx))
448
449	!- x axis
450	DO i = 0, nx
451	netcdf_data_1d(i) = i * delta
452	ENDDO
453	nc_stat = NF90_PUT_VAR(id_file, id_var_x, netcdf_data_1d, &
454	start = (/1/), count = (/nx+1/))
455	DO i = 0, 2
456	netcdf_data_1d(i) = i
457	ENDDO
458	nc_stat = NF90_PUT_VAR(id_file, id_var_xirr, netcdf_data_1d, &
459	start = (/1/), count = (/3/))
460
461	!- y axis
462	DO j = 0, ny
463	netcdf_data_1d(j) = j * delta
464	ENDDO
465	nc_stat = NF90_PUT_VAR(id_file, id_var_y, netcdf_data_1d, &
466	start = (/1/), count = (/ny+1/))
467	DEALLOCATE(netcdf_data_1d)
468
469	!- z axis
470	ALLOCATE(netcdf_data_1d(0:nsza))
471	DO k = 0, nsza
472	netcdf_data_1d(k) = k
473	ENDDO
474	nc_stat = NF90_PUT_VAR(id_file, id_var_z, netcdf_data_1d, &
475	start = (/1/), count = (/nsza+1/))
476	DEALLOCATE(netcdf_data_1d)
477
478
479
480	!- define axis for obstruction file
481	!- x axis
482	ALLOCATE(netcdf_data_1d(0:nx-1))
483	DO i = 0, nx-1
484	netcdf_data_1d(i) = i * dx
485	ENDDO
486	nc_stat = NF90_PUT_VAR(id_file, id_var_obsx, netcdf_data_1d, &
487	start = (/1/), count = (/nx-1+1/))
488	!- y axis
489	DO j = 0, ny-1
490	netcdf_data_1d(j) = j * dy
491	ENDDO
492	nc_stat = NF90_PUT_VAR(id_file, id_var_obsy, netcdf_data_1d, &
493	start = (/1/), count = (/ny-1+1/))
494	DEALLOCATE(netcdf_data_1d)
495
496	!- z axis
497	ALLOCATE(netcdf_data_1d(0:obs_nz))
498	DO k = 0, obs_nz
499	netcdf_data_1d(k) = k
500	ENDDO
501	nc_stat = NF90_PUT_VAR(id_file, id_var_obsz, netcdf_data_1d, &
502	start = (/1/), count = (/obs_nz+1/))
503	DEALLOCATE(netcdf_data_1d)
504	#endif
505
506
507	!- Close NetCDF file
508	CASE ('close')
509	#if defined ( __netcdf )
510	nc_stat = NF90_CLOSE(id_file)
511	#endif
512
513	!- Write data arrays to file
514	CASE ('write')
515	WRITE(*, '(A)') ' --- NetCDF Write! :-) ---'
516
517	! ALLOCATE(netcdf_data_2d(0:nx,0:ny))
518	ALLOCATE(netcdf_data_3d(0:nx,0:ny, 0:nsza))
519
520	!- Write geopotential
521	DO k = 0, nsza
522	DO i = 0, obs_nx
523	DO j = 0, obs_ny
524	netcdf_data_3d(i,j, k) = cube_rad(i,j, k) !normalweise j,i bei Feld2
525	ENDDO
526	ENDDO
527	ENDDO
528	#if defined ( __netcdf )
529	nc_stat = NF90_PUT_VAR(id_file, id_var_rad, netcdf_data_3d, &
530	start = (/1, 1, 1/), &
531	count = (/obs_nx+1, obs_ny+1, nsza+1/))
532
533	nc_stat = NF90_PUT_VAR(id_file, id_var_int, integration, &
534	start = (/1, 1, 1/), &
535	count = (/obs_nx+1, obs_ny+1/))
536
537	nc_stat = NF90_PUT_VAR(id_file, id_var_pfs, pfs, &
538	start = (/1, 1, 1/), &
539	count = (/obs_nx+1, obs_ny+1, 3/))
540
541	nc_stat = NF90_PUT_VAR(id_file, id_var_irr, irrad, &
542	start = (/1, 1, 1/), &
543	count = (/3, nsza+1/))
544
545	DEALLOCATE(netcdf_data_3d)
546
547
548	ALLOCATE(netcdf_data_3d(0:nx-1,0:ny-1, 0:obs_nz))
549
550	!- Write geopotential
551	DO k = 0, obs_nz
552	DO i = 0, nx-1
553	DO j = 0, ny-1
554	netcdf_data_3d(i,j, k) = obsi(i,j, k) !normalweise j,i bei Feld2
555	ENDDO
556	ENDDO
557	ENDDO
558	nc_stat = NF90_PUT_VAR(id_file, id_var_obs, netcdf_data_3d, &
559	start = (/1, 1, 1/), &
560	count = (/nx-1+1, ny-1+1, obs_nz+1/))
561
562	DEALLOCATE(netcdf_data_3d)
563	#endif
564
565
566	!- Print error message if unknown action selected
567	CASE DEFAULT
568
569	WRITE(, '(A)') ' * data_output: action "'// &
570	TRIM(action)//'"unknown!'
571
572	END SELECT
573
574	END SUBROUTINE data_output
575
576	END PROGRAM obstruction

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: palm/trunk/UTIL/uv_obstruction_preprocessor/calc_obstruction.f90

Download in other formats: