1 |
guez |
69 |
module writehist_m |
2 |
guez |
56 |
|
3 |
guez |
69 |
implicit none |
4 |
guez |
56 |
|
5 |
guez |
69 |
contains |
6 |
guez |
56 |
|
7 |
guez |
138 |
subroutine writehist(time, vcov, ucov, teta, phi, masse, ps) |
8 |
guez |
56 |
|
9 |
guez |
69 |
! From writehist.F 1403 2010-07-01 09:02:53Z |
10 |
|
|
! Écriture du fichier histoire au format IOIPSL |
11 |
|
|
! Appels successifs des routines histwrite |
12 |
|
|
! L. Fairhead, LMD, 03/99 |
13 |
guez |
56 |
|
14 |
guez |
69 |
use dimens_m, only: nqmx, llm, jjm |
15 |
|
|
use com_io_dyn, only: histid, histvid, histuid |
16 |
|
|
use paramet_m, only: ip1jm, ip1jmp1, iip1, jjp1 |
17 |
|
|
use temps, only: itau_dyn |
18 |
|
|
use histwrite_m, only: histwrite |
19 |
|
|
use histsync_m, only: histsync |
20 |
|
|
use covnat_m, only: covnat |
21 |
guez |
56 |
|
22 |
guez |
69 |
! Entree: |
23 |
|
|
! time: temps de l'ecriture |
24 |
|
|
! vcov: vents v covariants |
25 |
|
|
! ucov: vents u covariants |
26 |
|
|
! teta: temperature potentielle |
27 |
|
|
! phi : geopotentiel instantane |
28 |
|
|
! masse: masse |
29 |
|
|
! ps :pression au sol |
30 |
guez |
56 |
|
31 |
guez |
69 |
! Arguments |
32 |
guez |
56 |
|
33 |
guez |
69 |
REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) |
34 |
|
|
REAL teta(ip1jmp1, llm), phi(ip1jmp1, llm) |
35 |
guez |
83 |
REAL, intent(in):: ps(ip1jmp1), masse(ip1jmp1, llm) |
36 |
guez |
69 |
integer time |
37 |
guez |
56 |
|
38 |
guez |
69 |
! This routine needs IOIPSL to work |
39 |
|
|
! Variables locales |
40 |
|
|
|
41 |
|
|
integer ndexu(ip1jmp1*llm), ndexv(ip1jm*llm), ndex2d(ip1jmp1) |
42 |
|
|
logical ok_sync |
43 |
|
|
integer itau_w |
44 |
|
|
REAL vnat(ip1jm, llm), unat(ip1jmp1, llm) |
45 |
|
|
|
46 |
|
|
!--------------------------------------------------------------------- |
47 |
|
|
|
48 |
|
|
! Initialisations |
49 |
|
|
|
50 |
|
|
ndexu = 0 |
51 |
|
|
ndexv = 0 |
52 |
|
|
ndex2d = 0 |
53 |
|
|
ok_sync =.TRUE. |
54 |
|
|
itau_w = itau_dyn + time |
55 |
|
|
! Passage aux composantes naturelles du vent |
56 |
|
|
call covnat(llm, ucov, vcov, unat, vnat) |
57 |
|
|
|
58 |
|
|
! Appels a histwrite pour l'ecriture des variables a sauvegarder |
59 |
|
|
|
60 |
|
|
! Vents U |
61 |
|
|
|
62 |
|
|
call histwrite(histuid, 'u', itau_w, unat) |
63 |
|
|
|
64 |
|
|
! Vents V |
65 |
|
|
|
66 |
|
|
call histwrite(histvid, 'v', itau_w, vnat) |
67 |
|
|
|
68 |
|
|
! Temperature potentielle |
69 |
|
|
|
70 |
|
|
call histwrite(histid, 'teta', itau_w, teta) |
71 |
|
|
|
72 |
|
|
! Geopotentiel |
73 |
|
|
|
74 |
|
|
call histwrite(histid, 'phi', itau_w, phi) |
75 |
|
|
|
76 |
|
|
! Masse |
77 |
|
|
|
78 |
|
|
call histwrite(histid, 'masse', itau_w, masse) |
79 |
|
|
|
80 |
|
|
! Pression au sol |
81 |
|
|
call histwrite(histid, 'ps', itau_w, ps) |
82 |
|
|
|
83 |
|
|
if (ok_sync) then |
84 |
|
|
call histsync(histid) |
85 |
|
|
call histsync(histvid) |
86 |
|
|
call histsync(histuid) |
87 |
|
|
endif |
88 |
|
|
|
89 |
|
|
end subroutine writehist |
90 |
|
|
|
91 |
|
|
end module writehist_m |