1 function [E_CO2, E_CO, E_CH4, E_PM25] = FEMIS(fgi, burn_rate, U_F, V_F, RH4, mode)
2 %FEMIS function computing emission rates of CO2, CO, CH4 and particlate matter
4 % fgi - fuel load [kg/m2]
5 % burn_rate [kg/s] - mass of the fuel burnt in unit time.
6 % U_F, V_F [m/s] mid-flame wind speed components (used for computation of the WS_F)
7 % RH4 [%] -relative humidity 4 hours before (?)
8 % mode - 0 - default, computes consumption rates according to FEMIS
9 % allowing for the total consumption greater than the combustion rate
10 % - 1 - limits maximum consumption rate to the actual burn rate
11 % - 2 - computes emissions straight from burnt rate ignoring smoldering adjustment
13 %all the constants and equations come from the following US Forest Service
14 %technical report: Fire Emission Production Simulator (FEPS) User's Guide
15 %Version 1.0, Gary K. Anderson, David V. Sandberg and Robert A. Norheim
17 %available from: www.fs.fed.us/pnw/fera/feps/FEPS_users_guide.pdf
19 %Adam Kochanski 05.26.2012
23 %emission factor constants according to FEPS user's manual p. 63. [g/kg]
29 %emission factor coefficients according to FEPS user's manual p. 63. [g/kg]
35 %flam involvment sensitivity coeffcient (eq 5, p. 53)
38 % combustion effciency factors according to FEPS user's manual p. 62 [-]:
39 k_f=0.9; %flame phase (f)
40 k_sts=0.76; %short-term smoldering (sts)
41 k_lts=0.76; %long-term smoldering (lts)
43 %flaming phase consumption coeffcients from to FEPS user's manual p. 53 [-]
45 k_cbg=0.2; % not used since Cbg=0
47 % Duff bulk density [tons/acre-inch] p. 54
49 % STS depth sensitivity to consumption [tons/acre-inch] p.55
52 %flaming phase resindence time coeficients (eq. 8, p. 54)
56 %short-term moldering residence time coeffcients (eq.13 p. 55)
60 %flame phase diffusivity dependency
63 %short-term smoldering diffusivity dependency
66 %consumption threshold for flame involvment [tons/acre],
67 % 1 t/acre = 0.224 kg/m2
69 Cti=Cti * 0.224; % [kg/m2]
71 %wind speed benchmark for smoldering 1 mph = 0.44704 m/s
73 u_b= u_b * 0.44704; % [m/s]
75 %benchmark relative humidity deficit [%]
78 %Molecular weights [g/mol]
85 % by default computation from FEPS allowing greater and smaller consumption than the burn_rate
87 %==========================================================================
88 % computing flaming phase consumption, assuming no below ground consumption
89 % Cbg=0; and the above ground consumption (canopy, shrub, grass, wood
90 % litter) Cag equal to fuel combustion rate. eq. 6 p. 53
92 % flaming phase involvment
93 % flaming phase involvment
94 Cag=fgi; % we assume above the ground consuption equal to the fuel load [kg/m2]
95 %Cbg=0; % we assume no below ground consumption
97 Inv_f=100*(1-k_agi*exp(-Cag/Cti))
99 %==========================================================================
100 % computing short-term smoldering phase consumption [kg/m2]
103 % computing flaming phase consumption rate [kg/s]
104 CR_f= burn_rate * Inv_f/100 * k_cag
106 % computing smoldering adjustment from eq. 26, p 59
107 WS_F = sqrt (U_F^2 + V_F^2); % wind speed at midflame height [m/s]
108 %U_F midflame X wind speed component
109 %V_F midflame Y wind speed component
110 SA=sqrt(WS_F/u_b) * ((100-RH4)/RHb)
112 %computing sort-term smoldering consumption rate [kg/s]
113 CR_sts= SA* burn_rate * Inv_sts/100 * k_cag
115 % we assume no long term smoldering (>2h) since in most cases withing
116 % 2h all fuel within a grid cell will burn, and we have no data about the
120 %computing total consumption rate according for FEPS (eq. 30, p.61):
121 CR_total=CR_f + CR_sts + CR_lts; %[kg/s]
123 if (mode==1) % limit the consumption rate to the burn rate
124 if (CR_total>burn_rate)
127 elseif (mode==2) % compute the consumption rate straight from burn rate
132 %computing combustion efficiency:
133 CEff=(k_f*CR_f + k_sts*CR_sts + k_lts*CR_lts)/(CR_f + CR_sts + CR_lts + 1)
135 ER_CO2= (Kef_CO2 - EF_CO2*CEff)*CR_total/1000 %[kg/s]
136 ER_CO=(Kef_CO - EF_CO*CEff)*CR_total/1000 %[kg/s]
137 ER_CH4=(Kef_CH4 - EF_CH4*CEff)*CR_total/1000 %[kg/s]
138 ER_PM25=(Kef_PM25 - EF_PM25*CEff)*CR_total/1000 %[kg/s]
140 E_C02 = ER_CO2/(MW_CO2*3600) % CO2 emissions in mol/hr
141 E_C0 = ER_CO /(MW_CO *3600) % CO emissions in mol/hr
142 E_CH4 = ER_CH4/(MW_CH4*3600) % emissions in mol/hr
143 E_PM25 = ER_PM25*1000 % PM25 emissions in ug/s