-
Notifications
You must be signed in to change notification settings - Fork 312
/
MEGANFactorsMod.F90
307 lines (249 loc) · 10.1 KB
/
MEGANFactorsMod.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
module MEGANFactorsMod
!-----------------------------------------------------------------------
! !DESCRIPTION:
! Manages input of MEGAN emissions factors from netCDF file
!
! !USES:
use shr_kind_mod, only : r8 => shr_kind_r8
use abortutils, only : endrun
use clm_varctl, only : iulog
use shr_log_mod, only : errMsg => shr_log_errMsg
!
implicit none
private
!
! !PUBLIC MEMBERS:
public :: megan_factors_init
public :: megan_factors_get
public :: comp_names
!
! !PUBLIC DATA:
real(r8), public, allocatable :: LDF(:) ! light dependent fraction
real(r8), public, allocatable :: Agro(:) ! growing leaf age factor
real(r8), public, allocatable :: Amat(:) ! mature leaf age factor
real(r8), public, allocatable :: Anew(:) ! new leaf age factor
real(r8), public, allocatable :: Aold(:) ! old leaf age factor
real(r8), public, allocatable :: betaT(:)! temperature factor
real(r8), public, allocatable :: ct1(:) ! temperature coefficient 1
real(r8), public, allocatable :: ct2(:) ! temperature coefficient 2
real(r8), public, allocatable :: Ceo(:) ! Eopt coefficient
!
! !PRIVATE MEMBERS:
integer :: npfts ! number of plant function types
!
type emis_eff_t
real(r8) , pointer :: eff (:) ! [real(r8) (:)] emissions efficiency factor
real(r8) :: wght ! molecular weight
integer :: class_num ! MEGAN class number
endtype emis_eff_t
!
type(emis_eff_t) , pointer :: comp_factors_table (:) ! [type(emis_eff_t) (:)] hash table of MEGAN factors (points to an array of pointers)
integer , pointer :: hash_table_indices (:) ! [integer (:)] pointer to hash table indices
integer, parameter :: tbl_hash_sz = 2**16 ! hash table size
!
character(len=40), allocatable :: comp_names(:) ! MEGAN compound names
real(r8), allocatable :: comp_molecwghts(:)! MEGAN compound molecular weights
character(len=*), parameter, private :: sourcefile = &
__FILE__
!-----------------------------------------------------------------------
contains
!-----------------------------------------------------------------------
subroutine megan_factors_get( comp_name, factors, class_n, molecwght )
!
! !DESCRIPTION:
! Method for getting MEGAN information for a named compound
!
! !ARGUMENTS:
character(len=*),intent(in) :: comp_name ! MEGAN compound name
real(r8), intent(out) :: factors(npfts) ! vegitation type factors for the compound of intrest
integer, intent(out) :: class_n ! MEGAN class number for the compound of intrest
real(r8), intent(out) :: molecwght ! molecular weight of the compound of intrest
!
! LOCAL VARS:
integer :: hashkey, ndx
character(len=120) :: errmes
!-----------------------------------------------------------------------
hashkey = gen_hashkey(comp_name)
ndx = hash_table_indices(hashkey)
if (ndx<1) then
errmes = 'megan_factors_get: '//trim(comp_name)//' compound not found in MEGAN table'
write(iulog,*) trim(errmes)
call endrun(msg=errMsg(sourcefile, __LINE__))
endif
factors(:) = comp_factors_table( ndx )%eff(:)
class_n = comp_factors_table( ndx )%class_num
molecwght = comp_factors_table( ndx )%wght
end subroutine megan_factors_get
!-----------------------------------------------------------------------
subroutine megan_factors_init( filename )
!
! !DESCRIPTION:
! Initializes the MEGAN factors using data from input file
!
! !USES:
use ncdio_pio, only : ncd_pio_openfile,ncd_inqdlen
use pio, only : pio_inq_varid,pio_get_var,file_desc_t,pio_closefile
use fileutils , only : getfil
use clm_varpar , only : mxpft
!
! !ARGUMENTS:
character(len=*),intent(in) :: filename ! MEGAN factors input file
!-----------------------------------------------------------------------
character(len=256) :: locfn ! local file name
type(file_desc_t) :: ncid ! netcdf id
integer :: start(2), count(2)
integer :: ierr, i, vid
integer :: dimid, n_comps, n_classes, n_patchs
integer :: class_ef_vid,comp_ef_vid,comp_name_vid,class_num_vid
integer :: comp_mw_vid
integer, allocatable :: class_nums(:)
real(r8),allocatable :: factors(:)
real(r8),allocatable :: comp_factors(:)
real(r8),allocatable :: class_factors(:)
allocate(comp_factors_table(150))
allocate(hash_table_indices(tbl_hash_sz))
call getfil(filename, locfn, 0)
call ncd_pio_openfile (ncid, trim(locfn), 0)
call ncd_inqdlen( ncid, dimid, n_comps, name='Comp_Num')
call ncd_inqdlen( ncid, dimid, n_classes, name='Class_Num')
call ncd_inqdlen( ncid, dimid, n_patchs, name='PFT_Num')
if ( n_patchs /= mxpft )then
call endrun(msg='PFT_Num does NOT equal mxpft: '//errMsg(sourcefile, __LINE__))
end if
npfts = n_patchs
ierr = pio_inq_varid(ncid,'Class_EF', class_ef_vid)
ierr = pio_inq_varid(ncid,'Comp_EF', comp_ef_vid)
ierr = pio_inq_varid(ncid,'Comp_Name',comp_name_vid)
ierr = pio_inq_varid(ncid,'Class_Num',class_num_vid)
ierr = pio_inq_varid(ncid,'Comp_MW', comp_mw_vid)
allocate( factors(n_patchs) )
allocate( comp_factors(n_patchs) )
allocate( class_factors(n_patchs) )
factors(1:n_patchs) = 0._r8
comp_factors(1:n_patchs) = 0._r8
class_factors(1:n_patchs) = 0._r8
allocate( comp_names(n_comps) )
allocate( comp_molecwghts(n_comps) )
allocate( class_nums(n_comps) )
ierr = pio_get_var( ncid, comp_name_vid, comp_names )
ierr = pio_get_var( ncid, comp_mw_vid, comp_molecwghts )
ierr = pio_get_var( ncid, class_num_vid, class_nums )
! set up hash table where data is stored
call bld_hash_table_indices( comp_names )
do i=1,n_comps
start=(/i,1/)
count=(/1,mxpft/)
ierr = pio_get_var( ncid, comp_ef_vid, start, count, comp_factors )
start=(/class_nums(i),1/)
ierr = pio_get_var( ncid, class_ef_vid, start, count, class_factors )
factors(:) = comp_factors(:)*class_factors(:)
call enter_hash_data( trim(comp_names(i)), factors, class_nums(i), comp_molecwghts(i) )
enddo
allocate( LDF(n_classes) )
allocate( Agro(n_classes) )
allocate( Amat(n_classes) )
allocate( Anew(n_classes) )
allocate( Aold(n_classes) )
allocate( betaT(n_classes) )
allocate( ct1(n_classes) )
allocate( ct2(n_classes) )
allocate( Ceo(n_classes) )
ierr = pio_inq_varid(ncid,'LDF', vid)
ierr = pio_get_var( ncid, vid, LDF )
ierr = pio_inq_varid(ncid,'Agro', vid)
ierr = pio_get_var( ncid, vid, Agro )
ierr = pio_inq_varid(ncid,'Amat', vid)
ierr = pio_get_var( ncid, vid, Amat )
ierr = pio_inq_varid(ncid,'Anew', vid)
ierr = pio_get_var( ncid, vid, Anew )
ierr = pio_inq_varid(ncid,'Aold', vid)
ierr = pio_get_var( ncid, vid, Aold )
ierr = pio_inq_varid(ncid,'betaT', vid)
ierr = pio_get_var( ncid, vid, betaT )
ierr = pio_inq_varid(ncid,'ct1', vid)
ierr = pio_get_var( ncid, vid, ct1 )
ierr = pio_inq_varid(ncid,'ct2', vid)
ierr = pio_get_var( ncid, vid, ct2 )
ierr = pio_inq_varid(ncid,'Ceo', vid)
ierr = pio_get_var( ncid, vid, Ceo )
call pio_closefile(ncid)
deallocate( class_nums, comp_factors,class_factors,factors )
endsubroutine megan_factors_init
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
! Private methods...
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
subroutine bld_hash_table_indices( names )
character(len=*),intent(in) :: names(:)
integer :: n, i, hashkey
hash_table_indices(:) = 0
n = size(names)
do i=1,n
hashkey = gen_hashkey(names(i))
hash_table_indices(hashkey) = i
enddo
endsubroutine bld_hash_table_indices
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
subroutine enter_hash_data( name, data, class_n, molec_wght )
character(len=*), intent(in) :: name
real(r8), intent(in) :: data(:)
integer, intent(in) :: class_n
real(r8), intent(in) :: molec_wght
integer :: hashkey, ndx
integer :: nfactors
hashkey = gen_hashkey(name)
nfactors = size(data)
ndx = hash_table_indices(hashkey)
allocate (comp_factors_table(ndx)%eff(nfactors))
comp_factors_table(ndx)%eff(:) = data(:)
comp_factors_table(ndx)%class_num = class_n
comp_factors_table(ndx)%wght = molec_wght
end subroutine enter_hash_data
!-----------------------------------------------------------------------
!from cam_history
!
! Purpose: Generate a hash key on the interval [0 .. tbl_hash_sz-1]
! given a character string.
!
! Algorithm is a variant of perl's internal hashing function.
!
!-----------------------------------------------------------------------
integer function gen_hashkey(string)
implicit none
!
! Arguments:
!
character(len=*), intent(in) :: string
!
! Local vars
!
integer :: hash
integer :: i
integer, parameter :: tbl_max_idx = 15 ! 2**N - 1
integer, parameter :: gen_hash_key_offset = z'000053db'
integer, dimension(0:tbl_max_idx) :: tbl_gen_hash_key = (/61,59,53,47,43,41,37,31,29,23,17,13,11,7,3,1/)
hash = gen_hash_key_offset
if ( len_trim(string) /= 19 ) then
!
! Process arbitrary string length.
!
do i = 1, len_trim(string)
hash = ieor(hash , (ichar(string(i:i)) * tbl_gen_hash_key(iand(i-1,tbl_max_idx))))
end do
else
!
! Special case string length = 19
!
do i = 1, tbl_max_idx+1
hash = ieor(hash , ichar(string(i:i)) * tbl_gen_hash_key(i-1))
end do
do i = tbl_max_idx+2, len_trim(string)
hash = ieor(hash , ichar(string(i:i)) * tbl_gen_hash_key(i-tbl_max_idx-2))
end do
end if
gen_hashkey = iand(hash, tbl_hash_sz-1)
return
end function gen_hashkey
end module MEGANFactorsMod