1 | !> \file out_cptr_mod.f90 |
---|
2 | !! Module avec les routines d'ecriture et de lecture de fichiers cptr |
---|
3 | !< |
---|
4 | |
---|
5 | !> \namespace out_cptr |
---|
6 | !! This module gathers routines to read and write cptr files |
---|
7 | !! \author Vincent Peyaud |
---|
8 | !! \date janvier 2006 |
---|
9 | !! @note attention : en ecriture fichier runname.cptr |
---|
10 | !! @note en lecture fichier runname.CPTR |
---|
11 | !! @note il faut renommer le fichier ! |
---|
12 | !< |
---|
13 | !---------------------------------------------------------- |
---|
14 | |
---|
15 | |
---|
16 | !> subroutine: read_recovery |
---|
17 | |
---|
18 | !!Reprise d'un fichier compteur |
---|
19 | !! Used module: - use module3D_phy |
---|
20 | !! - use netcdf |
---|
21 | !! - use io_netcdf_grisli |
---|
22 | !! - use tracer_vars |
---|
23 | !! |
---|
24 | ! |
---|
25 | !> |
---|
26 | |
---|
27 | subroutine read_recovery(icpt) !---- reprise d'un fichier compteur ----- |
---|
28 | |
---|
29 | use module3D_phy, only: reprcptr,num_forc,time,S,H,B,Bsoc,ibase,bmelt,hwater,uxbar,uybar,T,ux,uy,hmx,hmy,& |
---|
30 | sdx,sdy,uxflgz,uyflgz,flot,debug_3D |
---|
31 | use geography, only: dx,dy,nx,ny,nz,nzm |
---|
32 | use runparam, only: itracebug,tbegin,num_tracebug |
---|
33 | use io_netcdf_grisli, only: ncdf_type,Read_ncdf_dim,Read_ncdf_var,lect_netcdf_type |
---|
34 | use tracer_vars, only: xdep_out,ydep_out,tdep_out |
---|
35 | implicit none |
---|
36 | |
---|
37 | character(len=20) :: titre |
---|
38 | integer :: nzzm |
---|
39 | integer :: icpt ! icompteur 1-> tout (topo, T, U) |
---|
40 | ! 2 -> tout sauf topo |
---|
41 | ! 3 -> tout sauf topo et eau basale |
---|
42 | |
---|
43 | |
---|
44 | real,dimension(nx,ny) :: bidon !< sert à sauter les lectures de Ss,H,B,... |
---|
45 | |
---|
46 | ! Variables specifiques pour les fichier .nc |
---|
47 | |
---|
48 | real*8, dimension(:,:), pointer :: tab => null() !< tableau 2d real ecrit dans le fichier |
---|
49 | real*8, dimension(:,:,:), pointer :: tab1 => null() !< tableau 3d real |
---|
50 | real*8, dimension(:,:,:), pointer :: tab1T => null() !< tableau 3d real pour la temperature |
---|
51 | real*8, dimension(:), pointer :: liste_time => null() |
---|
52 | |
---|
53 | integer :: ni,nzz,nxx,nyy |
---|
54 | integer :: i,j |
---|
55 | real,dimension(nx,ny) :: oldu |
---|
56 | |
---|
57 | ! rappel de f90 |
---|
58 | ! reprcptr est un string (ici le nom du fichier recovery) |
---|
59 | ! index (reprcptr,'.cptr') renvoie : |
---|
60 | ! 0, si 'cptr' n'est pas dans reprcptr |
---|
61 | ! la position de la premiere occurence de cptr sinon |
---|
62 | |
---|
63 | if (ncdf_type.eq.0) call lect_netcdf_type !< pour lire la valeur de netcdf_type (machine dependant) |
---|
64 | |
---|
65 | |
---|
66 | ascii_nc: if ((index(reprcptr,'.cptr').ne.0) .or.(index(reprcptr,'.CPTR').ne.0)) then !reprise fichier cptr |
---|
67 | |
---|
68 | if (itracebug.eq.1) call tracebug(' Entree dans routine read_recovery : ascii') |
---|
69 | |
---|
70 | open(num_forc,file=trim(reprcptr)) |
---|
71 | read(num_forc,*) time |
---|
72 | |
---|
73 | |
---|
74 | ! temps de la reprise |
---|
75 | !----------------------- |
---|
76 | if (tbegin.gt.1.d9) then ! si tbegin est tres grand , on reprend le temps du cptr |
---|
77 | tbegin=time |
---|
78 | else |
---|
79 | time=tbegin |
---|
80 | endif |
---|
81 | |
---|
82 | |
---|
83 | read(num_forc,*) titre |
---|
84 | read(num_forc,*) ni,nzz,nzzm,nxx,nyy |
---|
85 | read(num_forc,*) |
---|
86 | |
---|
87 | |
---|
88 | if ((nxx.ne.nx).or.(nyy.ne.ny).or.(nzz.ne.nz) & |
---|
89 | .or.(nzzm.ne.nzm)) write(6,*) 'attention', & |
---|
90 | 'les tailles de tableaux ne sont pas compatibles' |
---|
91 | |
---|
92 | ! lecture des tableaux 2D |
---|
93 | read(num_forc,*) |
---|
94 | read(num_forc,*) |
---|
95 | |
---|
96 | lect_topo: if (icpt.eq.1) then ! lecture topo seulement si icompteur = 1 |
---|
97 | read(num_forc,*) S |
---|
98 | read(num_forc,*) |
---|
99 | |
---|
100 | read(num_forc,*) H |
---|
101 | read(num_forc,*) |
---|
102 | |
---|
103 | read(num_forc,*) B |
---|
104 | read(num_forc,*) |
---|
105 | |
---|
106 | read(num_forc,*) Bsoc |
---|
107 | read(num_forc,*) |
---|
108 | |
---|
109 | else |
---|
110 | read(num_forc,*) bidon |
---|
111 | read(num_forc,*) |
---|
112 | |
---|
113 | read(num_forc,*) bidon |
---|
114 | read(num_forc,*) |
---|
115 | |
---|
116 | read(num_forc,*) bidon |
---|
117 | read(num_forc,*) |
---|
118 | |
---|
119 | read(num_forc,*) bidon |
---|
120 | read(num_forc,*) |
---|
121 | |
---|
122 | end if lect_topo |
---|
123 | |
---|
124 | ! lecture temperature dans tous les cas |
---|
125 | read(num_forc,*) Ibase |
---|
126 | read(num_forc,*) |
---|
127 | |
---|
128 | read(num_forc,*) bmelt |
---|
129 | read(num_forc,*) |
---|
130 | |
---|
131 | lect_Hwater: if (icpt.ne.3) then ! pas de lecture Hwater si icompteur = 3 |
---|
132 | read(num_forc,*) Hwater |
---|
133 | read(num_forc,*) |
---|
134 | else |
---|
135 | read(num_forc,*) bidon |
---|
136 | read(num_forc,*) |
---|
137 | Hwater(:,:)=0. |
---|
138 | |
---|
139 | end if lect_Hwater |
---|
140 | |
---|
141 | ! lecture vitesses dans tous les cas. |
---|
142 | read(num_forc,*) Uxbar |
---|
143 | read(num_forc,*) |
---|
144 | |
---|
145 | read(num_forc,*) Uybar |
---|
146 | |
---|
147 | |
---|
148 | ! tableaux 3d |
---|
149 | |
---|
150 | read(num_forc,*) |
---|
151 | read(num_forc,*) |
---|
152 | |
---|
153 | read(num_forc,*,end=22) T |
---|
154 | 22 continue |
---|
155 | |
---|
156 | read(num_forc,*,end=23) |
---|
157 | read(num_forc,*,end=23) |
---|
158 | |
---|
159 | read(num_forc,*,end=23) Ux |
---|
160 | 23 continue |
---|
161 | |
---|
162 | read(num_forc,*,end=24) |
---|
163 | read(num_forc,*,end=24) |
---|
164 | |
---|
165 | read(num_forc,*,end=24) Uy |
---|
166 | 24 continue |
---|
167 | |
---|
168 | ! aurel pour faire des cptr avec les traceurs : |
---|
169 | read(num_forc,*,end=25) |
---|
170 | read(num_forc,*,end=25) |
---|
171 | read(num_forc,*,end=25) Xdep_out !xdepk ou xdep ? a voir |
---|
172 | 25 continue |
---|
173 | |
---|
174 | read(num_forc,*,end=26) |
---|
175 | read(num_forc,*,end=26) |
---|
176 | read(num_forc,*,end=26) Ydep_out |
---|
177 | 26 continue |
---|
178 | |
---|
179 | read(num_forc,*,end=27) |
---|
180 | read(num_forc,*,end=27) |
---|
181 | read(num_forc,*,end=27) tdep_out |
---|
182 | 27 continue |
---|
183 | tdep_out(:,:,:)=tdep_out(:,:,:)+time |
---|
184 | |
---|
185 | close(num_forc) |
---|
186 | |
---|
187 | else if ((index(reprcptr,'.nc').ne.0)) then ! Modif hassine pour la reprise d'un fichier .nc |
---|
188 | |
---|
189 | if (itracebug.eq.1) call tracebug(' Entree dans routine read_recovery : netcdf') |
---|
190 | |
---|
191 | |
---|
192 | if (.not.associated(liste_time)) then |
---|
193 | allocate(liste_time(1)) |
---|
194 | end if |
---|
195 | if (.not.associated(tab)) then |
---|
196 | allocate(tab(nx,ny),tab1(nx,ny,nz),tab1T(nx,ny,nz+nzm)) |
---|
197 | end if |
---|
198 | |
---|
199 | call Read_ncdf_var('temps',TRIM(reprcptr),liste_time) |
---|
200 | |
---|
201 | ! temps de la reprise |
---|
202 | !----------------------- |
---|
203 | if (tbegin.gt.1.d9) then ! si tbegin est tres grand , on reprend le temps du .nc |
---|
204 | tbegin=liste_time(1) |
---|
205 | time = tbegin |
---|
206 | else |
---|
207 | time = tbegin |
---|
208 | endif |
---|
209 | |
---|
210 | if (itracebug.eq.1) write(num_tracebug,*) tbegin, liste_time(1) |
---|
211 | |
---|
212 | |
---|
213 | ! TBEGIN=liste_time(1)-liste_time(1) |
---|
214 | |
---|
215 | call Read_ncdf_dim('x' ,TRIM(reprcptr),nxx) |
---|
216 | call Read_ncdf_dim('y' ,TRIM(reprcptr),nyy) |
---|
217 | call Read_ncdf_dim('z' ,TRIM(reprcptr),nzz) |
---|
218 | call Read_ncdf_dim('zm',TRIM(reprcptr),nzzm) |
---|
219 | |
---|
220 | if ((nxx.ne.nx).or.(nyy.ne.ny).or.(nzz.ne.nz) & |
---|
221 | .or.(nzzm-nzz.ne.nzm)) write(6,*) 'attention', & |
---|
222 | 'les tailles de tableaux ne sont pas compatibles' |
---|
223 | |
---|
224 | ! lecture des tableaux 2D |
---|
225 | ! ----------------------------- |
---|
226 | |
---|
227 | if (icpt.eq.1) then ! lecture topo seulement si icompteur = 1 |
---|
228 | call Read_ncdf_var('S',TRIM(reprcptr),tab) |
---|
229 | S(:,:)=tab(:,:) |
---|
230 | call Read_ncdf_var('H',TRIM(reprcptr),tab) |
---|
231 | H(:,:)=tab(:,:) |
---|
232 | call Read_ncdf_var('B',TRIM(reprcptr),tab) |
---|
233 | B(:,:)=tab(:,:) |
---|
234 | call Read_ncdf_var('BSOC',TRIM(reprcptr),tab) |
---|
235 | Bsoc(:,:)=tab(:,:) |
---|
236 | end if |
---|
237 | |
---|
238 | ! lecture temperature dans tous les cas |
---|
239 | call Read_ncdf_var('IBASE',TRIM(reprcptr),tab) |
---|
240 | Ibase(:,:)=tab(:,:) |
---|
241 | call Read_ncdf_var('BMELT',TRIM(reprcptr),tab) |
---|
242 | Bmelt(:,:)=tab(:,:) |
---|
243 | |
---|
244 | if (icpt.ne.3) then ! pas de lecture Hwater si icompteur = 3 |
---|
245 | call Read_ncdf_var('HWATER',TRIM(reprcptr),tab) |
---|
246 | Hwater(:,:)=tab(:,:) |
---|
247 | else |
---|
248 | Hwater(:,:)=0. |
---|
249 | end if |
---|
250 | |
---|
251 | ! lecture vitesses dans tous les cas. |
---|
252 | call Read_ncdf_var('UXBAR',TRIM(reprcptr),tab) |
---|
253 | UXBAR(:,:)=tab(:,:) |
---|
254 | call Read_ncdf_var('UYBAR',TRIM(reprcptr),tab) |
---|
255 | UYBAR(:,:)=tab(:,:) |
---|
256 | |
---|
257 | do i=1,nx-1 |
---|
258 | do j=1,ny-1 |
---|
259 | oldu(i,j) = sqrt(((uxbar(i,j)+uxbar(i+1,j))/2)**2 & |
---|
260 | + ((uybar(i,j)+uybar(i,j+1))/2.)**2) |
---|
261 | end do |
---|
262 | end do |
---|
263 | |
---|
264 | ! tableaux 3D |
---|
265 | ! ----------------------------- |
---|
266 | |
---|
267 | call Read_ncdf_var('T',TRIM(reprcptr),tab1T) |
---|
268 | T(:,:,:)=tab1T(:,:,:) |
---|
269 | |
---|
270 | call Read_ncdf_var('UX',TRIM(reprcptr),tab1) |
---|
271 | Ux(:,:,:)=tab1(:,:,:) |
---|
272 | |
---|
273 | call Read_ncdf_var('UY',TRIM(reprcptr),tab1) |
---|
274 | Uy(:,:,:)=tab1(:,:,:) |
---|
275 | |
---|
276 | if (itracebug.eq.1) call tracebug('avant lectures variables traceur') |
---|
277 | |
---|
278 | |
---|
279 | call Read_ncdf_var('XDEP',TRIM(reprcptr),tab1) !aurel neem |
---|
280 | Xdep_out(:,:,:)=tab1(:,:,:) |
---|
281 | |
---|
282 | call Read_ncdf_var('YDEP',TRIM(reprcptr),tab1) |
---|
283 | Ydep_out(:,:,:)=tab1(:,:,:) |
---|
284 | |
---|
285 | call Read_ncdf_var('TDEP',TRIM(reprcptr),tab1) |
---|
286 | Tdep_out(:,:,:)=tab1(:,:,:)+time |
---|
287 | |
---|
288 | end if ascii_nc |
---|
289 | |
---|
290 | |
---|
291 | do i=1,nx-1 |
---|
292 | do j=1,ny-1 |
---|
293 | oldu(i,j) = sqrt(((Uxbar(i,j)+Uxbar(i+1,j))/2)**2 & |
---|
294 | + ((Uybar(i,j)+Uybar(i,j+1))/2.)**2) |
---|
295 | end do |
---|
296 | end do |
---|
297 | |
---|
298 | |
---|
299 | ! moyenne de l'epaisseur et pentes |
---|
300 | do j=2,ny |
---|
301 | do i=2,nx |
---|
302 | Hmy(i,j)=(H(i,j)+H(i,j-1))/2. |
---|
303 | Hmx(i,j)=(H(i,j)+H(i-1,j))/2. |
---|
304 | Sdx(i,j)=(S(i,j)-S(i-1,j))/dx |
---|
305 | Sdy(i,j)=(S(i,j)-S(i,j-1))/dy |
---|
306 | end do |
---|
307 | end do |
---|
308 | |
---|
309 | do i=2,nx |
---|
310 | Hmx(i,1)=(H(i,1)+H(i-1,1))/2. |
---|
311 | end do |
---|
312 | do j=2,ny |
---|
313 | Hmy(1,j)=(H(1,j)+H(1,j-1))/2. |
---|
314 | end do |
---|
315 | |
---|
316 | ! vitesse a la base |
---|
317 | Uxflgz(:,:) = Ux(:,:,nz) |
---|
318 | Uyflgz(:,:) = Uy(:,:,nz) |
---|
319 | |
---|
320 | where (flot(:,:)) |
---|
321 | debug_3D(:,:,122) = T(:,:,nz)+273.15 |
---|
322 | elsewhere |
---|
323 | debug_3D(:,:,120) = T(:,:,nz)+273.15 |
---|
324 | endwhere |
---|
325 | |
---|
326 | end subroutine read_recovery |
---|
327 | |
---|
328 | |
---|
329 | !> SUBROUTINE: read_no_recovery |
---|
330 | !!Pas de reprise d'un fichier compteur |
---|
331 | !> |
---|
332 | subroutine read_no_recovery |
---|
333 | |
---|
334 | use module3D_phy,only: T,ghf |
---|
335 | use geography, only: nx,ny,nz,nzm |
---|
336 | use param_phy_mod,only:dzm,cm |
---|
337 | |
---|
338 | implicit none |
---|
339 | |
---|
340 | integer :: i,j,k ! loop integers |
---|
341 | |
---|
342 | !!!pas reprise : il faut initier les temperatures ds le socle |
---|
343 | |
---|
344 | ! temperature dans le socle lineaire avec le gradient geothermique |
---|
345 | do k=nz+1,nz+nzm |
---|
346 | do j=1,ny |
---|
347 | do i=1,nx |
---|
348 | T(i,j,k)=T(i,j,nz)-dzm*(k-nz)*ghf(i,j)/cm |
---|
349 | end do |
---|
350 | end do |
---|
351 | end do |
---|
352 | |
---|
353 | |
---|
354 | end subroutine read_no_recovery |
---|
355 | |
---|
356 | !---------------------------------------------------------- |
---|
357 | |
---|
358 | !---------------------------------------------------------- |
---|
359 | |
---|
360 | !> SUBROUTINE: out_recovery |
---|
361 | !! Sortie d'un fichier de reprise |
---|
362 | !! Attention : fichier en sortie .cptr ou fichier .nc |
---|
363 | !! used modules: - use module3D_phy |
---|
364 | !! - use netcdf |
---|
365 | !! - use io_netcdf_grisli |
---|
366 | !! - use util_recovery |
---|
367 | !! - use tracer_vars |
---|
368 | !> |
---|
369 | |
---|
370 | subroutine out_recovery(isortie) |
---|
371 | |
---|
372 | ! --------------------------------------------------------------- |
---|
373 | ! SORTIE COMPTEUR |
---|
374 | ! |
---|
375 | ! |
---|
376 | ! sortie qui stocke dans le fichier runname.cptr les variables |
---|
377 | ! S, H, B, HDOT, BDOT, et T plusieurs fois dans le programme |
---|
378 | ! (chaque nouvelle sortie ecrase la precedente) => |
---|
379 | ! on peut faire repartir le programme a partir de ce fichier, |
---|
380 | ! et donc du dernier pas de temps stocke |
---|
381 | ! |
---|
382 | ! Maintenant, il y a deux formats de fichiers possibles : |
---|
383 | ! ascii (isortie = 1) et netcdf (isortie = 2). |
---|
384 | ! --------------------------------------------------------------- |
---|
385 | |
---|
386 | use module3D_phy, only: num_file1,s,h,b,bsoc,ibase,bmelt,hwater,ux,uy,uxbar,uybar,t,time |
---|
387 | use runparam, only: itracebug |
---|
388 | use geography, only: nx,ny,nz,nzm,geoplace |
---|
389 | use netcdf, only: nf90_64bit_offset,nf90_create,nf90_close,nf90_write |
---|
390 | use io_netcdf_grisli, only: ncdf_type,write_ncdf_var |
---|
391 | use util_recovery,only: logic_out,time_out,testout_recovery |
---|
392 | use tracer_vars,only: xdep,ydep,tdep,xdep_out,ydep_out,tdep_out |
---|
393 | |
---|
394 | implicit none |
---|
395 | |
---|
396 | character(len=80) :: filin |
---|
397 | integer ncid,status |
---|
398 | real*8, dimension(:,:), pointer :: tab => null() !< tableau 2d real ecrit dans le fichier |
---|
399 | real*8, dimension(:,:,:), pointer :: tab1 => null() !< tableau 3d real |
---|
400 | real*8, dimension(:,:,:), pointer :: tab1T => null() !< tableau 3d real pour la temperature |
---|
401 | character(len=20),dimension(3) :: dimnames2d !< dimensions pour netcdf pour tableau 2d pour les noeud majeur |
---|
402 | character(len=20),dimension(4) :: dimnames3d !< pour 3d troisieme dim est nz |
---|
403 | character(len=20),dimension(4) :: dimnames3dT !< pour 3d troisieme dim est nz+nzm |
---|
404 | ! integer :: nrecs=1 !< compteur pour les enregistrements temps des variables |
---|
405 | ! integer :: idef=0 !< pour savoir si la variable a ete definie ou non |
---|
406 | integer :: isortie !< pour le choix de type de fichier en sortie |
---|
407 | |
---|
408 | if (itracebug.eq.1) call tracebug(' Entree dans routine out_recovery') |
---|
409 | |
---|
410 | call testout_recovery(filin) |
---|
411 | if (itracebug.eq.1) call tracebug(' Entree dans routine out_recovery 2') |
---|
412 | |
---|
413 | if (logic_out) then |
---|
414 | Select Case (Isortie) |
---|
415 | Case (1) ! sortie de type ascii |
---|
416 | |
---|
417 | filin=trim(filin)//'.cptr' |
---|
418 | open(num_file1,file=trim(filin)) |
---|
419 | |
---|
420 | write(num_file1,*) time_out(1), ' TIME ' |
---|
421 | write(num_file1,*) geoplace |
---|
422 | write(num_file1,*) nx*ny,nz,nzm,nx,ny |
---|
423 | write(num_file1,*) |
---|
424 | |
---|
425 | ! ecriture des tableaux 2D |
---|
426 | write(num_file1,*)' Tableaux 2D' |
---|
427 | write(num_file1,*)' Surface ' |
---|
428 | write(num_file1,*) S |
---|
429 | |
---|
430 | write(num_file1,*)' Epaisseur ' |
---|
431 | write(num_file1,*) H |
---|
432 | |
---|
433 | write(num_file1,*)' Base glace ' |
---|
434 | write(num_file1,*) B |
---|
435 | |
---|
436 | write(num_file1,*)' Socle ' |
---|
437 | write(num_file1,*) Bsoc |
---|
438 | |
---|
439 | write(num_file1,*)' Ibase ' |
---|
440 | write(num_file1,*) Ibase |
---|
441 | |
---|
442 | write(num_file1,*)' Bmelt ' |
---|
443 | write(num_file1,*) Bmelt |
---|
444 | |
---|
445 | write(num_file1,*)' Hwater ' |
---|
446 | write(num_file1,*) Hwater |
---|
447 | |
---|
448 | write(num_file1,*)' Uxbar ' |
---|
449 | write(num_file1,*) Uxbar |
---|
450 | |
---|
451 | write(num_file1,*)' Uybar ' |
---|
452 | write(num_file1,*) Uybar |
---|
453 | |
---|
454 | ! tableau 3D |
---|
455 | |
---|
456 | write(num_file1,*) |
---|
457 | write(num_file1,*) 'Temperature (3D y compris le socle)' |
---|
458 | write(num_file1,*) T(:,:,:) |
---|
459 | write(num_file1,*) |
---|
460 | |
---|
461 | write(num_file1,*) 'vitesse selon x' |
---|
462 | write(num_file1,*) Ux(:,:,:) |
---|
463 | write(num_file1,*) |
---|
464 | |
---|
465 | write(num_file1,*) 'vitesse selon y' |
---|
466 | write(num_file1,*) Uy(:,:,:) |
---|
467 | |
---|
468 | !aurel neem : traceurs dans le cptr |
---|
469 | |
---|
470 | write(num_file1,*) |
---|
471 | write(num_file1,*) 'origine de la glace en x' |
---|
472 | write(num_file1,*) xdep(:,:,:) |
---|
473 | |
---|
474 | write(num_file1,*) |
---|
475 | write(num_file1,*) 'origine de la glace en y' |
---|
476 | write(num_file1,*) ydep(:,:,:) |
---|
477 | |
---|
478 | write(num_file1,*) |
---|
479 | write(num_file1,*) 'age de la glace' |
---|
480 | tab1(:,:,:) = tdep(:,:,:) - time |
---|
481 | write(num_file1,*) tab1(:,:,:) ! note: on decide d'ecrire le temps 'relatif' |
---|
482 | |
---|
483 | close(num_file1) |
---|
484 | |
---|
485 | |
---|
486 | Case (2) ! sortie de type netcdf |
---|
487 | |
---|
488 | if (.not.associated(tab)) then |
---|
489 | allocate(tab(nx,ny),tab1(nx,ny,nz),tab1T(nx,ny,nz+nzm)) |
---|
490 | end if |
---|
491 | filin=trim(filin)//'.nc' |
---|
492 | |
---|
493 | ! creation du fichier |
---|
494 | !------------------------ |
---|
495 | |
---|
496 | if (ncdf_type.eq.32) then |
---|
497 | status = nf90_create(TRIM(filin),NF90_WRITE,ncid) ! en 32 bits |
---|
498 | else if (ncdf_type.eq.64) then |
---|
499 | status = nf90_create(trim(filin),and(nf90_write,nf90_64bit_offset),ncid) ! sur r2d2 |
---|
500 | else |
---|
501 | write(6,*)'pb de lecture de netcdf_type dans out_cptr:', ncdf_type |
---|
502 | endif |
---|
503 | |
---|
504 | status = nf90_close(ncid) ! fermeture |
---|
505 | |
---|
506 | |
---|
507 | call initfile(nx,ny,nz,nzm,dimnames2d,dimnames3d,dimnames3dT,filin) ! initialise le fichier nc |
---|
508 | |
---|
509 | call write_ncdf_var('temps','time',trim(filin),time_out,'double') ! ecrit le temps du snapshot |
---|
510 | |
---|
511 | ! variables 2D |
---|
512 | |
---|
513 | tab(:,:)= S(:,:) |
---|
514 | call write_ncdf_var('S',dimnames2d,trim(filin),tab,'double') ! ecrit S |
---|
515 | |
---|
516 | tab(:,:)= H(:,:) |
---|
517 | call write_ncdf_var('H',dimnames2d,trim(filin),tab,'double') |
---|
518 | |
---|
519 | tab(:,:)= B(:,:) |
---|
520 | call write_ncdf_var('B',dimnames2d,trim(filin),tab,'double') |
---|
521 | |
---|
522 | tab(:,:)= Bsoc(:,:) |
---|
523 | call write_ncdf_var('BSOC',dimnames2d,trim(filin),tab,'double') |
---|
524 | |
---|
525 | tab(:,:)= IBASE(:,:) |
---|
526 | call write_ncdf_var('IBASE',dimnames2d,trim(filin),tab,'double') |
---|
527 | |
---|
528 | tab(:,:)= Bmelt(:,:) |
---|
529 | call write_ncdf_var('BMELT',dimnames2d,trim(filin),tab,'double') |
---|
530 | |
---|
531 | tab(:,:)= Hwater(:,:) |
---|
532 | call write_ncdf_var('HWATER',dimnames2d,trim(filin),tab,'double') |
---|
533 | |
---|
534 | tab(:,:)= UXBAR(:,:) |
---|
535 | call write_ncdf_var('UXBAR',dimnames2d,trim(filin),tab,'double') |
---|
536 | |
---|
537 | tab(:,:)= UYBAR(:,:) |
---|
538 | call write_ncdf_var('UYBAR',dimnames2d,trim(filin),tab,'double') |
---|
539 | |
---|
540 | ! variables 3D |
---|
541 | |
---|
542 | tab1T(:,:,:)= T(:,:,:) |
---|
543 | call write_ncdf_var('T',dimnames3dT,trim(filin),tab1T,'double') ! ecrit T |
---|
544 | |
---|
545 | tab1(:,:,:)= UX(:,:,:) |
---|
546 | call write_ncdf_var('UX',dimnames3d,trim(filin),tab1,'double') |
---|
547 | |
---|
548 | tab1(:,:,:)= UY(:,:,:) |
---|
549 | call write_ncdf_var('UY',dimnames3d,trim(filin),tab1,'double') |
---|
550 | |
---|
551 | ! variables 3D pour les traceurs (Aurel NEEM) |
---|
552 | ! attention ce sont les tabelaux _out qui ont les bonnes dimensions. |
---|
553 | |
---|
554 | tab1(:,:,:)=XDEP_out(:,:,:) |
---|
555 | call write_ncdf_var('XDEP',dimnames3d,trim(filin),tab1,'double') !aurel neem |
---|
556 | |
---|
557 | tab1(:,:,:)=YDEP_out(:,:,:) |
---|
558 | call write_ncdf_var('YDEP',dimnames3d,trim(filin),tab1,'double') |
---|
559 | |
---|
560 | tab1(:,:,:)=TDEP_out(:,:,:) - time |
---|
561 | call write_ncdf_var('TDEP',dimnames3d,trim(filin),tab1,'double') |
---|
562 | |
---|
563 | |
---|
564 | End Select |
---|
565 | End if |
---|
566 | end subroutine out_recovery |
---|
567 | |
---|
568 | !> subroutine initfile |
---|
569 | !! Initialise the netcdf file |
---|
570 | !! Used module: - use netcdf |
---|
571 | !! - use io_netcdf_grisli |
---|
572 | !! @param nxx dimension along x |
---|
573 | !! @param nyy dimension along y |
---|
574 | !! @param nzz dimension along z |
---|
575 | !! @param nzmm dimension along z for T |
---|
576 | !! @param fil_sortie name of the file to initialise |
---|
577 | !! @param dimnames2d dimensions for netcdf |
---|
578 | !> |
---|
579 | subroutine initfile(nxx,nyy,nzz,nzmm,& |
---|
580 | dimnames2d,dimnames3d,dimnames3dT,file) |
---|
581 | |
---|
582 | use io_netcdf_grisli, only: write_ncdf_dim |
---|
583 | implicit none |
---|
584 | integer,intent(in) :: nxx !< dimension along x |
---|
585 | integer,intent(in) :: nyy !< dimension along y |
---|
586 | integer,intent(in) :: nzz !< dimension along z |
---|
587 | integer,intent(in) :: nzmm !< dimension along z in socle |
---|
588 | character(len=20),dimension(3) :: dimnames2d !< dimensions pour netcdf pour tableau 2d pour les noeud majeur |
---|
589 | character(len=20),dimension(4) :: dimnames3d !< pour 3d troisieme dim est nz |
---|
590 | character(len=20),dimension(4) :: dimnames3dT !< pour 3d troisieme dim est nz+nzm |
---|
591 | character(len=*) ::file !< name of the file to init |
---|
592 | ! initialisation |
---|
593 | call write_ncdf_dim('x',trim(file),nxx) !< dimensions des tableaux |
---|
594 | call write_ncdf_dim('y',trim(file),nyy) |
---|
595 | call write_ncdf_dim('z',trim(file),nzz) |
---|
596 | call write_ncdf_dim('zm',trim(file),nzz+nzmm) |
---|
597 | call write_ncdf_dim('time',trim(file),0) |
---|
598 | dimnames2d(1)='x' |
---|
599 | dimnames2d(2)='y' |
---|
600 | dimnames2d(3)='time' |
---|
601 | |
---|
602 | dimnames3d(1)='x' |
---|
603 | dimnames3d(2)='y' |
---|
604 | dimnames3d(3)='z' |
---|
605 | dimnames3d(4)='time' |
---|
606 | |
---|
607 | dimnames3dT(1)='x' |
---|
608 | dimnames3dT(2)='y' |
---|
609 | dimnames3dT(3)='zm' |
---|
610 | dimnames3dT(4)='time' |
---|
611 | |
---|
612 | end subroutine initfile |
---|
613 | |
---|
614 | |
---|
615 | !> SUBROUTINE: symetry_cptr |
---|
616 | !! Subroutine utilisee dans les experiences axysymetriques pour repartir d'un cptr dans lequel |
---|
617 | !! tous les champs sont axysmetriques ou symetriques par rapport a un axe |
---|
618 | !! Used module: - use module3D_phy |
---|
619 | !! @param iaxe symetrie par rapport a iaxe |
---|
620 | !! @param jaxe symetrie par rapport a jaxe |
---|
621 | !> |
---|
622 | |
---|
623 | subroutine symetry_cptr(jaxe) |
---|
624 | |
---|
625 | use module3D_phy,only: s,h,b,bsoc,ibase,bmelt,hwater,uxbar,uybar,t,ux,uy |
---|
626 | use geography,only: nx,ny |
---|
627 | use tracer_vars,only: xdep,ydep,tdep |
---|
628 | |
---|
629 | implicit none |
---|
630 | |
---|
631 | integer :: i,j |
---|
632 | integer :: jaxe |
---|
633 | integer :: jsym |
---|
634 | |
---|
635 | !symetrique par rapport à jaxe la référence est le bas. Pour l'instant seulement lui |
---|
636 | jsym=0 |
---|
637 | do j=jaxe+1,ny |
---|
638 | jsym=jsym+1 |
---|
639 | |
---|
640 | do i=1,nx |
---|
641 | S(i,j)=S(i,jaxe-jsym) |
---|
642 | H(i,j)=H(i,jaxe-jsym) |
---|
643 | B(i,j)=B(i,jaxe-jsym) |
---|
644 | Bsoc(i,j)=Bsoc(i,jaxe-jsym) |
---|
645 | Ibase(i,j)=Ibase(i,jaxe-jsym) |
---|
646 | Bmelt(i,j)=bmelt(i,jaxe-jsym) |
---|
647 | hwater(i,j)=hwater(i,jaxe-jsym) |
---|
648 | uxbar(i,j)=uxbar(i,jaxe-jsym) |
---|
649 | uybar(i,j)=-uybar(i,jaxe-jsym+1) ! attention different des autres |
---|
650 | ! a cause des staggered grids |
---|
651 | |
---|
652 | T(i,j,:)=T(i,jaxe-jsym,:) |
---|
653 | ux(i,j,:)=ux(i,jaxe-jsym,:) |
---|
654 | uy(i,j,:)=-uy(i,jaxe-jsym+1,:) ! different des autres |
---|
655 | |
---|
656 | XDEP(i,j,:)=XDEP(i,jaxe-jsym,:) |
---|
657 | YDEP(i,j,:)=YDEP(i,jaxe-jsym,:) |
---|
658 | TDEP(i,j,:)=TDEP(i,jaxe-jsym,:) |
---|
659 | |
---|
660 | end do |
---|
661 | end do |
---|
662 | return |
---|
663 | end subroutine symetry_cptr |
---|