/[lmdze]/trunk/Sources/phylmd/CV30_routines/cv30_uncompress.f
ViewVC logotype

Annotation of /trunk/Sources/phylmd/CV30_routines/cv30_uncompress.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 201 - (hide annotations)
Mon Jun 6 17:42:15 2016 UTC (7 years, 11 months ago) by guez
File size: 2856 byte(s)
Removed intermediary objects of cv_thermo_m, access suphec_m
directly. Procedure cv_thermo disappeared, all objects are named
constants.

In cv_driver and below, limited extents of arrays to what is needed.

lv, cpn and th in cv30_compress were set at level nl + 1 but lv1, cpn1
and th1 are not defined at this level. This did not lead to an error
because values at nl + 1 were not used.

Removed test on ok_sync in phystokenc because it is not read at run
time. Printing min and max of output NetCDF variables is heavy and
archaic.

Used histwrite_phy in phytrac.

1 guez 185 module cv30_uncompress_m
2 guez 47
3     implicit none
4    
5 guez 97 contains
6 guez 47
7 guez 185 SUBROUTINE cv30_uncompress(idcum, iflag, precip, VPrecip, sig, w0, ft, fq, &
8 guez 189 fu, fv, inb, Ma, upwd, dnwd, dnwd0, qcondc, cape, da, phi, mp, iflag1, &
9     precip1, VPrecip1, sig1, w01, ft1, fq1, fu1, fv1, inb1, Ma1, upwd1, &
10     dnwd1, dnwd01, qcondc1, cape1, da1, phi1, mp1)
11 guez 47
12 guez 196 ! UNCOMPRESS THE FIELDS
13    
14 guez 185 USE cv30_param_m, ONLY: nl
15 guez 103 use dimphy, only: klon, klev
16 guez 47
17 guez 103 integer, intent(in):: idcum(:) ! (ncum)
18 guez 196 integer, intent(in):: iflag(:) ! (ncum)
19 guez 103 real, intent(in):: precip(klon)
20     real, intent(in):: VPrecip(klon, klev+1)
21     real, intent(in):: sig(klon, klev), w0(klon, klev)
22     real, intent(in), dimension(klon, klev):: ft, fq, fu, fv
23 guez 201 integer, intent(in):: inb(:) ! (ncum)
24 guez 103 real, intent(in):: Ma(klon, klev)
25     real, intent(in):: upwd(klon, klev), dnwd(klon, klev), dnwd0(klon, klev)
26     real, intent(in):: qcondc(klon, klev)
27 guez 189 real, intent(in):: cape(klon)
28 guez 201 real, intent(in):: da(:, :) ! (klon, klev)
29     real, intent(in):: phi(:, :, :) ! (klon, klev, klev)
30     real, intent(in):: mp(:, :) ! (ncum, nl)
31 guez 47
32 guez 97 ! outputs:
33 guez 196 integer, intent(out):: iflag1(:) ! (klon)
34 guez 103 real precip1(klon)
35     real VPrecip1(klon, klev+1)
36     real sig1(klon, klev), w01(klon, klev)
37     real ft1(klon, klev), fq1(klon, klev), fu1(klon, klev), fv1(klon, klev)
38     integer, intent(inout):: inb1(klon)
39     real Ma1(klon, klev)
40     real upwd1(klon, klev), dnwd1(klon, klev), dnwd01(klon, klev)
41     real qcondc1(klon, klev)
42 guez 189 real cape1(klon)
43 guez 103 real, intent(inout):: da1(klon, klev), phi1(klon, klev, klev)
44     real, intent(inout):: mp1(klon, klev)
45 guez 47
46 guez 180 ! Local:
47 guez 103 integer ncum, i, k, j
48 guez 47
49 guez 97 !-------------------------------------------------------------------
50 guez 47
51 guez 103 ncum = size(idcum)
52 guez 196 iflag1 = 42 ! for non convective points
53 guez 103
54 guez 97 do i=1, ncum
55     precip1(idcum(i))=precip(i)
56     iflag1(idcum(i))=iflag(i)
57     inb1(idcum(i))=inb(i)
58     cape1(idcum(i))=cape(i)
59     end do
60 guez 47
61 guez 97 do k=1, nl
62     do i=1, ncum
63     VPrecip1(idcum(i), k)=VPrecip(i, k)
64     sig1(idcum(i), k)=sig(i, k)
65     w01(idcum(i), k)=w0(i, k)
66     ft1(idcum(i), k)=ft(i, k)
67     fq1(idcum(i), k)=fq(i, k)
68     fu1(idcum(i), k)=fu(i, k)
69     fv1(idcum(i), k)=fv(i, k)
70     Ma1(idcum(i), k)=Ma(i, k)
71     upwd1(idcum(i), k)=upwd(i, k)
72     dnwd1(idcum(i), k)=dnwd(i, k)
73     dnwd01(idcum(i), k)=dnwd0(i, k)
74     qcondc1(idcum(i), k)=qcondc(i, k)
75     da1(idcum(i), k)=da(i, k)
76     mp1(idcum(i), k)=mp(i, k)
77     end do
78     end do
79    
80     do i=1, ncum
81 guez 103 sig1(idcum(i), klev)=sig(i, klev)
82 guez 97 end do
83    
84 guez 103 do j=1, klev
85     do k=1, klev
86 guez 97 do i=1, ncum
87     phi1(idcum(i), k, j)=phi(i, k, j)
88     end do
89     end do
90     end do
91    
92 guez 185 end SUBROUTINE cv30_uncompress
93 guez 97
94 guez 185 end module cv30_uncompress_m

  ViewVC Help
Powered by ViewVC 1.1.21