source: codes/icosagcm/trunk/src/disvert.f90 @ 78

Last change on this file since 78 was 75, checked in by ymipsl, 12 years ago

write file apbp.nc parameters for vertical discretisation

YM
File size: 3.9 KB
Line 
1MODULE disvert_mod
2  USE icosa
3  REAL(rstd), SAVE, POINTER :: ap(:)
4  REAL(rstd), SAVE, POINTER :: bp(:)
5  REAL(rstd), SAVE, POINTER :: presnivs(:)
6
7CONTAINS
8
9  SUBROUTINE init_disvert
10  USE disvert_std_mod, ONLY: ap_std=>ap, bp_std=>bp, presnivs_std=>presnivs, init_disvert_std=>init_disvert
11  USE disvert_ncar_mod, ONLY: ap_ncar=>ap, bp_ncar=>bp, presnivs_ncar=>presnivs, init_disvert_ncar=>init_disvert
12  USE disvert_ncarl30_mod, ONLY: ap_ncarl30=>ap, bp_ncarl30=>bp, presnivs_ncarl30=>presnivs, init_disvert_ncarl30=>init_disvert
13  USE icosa
14  IMPLICIT NONE
15    CHARACTER(LEN=255) :: disvert_type = 'std'
16   
17    CALL getin("disvert",disvert_type)
18   
19    SELECT CASE (TRIM(disvert_type))
20      CASE('std')
21   
22        CALL init_disvert_std
23        ap=>ap_std
24        bp=>bp_std
25        presnivs=>presnivs_std
26     
27      CASE ('ncar')
28
29        CALL init_disvert_ncar
30        ap=>ap_ncar
31        bp=>bp_ncar
32        presnivs=>presnivs_ncar
33
34      CASE ('ncarl30')
35
36        CALL init_disvert_ncarl30
37        ap=>ap_ncarl30
38        bp=>bp_ncarl30
39        presnivs=>presnivs_ncarl30
40       
41      CASE default
42        PRINT*,'Bad selector for variable disvert : <', TRIM(disvert_type),"> options are <std>, <ncar>, <ncarl30>" 
43        STOP
44       
45    END SELECT
46
47  END SUBROUTINE init_disvert 
48 
49  SUBROUTINE write_apbp
50  USE icosa
51  USE netcdf_mod
52  IMPLICIT NONE
53    REAL(rstd) :: val(llm)
54    INTEGER :: status
55    INTEGER :: lev,ilev
56    INTEGER :: ncid,levid,ilevid,hyaiid,hybiid,hyamid,hybmid,P0id
57    INTEGER :: l
58   
59    status = NF90_CREATE('apbp.nc', NF90_CLOBBER, ncid)
60    status = NF90_DEF_DIM(ncid,'lev',llm,lev)
61    status = NF90_DEF_DIM(ncid,'ilev',llm+1,ilev)
62   
63    status = NF90_DEF_VAR(ncid,'lev',NF90_DOUBLE,(/ lev /),levid)
64    status = NF90_PUT_ATT(ncid,levid,"long_name","hybrid level at midpoints (1000*(A+B))")
65    status = NF90_PUT_ATT(ncid,levid,"units","Pa")
66    status = NF90_PUT_ATT(ncid,levid,"positive","up")
67    status = NF90_PUT_ATT(ncid,levid,"standard_name","atmosphere_hybrid_sigma_pressure_coordinate")
68    status = NF90_PUT_ATT(ncid,levid,"formula_terms","a: hyam b: hybm p0: P0 ps: PS")
69
70    status = NF90_DEF_VAR(ncid,'ilev',NF90_DOUBLE,(/ ilev /),ilevid)
71    status = NF90_PUT_ATT(ncid,ilevid,"long_name","hybrid level at interfaces (1000*(A+B))")
72    status = NF90_PUT_ATT(ncid,ilevid,"units","Pa")
73    status = NF90_PUT_ATT(ncid,ilevid,"positive","up")
74    status = NF90_PUT_ATT(ncid,ilevid,"standard_name","atmosphere_hybrid_sigma_pressure_coordinate")
75    status = NF90_PUT_ATT(ncid,ilevid,"formula_terms","a: hyai b: hybi p0: P0 ps: PS")
76   
77    status = NF90_DEF_VAR(ncid,'hyai',NF90_DOUBLE,(/ ilev /),hyaiid)
78    status = NF90_PUT_ATT(ncid,hyaiid,"long_name","hybrid A coefficient at layer interfaces")
79
80    status = NF90_DEF_VAR(ncid,'hybi',NF90_DOUBLE,(/ ilev /),hybiid)
81    status = NF90_PUT_ATT(ncid,hybiid,"long_name","hybrid B coefficient at layer interfaces")
82
83    status = NF90_DEF_VAR(ncid,'hyam',NF90_DOUBLE,(/ lev /),hyamid)
84    status = NF90_PUT_ATT(ncid,hyamid,"long_name","hybrid A coefficient at midpoint interfaces")
85
86    status = NF90_DEF_VAR(ncid,'hybm',NF90_DOUBLE,(/ lev /),hybmid)
87    status = NF90_PUT_ATT(ncid,hybmid,"long_name","hybrid B coefficient at midpoint interfaces")
88   
89    status = NF90_DEF_VAR(ncid,'P0',NF90_DOUBLE,varid=P0id)
90
91    status = NF90_ENDDEF(ncid)   
92   
93    PRINT*,ap
94    PRINT*,bp
95    status=NF90_PUT_VAR(ncid,ilevid, ap(:)+bp(:)*Preff)
96   
97    DO l=1,llm
98      val(l)= 0.5*(ap(l+1)+ap(l)+Preff*(bp(l)+bp(l+1)))
99    ENDDO
100   
101    status=NF90_PUT_VAR(ncid,levid, val)
102
103    status=NF90_PUT_VAR(ncid,hyaiid, ap(:)/Preff)
104    status=NF90_PUT_VAR(ncid,hybiid, bp(:)/Preff)
105   
106    DO l=1,llm
107      val(l)= 0.5*(ap(l+1)+ap(l))
108    ENDDO
109    status=NF90_PUT_VAR(ncid,hyamid, val(:)/Preff)
110   
111     DO l=1,llm
112      val(l)= 0.5*(bp(l+1)+bp(l))
113    ENDDO
114    status=NF90_PUT_VAR(ncid,hybmid, val(:))
115 
116    status=NF90_PUT_VAR(ncid,P0id, Preff)
117   
118   status=NF90_CLOSE(ncid)
119  END SUBROUTINE write_apbp
120
121END MODULE disvert_mod
Note: See TracBrowser for help on using the repository browser.