UBC Faculty Research and Publications

Evapotranspiration-Interception Model for Urban Areas. Grimmond, C. S. B.; Oke, Timothy R. 1991

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


52383-Oke_AGU_1991WR00557.pdf [ 1.45MB ]
JSON: 52383-1.0041943.json
JSON-LD: 52383-1.0041943-ld.json
RDF/XML (Pretty): 52383-1.0041943-rdf.xml
RDF/JSON: 52383-1.0041943-rdf.json
Turtle: 52383-1.0041943-turtle.txt
N-Triples: 52383-1.0041943-rdf-ntriples.txt
Original Record: 52383-1.0041943-source.json
Full Text

Full Text

WATER RESOURCES RESEARCH, VOL. 27, NO. 7, PAGES 1739-1755, JULY 1991 An Evapotranspiration-Interception Model for Urban Areas C. S. B. GRIMMOND • AND T. R. 0KE Atmospheric Science Programtne, Departmen! of Geography, University qf British Columbia, Vancouver, Canada A model to calculate evapotranspiration from urban areas over a wide range of meteorological conditions is presented. An evapotranspiration-interception approach is used because it is necessary to cope with the changing water availability on the surface, during and following rainfall or irrigation. The model is applicable to areas ranging from the size of city blocks to land use zones and time periods of one hour and longer. The modeled evaporation is compared with that from micrometeorological measurements conducted from January to June 1987 in a suburb of Vancouver, British Columbia, Canada. The results show that this approach to modeling urban evapotranspiration provides realistic hourly and daily estimates of the areally averaged latent heat flux and surface water state. 1. INTRODUCTION Evapotranspiration (E), the mass equivalent of the latent heat energy flux, is the term that links the energy and the water balances. For urban areas the water balance is written P + I + F = E + r + ZXS + •XA mm hr -• (1) where P is precipitation, I is piped water supply, F is water released due to anthropogenic activities, r is runoff, •XS is the change in water storage for the period of interest, and •XA is net moisture advection; the urban energy balance is written Q* + QF = Qœ + QH + zXQs + •QA W m -2 (2) where Q* is net all-wave radiation, Qe is the energy released due to anthropogenic activities, Q• is the sensible heat flux, •Qs is the storage heat flux, and ZXQA is the net heat advection; Qe = L•,E where L,, is the latent heat of vaporization. These balances apply to the top of a volume which extends to sufficient depth that vertical heat and water exchange is negligible. Fluxes are determined per unit area of the top of the volume. A fuller explanation of this concept is given by Oke [1988]. Evaporation from urban areas has been assumed to be considerably ess than that from neighboring rural areas because of the supposed contrast between the hydrologic properties ofbuilding materials and vegetation-covered soils [e.g., Chandler, 1976]. Recent work by Kalanda et al. [1980], Cleugh and Oke [1986], and Grimmond and Oke [1986] shows that evaporation has a larger magnitude than originally thought, for example, in a temperate city (Vancou- ver, British Columbia (B.C.)) evaporation has been shown to constitute 38% of the losses of the annual external water balance and 81% of the losses of the summer water balance [Grirnmond and Oke, 1986]. Urban evaporation has received little attention in the fields of both urban climatology and urban hydrology. Indeed, it hardly rates mention in textbooks on the subject [e.g., Lazaro, 1979; Landsberg, 1981; Hall, 1984]. In most urban 1Now atClimate and Meteorology Program, Department of Geography, Indiana University, Bloomington. Copyright 1991 by the American Geophysical Union. Paper number 91WR00557. 0043-1397/91/91 WR-00557505.00 hydrologic work, evapotranspiration is either ignored or dealt with via evaporation pan data [e.g., Alley et al., 1980; Wenzell and Voorhees, 1980; Smart et al., 1988] so that van den Ven [ 1988] concludes that the weakest point in the study of the urban water balance remains the estimation of evapo- transpiration. A model for use in conjunction with available runoff models to determine continuous hourly evapotranspiration would be useful to model water quantity, water quality, and urban climate. For water quantity it is of interest to be able to predict the surface water storage status. A recent rend in urban runoff work favors continuous simulation rather than modeling the runoff from individual rainfall events. Be- tween-event simulation makes it necessary to account not only for the runoff but also for the other components of the water balance, including evaporation. For water quality modeling it is necessary to know the volumes of water present in order to calculate pollutant concentrations. In urban climatology there is a greater awareness of the evap- orative term because of its role in the amelioration of thermal climate, even so there are only a few studies. Models to determine the evapotranspiration from urban areas can be divided into two groups. First, there are those which consider Q• as one flux in the total surface energy balance. These are one-, two-, and quasi-three-dimensional numerical boundary layer models. The convective Q•v and Qe fluxes are often partitioned using a surface moisture availability factor or an assumed Bowen ratio (/3). After attempting to validate three such urban climate models against measured energy balance fluxes (for dry summertime conditions), Ross and Oke [!988] conclude that Qœ is the poorest of the fluxes imulated. The second group includes simple statistical algorithms and physically based models utilizing hydrometeorological d ta as input. The objective of this paper is to present a physically based model to estimate hourly, and longer-period, evapotranspi- ration associated with urban units ranging in size from city blocks to land use zones. 2. THE MODEL If evapotranspiration s to be calculated under diverse meteorological conditions, it is necessary for a model to cope with changing water availability on the surface during and following rainfall. Further, since the physical structure of the city has been compared to that of a forest with its 1739 1740 GRIMMOND AND OKE: EVAPOTRANSPIRAT!O N~INTERCEPTION MODEL ENERGY WATER Net Sensible radiation heat I '" Fig. 1. Schematic of the conceptual framework of the model. canopy, an evapotranspiration-interception approach similar to that used for forests is adopted. The urban canopy layer consists of the air contained between the buildings and trees. We propose an appropriate framework for urban evapotrans- piration to be the Penman-Monteith-Rutter-Shuttleworth evapotranspiration-interception model which was originally developed for forests [Penman, 1948; Monteith, 1965; Rutter et al., 1971; Shuttleworth, 1978]. It is the most successful, rigorous, robust, physically based evapotranspiration- interception model currently available [Gash, 1979; Shuttle- worth, 1983, 1989]. In essence, the model calculates a running water balance of the "canopy" and "trunks." It requires inputs of hourly rainfall and meteorological param- eters and provides outputs of throughfall, stemflow, and evaporation of intercepted water. Evaporation is calculated from the Penman-Monteith equation. The model allows for a continuous treatment of surface resistance for the transition between wet and dry surfaces. The conceptual framework of this model, as adapted for urban areas, is shown in Figure 1. The following is an outline of the model and the adaptations incorporated so that it can be applied to urban areas. The process of evapotranspiration from surfaces which have some vegetative cover can be considered to consist of three stages: (1) evaporation from a totally wet surface (due to intercepted water), (2) evapotranspiration from a partially wet surface, and (3) transpiration from a dry surface (water loss from outwardly dry vegetation via the stoma). To calculate E at any time it is necessary to know the surface water state to determine which stage is applicable. The Rutter model can be used to calculate the state of the running water balance of the canopy (C), summarized by the conti- nuity equation dC =P(1- TF-Ts)-D-E mmhr -1 (3) 'dr where t is time, TF is the throughfall or proportion of rain falling directly to the ground through the gaps in the canopy, T s is the proportion of rain diverted to the "trunks," and D is the drainage rate. In the original model the canopy is represented as a single-layer store of moisture filled by incoming rainfall. It is capable of holding a preset amount 0f water S and drains at a rate D expressed as an empirical function in proportion to the amount of water in the store. Throughfall represents water not available for evaporation from the canopy. Subsidiary evaporation may arise from trunks wetted by stemflow. Adapting this to the urban environment, the canopy is treated as a single-layer moisture store, with parallel stores which allow for the range of different surface types found in urban areas. As a first approximation, the urban system consists of two surface types: built/impervious and vege- tated. There is little relevance to the concept of throughfall for impervious surfaces (buildings, paved road, paths), and since urban trees usually have vegetation (grass) below, the water which falls between them is caught by the vegetation below. Therefore throughfall can be set to zero or a very small number for both the vegetated and built surface types. In urban areas the subsidiary loss, due to "stemflow", from wetted trunks can be interpreted as that from walls of buildings and from trunks of the tall vegetation. Generally, walls are not wetted by stemflow from the roof because of roof guttering. It is possible for them to be wetted directly b driving rain, but they are not extensively wetted under normal rain [Lacy, 1977]. Further, since the surface area of actual tree trunks in the urban system is not large [Grim. mond, 1988], urban stemflow is also set to zero. The piped water supply (I) is an additional source of water to wet he urban canopy. Therefore by neglecting TF and Ts but adding I the model can be reformulated for urban areas to read dCi - (P+Ii)-Di-E (4/ dt where subscript  is the ith surface type. In this study, six surface types are used to characterize th  study area: paved. built, coniferous, deciduous, irrigated grass, and unirrigated grass. A running water balance is maintained for each surface, and the evaporative flux calculated foran individual GRIMMOND AND OKE: EVAPOTRANSPIRATION-INTERCEPTION MODEL 1741 TABLE 1. Hourly Input Data Requirements and Outputs Generated by the Model , Units Input Rainfall (P) mm hr -1 Net radiation (Q*) W m -2 Wind speed (u) m s -1 Wind direction (q,) deg, rad Standard eviation of q>(c%) deg, rad Air temperature (Ta) øC Wet-bulb temperature (Tw) øC or relative humidity % Pressure (p) Pa Sensible heat flux (QH) or W m 2 temperature difference (ATa) External water use (wu) mm hr -1 Leaf area index (LAD Soil moisture deficit (80) mm Storage heat flux (AQs) W m -2 Anthropogenic heat flux (QF) W m -2 Output Latent heat flux (Qœ) W m -2 Evapotranspiration (E) mm hr-1 Cumulative evapotranspiration mm Drainage (D) mm hr - 1 Surface state* mm Surface storage change mm hr -1 -1 Aerodynamic resistance (ra) s m -1 Surface resistances (rs) s m *For total area and each surface type. hour is weighted for the abundance of the i types occurring in the source area of the measurements (see section 3.1). The change in storage for the whole system is determined from the summation of all the surface types. If it is assumed that the components of the urban surface are exposed to the same local scale climate, then a single- source (layer) model can be used to determine the latent heat flux (QE). The best tested approach is the Monteith [1965] version of the Penman equation [Stewart, !984; Dolman et al., 1988]. In the urban case it is modified to read s(Q* + QF- AQs) + (CaV)/ra Q•r = s + T(1 + rs/ra) where s is the slope of the saturation vapor pressure versus temperature curve (Pa øC-•), •, is the psychrometric " on- stant" (Pa øC-•), Ca is the heat capacity of the air (J m -3 øC-i), Vis the vapor pressure d ficit of the air (Pa), and ra and rs are the aerodynamic and surface resistances ( m-l). In this form the resistances for sensible and latent heat are equal nd related to the equivalent resistance for momentum (often called the aerodynamic resistance ra). The anthropo- genic heat flux (QF) is included as an additional source term. The quation is applicable.to dry surfaces ( tage 3); when the surface is completely wet (stage 1), the rs term is set o zero. The overlap between the evaporation f intercepted water and transpiration s poorly defined in the Rutter model [Shuttleworth, 1983]. To link the two stages, a physically continuous transition between wet and dry canopies is required. Shuttleworth [1978] proposed such a model in which r s is replaced bya redefined surface r sistance r ss, equal to rs in dry conditions and set to zero in wet condi- tions, with asmooth transition between the two depending on the fractional surface wetness. Shuttleworth's theoretical analysis has been tested over tall vegetation [Shuttleworth, 1978] and takes the form ß • rss = rb(s/• + 1) + 'r s + rt,(s/•/ + 1) rt,(s/•/ + 1) where W=I C_>S R-1 W= C<S R- s/c (rs/ra)(ra- rb) rs + r,(s/7 + 1) where rb = 1.1u/-1 + 5.6u. •/3 is the mean boundary layer resistance [Shuttleworth, 1983], and u. is frictional velocity (m s -1) (see section 3.4. for method of calculation). The complete urban evapotranspiration-interception model uses submodels (see section 3) to identify and de- scribe turbulent source areas and to calculate anthropogenic heat and storage heat fluxes, aerodynamic and surface resistances, drainage, and storage capacities. The model requires inputs to describe the surface and hourly meteoro- logical data and provides outputs of surface water state and evapotranspiration (Table 1). 3. SUBMoDELS 3.1. Source Area Identification The urban atmosphere interface is extremely complex, defying simple description. For example, in a suburban land 1742 GRIMMOND AND OKE: EVAPOTRANSPIRATION-INTERCEPTION MODEL use zone it is possible at one scale to identify an area which could be considered homogeneous (e.g., a residential area), within which there is a large degree of smaller-scale hetero- geneity (e.g., buildings, vegetation, roads). The collection of surface elements that are sampled by a micrometeorological sensor are termed its "source area." In urban areas the surface is not necessarily homogenous at the scales of the source areas of the turbulent fluxes. The turbulent source area is upwind of the measurement site and in the direction of the prevailing wind. Its upwind, downwind, and lateral boundaries are dependent on the characteristics of the flow and on the boundary layer development in the atmospheric layer between the surface and sensor level. Further, the turbulent fluxes measured at a point have a continuously changing source area due to the variations in the flow and stability of the boundary layer. In order to model evapotranspiration it is necessary to have surface descriptive information for the calculation of QF, lXQs, D, surface resistance, and the overall QE- There- fore it is necessary to have a method to identify the appro- priate source area and to retrieve the surface characteristics for the same area [Schmid et al., 1991]. To identify these areas, a modified version of the Schmid and Oke [1990] source area model was used. The inputs necessary to deter- mine the weighted source areas are the site conditions (roughness and zero-plane displacement) and meteorological scaling parameters (Obukhov stability length, friction veloc- ity, and either the standard deviation of wind direction and mean wind speed or the mixed layer depth). The output of the model consists of nine source area ellipses for each hour (Figure 2). The area between each ellipse has equal weight, that is, the area labeled "a" contributes the same amount to the measured flux as that labeled "b". In the current work the Schmid and Oke [1990] model, which was developed for unstable conditions, was extended using the relations of Grynning et al. [1987] to incorporate stable conditions. The dimensions increase as conditions change from unstable to stable (Figure 2). Combining the source area calculations with a surface description data base and a data accessing system allows a parameter to be assessed in response to different source areas shapes as meteorological conditions and scenarios change. Surface data for 100 x 100 m squares were collected for a 5-km radius circle centered on the tower (i.e., approx- imately 8000 squares). 3.2. Anthropogenic Heat Flux (QF) The energy released due to human activities cannot be measured directly. C.S.B. Grimmond and T.R. Oke (manu- script in preparation, 1991) present the details of the method used to determine Q F on an hourly basis. Q v can be subdivided into three components: QI* = Q•'v + QF• + QI*M W m -2 Qvv, the heat released by vehicles, is a function of the amount of fuel used and its energy content. It can be related to the number of vehicles traveling within the source areas, the distance they travel, the type of fuel, and their fuel efficiency. Q FH, the heat from stationary sources (primarily buildings), is determined from an inventory of the types of consumers and their consumption patterns within the source areas. Q•u, the heat of metabolism, is small and a function N 000 LAT irection = 246 ---44 m u* = 0.17 ms-1 •i•_•" '/.• , N •300 LAT Wind Direction = 310 L = 462 m Fig. 2. Examples of the shape of the source areas calculated using the plume model. The nine weighted bands within the source area are also shown (radius of inner circle, 2 kin; radius of outer circle, 5 kin). of the number of people and animals within the source areas and their metabolic rates. Hence each component requires information in the surface data base which can be accessed to fit the source area for that hour. 3.3. Storage Heat Flux (AQs) In urban areas the subsurface or storage heat flux is the net uptake or release of energy from the urban system. It includes latent and sensible heat changes in the air, build- ings, vegetation, and ground, extending from above roof level to a depth in the ground where the net heat exchange over the period of study is negligible [Oke and Cleugh, 1987]. It is not possible to directly measure the integrated heat storage of the urban system. Here the parameterizati0n scheme of Grimmond et al. [1991] is used. The daytime storage heat flux is given: OQ* AQs=a•Q* + a 2•+ a 3 ot GRIMMOND AND OKE: EVAPOTRANSPIRATION-INTERCEPTION MODEL 1743 where a•, a2, and a3 are empirical coefficients. The second term gives the hysteresis loop departures from the linear relationship. Theempirical coefficients are derived from data for individual surface types and weighted according to their abundance in the source area. The active surface area in the source area involves determination ofthe plan area of green space, the three-dimensional areaof buildings (roofs plus walls), and the plan area of paved surfaces. The surface coefficients are determined each hour for the appropriate source areas using the source area model and the surface data base. At night the storage heat flux is found to be approximately equal to the sum of the net all-wave radiation and anthropogenic heat flux. 3.4. Aerodynamic Resistance (%) For a given wind speed and vapor pressure the rate of removal of water vapor depends on the atmospheric turbu- lence created by the wind blowing over the surface rough- ness elements. The integrated transfer coefficient for water vapor between the evaporating surface and some reference height in the free atmosphere is the aerodynamic resistance (r•) or its reciprocal the aerodynamic onductance (.q,). Under neutral atmospheric conditions, r, is determined using the logarithmic wind profile. In nonneutral conditions, and for water vapor exchange, the following expression is used [Thom, 1972]: % = {In [(z - d)/zo] - ½}{In [(z - d)/zovl - ½7,}/(k2u) where z is the height measurement (m), d is the displacement length (m), z0 is the momentum roughness length (m), k is von K•rmfin's constant (0.41), u is the horizontal wind speed (m s-l); ½ and ½•, are the stability functions for momentum exchange and water vapor, respectively, and z0•,, the water vapor roughness length, is taken to be 0.1z0 [Brutsaert, 1982]. For unstable conditions the stability function for momentum is [van Ulden and Holtslag, 1985]: .6 = 2 In [(1 + X)/2] + In [(1 + X2)/2] - 2 tan -i (X) + rr/2 where X is defined by Dyer [1974] as X = (1 - 16 z/L) 0'25, L is the Monin Obukhov stability length (m), and the stability function for heat and water vapor is [van Ulden and Holtslag, 1985] ½H = 2 In {-(1 + y2)/2] where y is defined by Dyer [1974] as y = (! - 16 z/L) ø'25 For the stable case the functions of Dyer [1974] are used for heat and water vapor: ½•v = -5(z - d)/L and the functions of van Ulden and Holtslag [1985] for momentum: • = -17{! - exp [-0.29(z - d)/L]} L is calculated iteratively with u, using Randerson [ 1984] when sensible h at flux data are available: L = -(u ,3pcpTa)/(k#QH) where # is acceleration due to gravity (9.81 m s-2), c p is specific heat capacity of air at constant pressure (J kg -• K-•), p is density (kg m-3), and T a is dry bulb temperature (øC). L is calculated using van Ulden and Holtslag [1985] when temperature differences are used: œ = where ©,, the temperature scale for turbulent heat transfer, is o, = [(z2 - - 0[(z2 - + - 6.0 is the potential temperature difference [Mcintosh and Thom, 1972] zX© = Ar(lOO/p) g K = R/cp, R is the gas constant of dry air (287.04 J kg -• K-•), and p is pressure (Pa). Here u, is calculated using Van Ulden and Holtslag [1985]: u, = [ku(z3 - d)]/{ln [(z3 - d)/zo] - ½[(z3 - d)/L] + ½(z0/L)} where u(z3) is the wind speed at the height of the anemom- eter (m s -•). When information to calculate the stability correction is not available, the neutral form of r•, is used. This has been found satisfactory over forests [Shuttleworth, 1989]. 3.5. SurJbce Resistance (rs) When there is no rain or dew on vegetation, water loss occurs through the stoma. This process is biologically con- trolled and changes with micrometeorological conditions. The biophysical and biochemical mechanisms of stomatal resistance (rsr) are poorly understood [Shuttleworth, 1989]. The form of stomatal response is qualitatively similar from one species to the next but exhibits marked variations in quantitative terms both between and within species [Jarvis et al., 1976]. At the whole canopy scale, Shuttleworth [1976, 1978, 1979] shows that the bulk stomatal resistance is closely related to the surface resistance (rs), as defined by Monteith [!965]. In this approach the complex array of r st for the individual leaves is replaced by an equivalent system con- taining a single hypothetical leaf with a single rs. Thom [1972] showed that if the source of evaporation is a forest canopy, and the aerodynamic resistance correctly describes the transfer of water, then rs is equivalent to the bulk rsr. In the more general case where some of the water vapor comes from other layers (e.g., litter or understory), changes in rs also reflect changes in their rates of evaporation [Stewart, 1988]. In urban areas the diversity of vegetation makes the determination of a bulk stomatal resistance virtually impos- sible because of sampling problems and the added compli- cation posed by the nonvegetated areas. Here rs is modeled as a single integrated resistance for the whole system. At present there is insufficient knowledge to interpret the effects of environmental variables on rst using a mechanistic model [Stewart, 1988], but a physically based model pro- posed by Jarvis [!976] can be used to predict r st from empirical relationships between rsr and controlling environ- 1744 GRIMMOND AND OKE' EVAPOTRANSPIRATION-INTERCEPTION M DEL mental variables. This form of model is promising for use in urban areas. It has the advantage of being dynamically responsive to both meteorological and seasonal conditions. The general form of whole-canopy models for forests [e.g., Lindroth, 1985; Gash et al., 1989; Sellers et al., 1986; Shuttleworth, 1988, !989; Stewart, 1988], expressed as a conductance, (qs) is gs = G•L#(variables) where G1 is the maximum value of the surface conductance (m s-l), L is the leaf area index (LAI), and #(variables) are the functions of the environmental variables which have values between zero and unity. This general form was used to develop a model for the Sunset study site. It assumes that surface conductance is a function of the nonsynergistic interaction between stomatal aperture, LAI, and environmental variables. It operates on an hourly time scale. The environmental variables identified by Jarvis [1976] were not available and were replaced by closely related ones: quantum flux density was replaced with net radiation (Q*), leaf-air specific humidity was replaced with specific humidity deficit ($q) [c.f. Stewart, 1988], leaf temperature was replaced by air temperature (T), and leaf water status was replaced by soil moisture deficit ($0), so that gs = Gjg(Q*)g(Sq)g(T)g(80)g(L) The dependence on Q* is expressed using the same mathe- matical form of equation used for solar radiation [e.g., Shuttleworth, 1989]' g(Q*) = [Qm/(G2 + Q*)]/[Qm/(Qm + G2) ] where Qm is the annual maximum hourly net radiation (W m-2), and G2 is a fitted parameter (W m-2). The dependence on $q is expressed as [Gash et al., 1989; Stewart, 1988] g($q) = 1 - G3t5 q O< $q< G 4 g( •q) = 1 - G3G 4 •q -> G4 The dependence on T is given [Stewart, 1988]' g(r) = (r- r)(rn- T)rc/[(r- rD(rn- r)rq where Tc = (T•q- Gs)/(G 5 - Tœ) and TL and Tn are minimum and maximum temperature limits (øC). The dependence on $0 is described by [Stewart, •988] g(80) = 1 - exp {G6[t50 - (S1/G 6 + S2)]} where S1 and S2 are fitted parameters which relate to the maximum •0. The dependence on LAI follows Dolman et aI. [ 1988] but allows for the greater diversity of vegetation found in a suburban area compared to a temperate forest: g(L) = [(LA u)/Lm) + At]/(A u + At) where A u is the area unirfigated (m2), Lm is the maximum leaf area index in A u, and At is the area irrigated (m2). The values of the parameters (G1-G6) for the study area were optimized using nonlinear least squares regression [Denis et al. , 198!; Gay, 1983; Moore, 1986] between"rnea. sured" values of r s and the selected environmental v ri. ables. "Measured" surface resistance (rs) is determined when the surface is dry by rearranging the Penman-Monteith equation and using measured Qlc and/3 [Monteith, 1965] rs = - 1 r a + yQœ 3.6. Drainage (D) Determination of the drainage rate for urban areas is a problem because it does not correlate directly with urban runoff. It has to be envisaged as that loss from the "storage layer" which is no longer available for evaporation, for example, the water may enter a pipe system or infiltrate into the soil. Given the range of vegetation and other surface types in urban areas, the measurement of the drainage rate involves a vast sampling problem. Even in forest areas there is a significant error involved in the measurement of this term [Shuttleworth, !989]. The urban drainage function developed here describes the rate at which water drains from a single-layer moisture store (Figure 1) and is set proportional to the current water status of the store. The original Rutter et al. [ 1971, 1975] formula. tion had the inherent difficulty of predicting a finite drainage when the canopy is dry [Calder, 1977]. It has been rewritten [e.g., Halldin et al., 1979] to rectify this so that it reads D = D0 exp [(bC) - 1] (5) where D 0 is the drainage rate when C = S, b is an empirical coefficient, and S is the amount of water the canopy retains after the rainfall and throughfall cease (the storage capacity). This equation is used in the model for the proportion of the source area which is vegetated but not irrigated. For the impervious area, values from the urban runoff literature are used. They usually relate to the runoff from a single surface type rather than a land use area. Falk and Niemczynowicz [1978] propose the following equations for the relationship between drainage and storage: D = Do(C t_ 1 - S) b (6) D = Do(C t_ i) b (7) Equation (6) allows no drainage when C < S. This assump- tion has also been used in forested areas, for example by Gash and Morton [1978] and Mulder [1985]. Pratt et al. [ 1984] suggest that the rate of drainage from flat roofs can be assumed to be the same as paved surfaces of 1% slope. Equation (7) is used for the paved, built surfaces and for irrigated grass. The values of the coefficients for (5)-(7)are taken from Calder and Wright [1985] and Falk and Niemc. zynowicz [1978] and are listed in Table 2. 3.7. Storage Capacity (S) The values of the surface storage capacity, like the drain- age functions, are based on values taken from the literature (Table 2). The capacity for a surface type is assumed to be constant through time for all types except deciduous vege- tation, where allowance is made for seasonal change. The values used are as follows: paved surfaces, 0.48 mm [Falk and Niemczynowicz, 1978]; roofs, 0.25 mm [Davies and GRIMMOND AND OKE' EVAPOTRANSPIRATION-INTERCEPTION MODEL 1745 TABLE 2. Values Assigned to Parameters for the "Base" Run of the model in Vancouver, 1987 Surface Characteristics Unirrigated Irrigated Pavement Building Coniferous Deciduous Grass Grass S, mm 0.48 0.25 1.2 0.3 1.3 1.3 D function 7 7 5 5 5 7 Coefficient Do 10 10 0.013 0.013 0.013 10 b 3 3 1.71 1.71 1.71 3 C(t = O) 0.0 0.0 0.0 0.0 0.0 0.0 , Parameter Value Site Characteristics d, m 3.50 Z0, m 0.52 Deciduous Transitions T• YD 87/65 Te YD 87/115 Sas, mm 0.80 Surface Conductance Lm 3.10 Tt•, øC 40.00 Tœ, øC 0.00 S 1, mm 0.45 S2, mm 15.00 Q max, W m -2 725.00.00 Gi, mm s -l 53.95 G2, W m -2 633.95 G3, kg g-• 0.0821 G 4, g kg -1 8.91 Gs, øC 18.88 G6, mm -• 0.0107 rs max, s m-l 9999.00 Hollis, 1982]; coniferous trees, 1.2 mm [Shuttleworth, 1989]; deciduous trees when leafless, 0.3 mm; deciduous trees when in leaf, 0.8 mm [Shuttleworth, 1989]; and grass, 1.3 mm [Zinke, 1967]. 4. FIELD MEASUREMENTS In order to test the evapotranspiration model and its submodels, data were collected in the Sunset suburb of Vancouver, B.C. (49 ø 15'N, 123 ø 18'W) (Figure 3).' The city is situated atthe west end of the Lower Fraser Valley and has low to moderate relief. The Sunset area has predomi- nantly one- and two-story single family dwellings and slopes slightly southwest toward the Fraser River. The objective was to determine spatially averaged fluxes from an area of suburban land use at the local scale [Oke, 1984], rather than from individual surfaces such as lawn or pavement. Therefore measurements were conducted from a tower well above the height of the roughness elements to ensure that observations were conducted in the surface layer. On the basis of the results of Roth et al. [1989] all sensors were mounted at heights >20 m above zero-plane displacement, and the averaging time was 60 min. The environmental variables measured were net radiation 1tux density (with a net pyrradiometer), s nsible h at flux density (with asonic anemometer-thermometer), wind speed and irection, wet- and dry-bulb air temperature (s e below), relative humidity, surface wetness ( ee below), precipita- tion, and soil moisture deficit (see below). Full details of the measurements are available from Grimmond [1988]. All times referred to are local apparent time (LAT) (plus 0800 UT) and year day (YD). The reversing temperature difference measurement sys- tem used to determine/3, was similar to that described by McCaughey et at. [1987]. The wet- and dry-bulb thermocou- ples (10 junction copper/constantan thermopiles) were mounted within shields on two carts which move up and down on a trackway. The carts, separated by 7.1 m, are reversed twice each hour. A ten-min period allows for travel and equilibration at the new level. The remaining 20 min at each level allowed for two 10-min averages. The surface wetness sensors followed the design of Weiss and Lukens [1981]. The sensors consist of a frame (0.1 x 0.1 m) with 30 AWG wire threaded across the center and muslin material weaved through the wire. The wire has 2 V ac passed through it, and the output voltage varies depending on the wetness of the material (a surrogate for the surface). The output does not give a numerical value of degree of wetness but an indication of whether the surface is wet or dry. Two such sensors were located on unirrigated grass and a third in a coniferous hedge, starting on YD 87/38. Visual observations of the drying of different types were also made throughout the measurement period. Gravimetric soil moisture samples were taken approxi- mately three times a week from two sites near the tower. 1746 GRIMMOND AND OKE: EVAPOTRANSPIRATION-INTERCEPTION M DEL parks. recreational, agricultural and undeveloped are.s industrial. commercial and institutional areas English B•y o 1 2 Kilometres Fig. 3. The location of the Sunset study site in the Vancouver area and the surrounding land uses. From mid-May (YD 87/133), samples were gathered at an additional three sites, both irrigated and nonirrigated, around the tower area. The average dry weight moisture content for the profile from the surface to a depth of 200 mm was determined and used with dry bulk density measure- ments to determine soil moisture deficit. Water use was monitored from YD 87/139 by the Engi- neering Department, City of Vancouver, in the nearby residential area of Oakridge (Figure 3). The water use for the 21-ha catchment was recorded at 15-min intervals with a resolution of 2.832 m 3. The external water use (i.e., mainly for irrigation) was separated from the total use on the basis of experience in the same catchment [Grimmond and Oke, 1986]. It was assumed that the same depth of water was applied to the irrigated area in the Sunset suburb as was applied in the Oakridge catchment. From visual surveys at both sites it was determined that different proportions of the residential vegetation were irrigated (Oakridge, 95%; Sun- set, 60%). The external water use data set was extended back to YD 87/121, selected as the day of the onset of irrigation (based on observations and Grimmond [1983]), by using a multiple regression equation developed from the hourly 1987 data set incorporating time of day, hours since rain, available energy, temperature, and vapor pressure deficit. 5. DEVELOPMENT OF THE SURFACE CONDUCTANCE SUBMODEL Figure 4 shows the ensemble mean and standard deviation of g s at the Sunset suburban site for 543 hours. The diurnal trend is similar to that reported for forests [e.g., Shuttle- worth, 1989], but the values are approximately 75% smaller and have approximately double the hourly standard evia- tion [e.g., Stewart and de Bruin, 1984; Shuttleworth, 1988: Stewart, 1988]. The smaller size of #s in urban areas i undoubtedly due to the presence of nonvegetated surfaces. The surface resistance submodel was developed and tested using data for 543 daylight hours when the surface was dry. To be dry, a surface must meet the following criteria [Shuttleworth, 1988]: rain is not recorded in the current or previous 3 hours; 3 hours with positive Q* have occurred since rainfall; the hour preceding recorded rain was regarded as not dry; in the morning, 2 hours have elapsed since Q* first became positive; the wetness ensors record the surface as dry; and there are no visible puddles or other traces of rain. To develop and test the model, the Sunset data set was split into two independent sets by arbitrarily putting alter. nate days into separate s ts: one of 300 hours (AD1) and the other of 243 hours (AD2). This approach was adopted because no independent data set exists for a suburban rea GRIMMOND AND OKE' EVAPOTRANSPIRATION-INTERCEPTION MODEL 1747 16 14 12 ! !•!•!ii:•.. i!!i:::: i::: ?•?i :!i!i•!!•!i!i !i•i5 ! i !ii i i!! !:•i::: :5:5::!!: : •:i::!: i :!::i !:::i : : i i :' !!:!i• !!!:•:i!i:i. : ::::::::::::::::::::::::::::: 0 6 8 lb 12 1• '1•5 18 Time of Day Fig. 4. Diurnal trend of ensemble mean and standard deviation of measured surface conductance for Sunset suburban site, Vancou- ver, British Columbia, for dry, daylight hours. The shaded area shows the range of values from eight forests [from Shuttleworth, 1989, Figure 10]. g • 8 which covers a sufficiently wide range of conditions. The environmental variables controlling conductance were determined as follows. Q* and T were measured directly, $q was calculated from measured T and pressure p. •0 was determined for unirrigated areas using profiles of gravimetric soil moisture and the dry bulk density of the soil. Irrigated areas were found to have no $0 at any time. The surface description data base and the source area model enabled the areas with coniferous, deciduous, and irrigated and unirrigated grass to be identified for each hour. Thus an areally weighted 80 was determined for each hour. It was assumed that all of the irrigated area is grass and that mowing keeps the LAI of the grass constant at 1.6 [Ripley and Redmann, 1975]. The LAI for the unirrigated area was determined as an areally weighted average of the coniferous, deciduous, and grass values for each hour. The LAI of coniferous and deciduous vegetation were set to lie in the ranges of 4-6 and 0-4, respectively [Kramer and Kozlowski, 1979]. A full circle panorama of photographs was taken from the tower weekly to identify the stages of the growing season. They were used to model LAI for each day for the two vegetation types using [Street and Opik, 1984] L+I log = k[(d- do) - (dl a-(L+l) - ao)] where a is the maximum LAI minus the winter minimum, k is a fitted parameter (0.05), d is the day of interest, d o is the starting day of the curve from winter minimum of LAI (YD 87/45), and d l is the day of inflection of the S curve (YD 87/95). It was possible to assign a daily value of LAI for each vegetation type and, depending on the source area being sampled, to determine the weighted average of LAI for each hour. The proportion of the source area irrigated was deter- mined from visual surveys and water use data in conjunction with the data base and the source area model. The other model parameters were assigned as follows: maximum areally weighted Lm is 3.1, based on the calcu- lated LAI; Qm was set to 725 W m -2 using measured Q* during June; T• was set to 40øC and TL to 0øC following Stewart [1988] and Gash et al. [1989]; and Sl and S2 were TABLE 3. Fitted Parameters and Test Statistics for the Surface Resistance Model G1, G 2 , G 3 , G4 Gs, G6, Model n mm s -I W m -2 kg g-I g kg 'I øC mrrr -t (a) Fitted Model Parameters AD1 300 56.35 694.0 0.0857 9.44 20.07 0.0105 AD2 243 53.95 634.0 0.0821 8.91 18.88 0.0107 All 543 55.20 669.6 0.0840 9.02 19.56 0.0106 c, RMSE, RMSEs 1, RMSE,,•, Model Test r 2 m0 m mms -• mm s -• mm s mm s- d n AD I AD 1 0.61 0.90 AD 2 AD 2 0.71 0.93 AII AII 0.65 0.91 AD 1 AD 2 0.71 0.89 AD 2 AD 1 0.61 0.92 May-June AD 1 AD 2 0.72 0.89 AD 2 AD 1 0.65 0.92 (b) Model Developtnent Statistics 0.62 2.37 2.35 1.45 1.85 0.87 0.73 1.47 1.75 0.88 1.52 0.91 0.66 1.99 2.11 1.21 1.73 0.89 (c) Model Test Statistics 0.74 1.20 1.79 0.91 1.54 0.91 0.61 2.60 2.36 1.47 1.85 0.87 0.74 1.14 1.75 0.89 1.50 0.91 184 0.65 2.49 2.26 1.38 1.79 0.88 202 Statistics are n, number of hours; m0, slope forced through t e origin; m and c, slope and intercept for ordinary linear regression; RMSE, root-mean-square error;RMSEs systematic error; RMSE,,, unsystematic error; d, index of agreement [Willmott, 1984]. 1748 GRIMMOND AND OKE: EVAPOTRANSPIRATION-[NTERCEPTION MODEL 14 •12 '610 , 8 • 6 g4 •o 2 Fig. 5. Comparison between hourly measured and modeled sur- face conductance using data set AD2 for modeled AD1. TABLE 4. Statistics forthe Comparison Between Measured an Modeled (Using "Base" Values) Q•z for Vancouver, British Columbia Statistic Hourly Daily Mean, W m -2 Measured Modeled Standard deviation W m Measured Modeled m 2 I' d N&S* RMSE, W m -2 RMSE s, W m -2 RMSEt,, W m -2 -2 2944 125 40.10 40.49 40.38 41.87 63.03 21.80 58.20 23.87 0.83 0.93 0.81 0.71 0.95 0.91 0.81 0.64 27.65 13.00 !0.72 2.17 25.49 12.82 *Nash and Sutcliff? [1970] goodness of fit. determined by trial and error [Doltnan et al., 1988] to be 0.45 and 15 mm, respectively. Table 3 lists the fitted parameters and the test statistics for the separate and combined data sets. The values of the model parameters G3 to G 6 are similar to those reported by Stewart [1988] and Gash et al. [1989] for temperate forests (Thetford, U.K., and Les Landes, France, respectively). G i and G2 are not the same as those of these forested sites. This is expected for G1 because the urban canopy vegetation is more heterogeneous and sparse. Second, G2 is associated with net all-wave radiation in this study rather than net solar radiation, as in the other studies. This is the best form of the model after several attempts to incorporate urban features. The statistics pertaining to model development and per- formance are presented in Tables 3b and 3c, respectively, and in Figure 5. The model was first developed using data set AD2 and tested against AD1, then the process was reversed. The r 2 of 0.71 for model AD2 is similar to that reported by Stewart [1988]. AD1 always gave a lower degree of expla- nation. This is attributed to a greater number of outliers associated with cloudy and pa_rtly cloudy hours, suggesting that the net radiation function is inadequate for such condi- tions. All three models show larger unsystematic than sys- 8 3 2 1 . measured ø modelled 1•2 1'4 1• 1'8 ..... Time of Day Fig. 6. Ensemble mean of measured and modeled surface conduc- tance using data set AD2 for modeled AD1. (a) 300. '• 200- o ! • 100. CY 0 -lOO -lOO 1!}o 200 300 (b) QE (W m-2)- measured ,, x x xx.,•:•x x X • x • , 150- ' 50- CY 0- -50 -50 0 50 150 QE (W m-2) - measured Fig. 7. Comparison of measured and modeled Q•r for the "base" model run for Vancouver, British Columbia: (a) hourly (b) daily. GRIMMOND AND OKE' EVAPOTRANSPIRATION-INTERCEPTION MODEL 1749 tematic errors, which indicates there is only a small bias in the model estimates. The slope of the best fit line indicates an overprediction for low -qs (high r s) and an underpredic- tion for high #s (low rs). The slopes forced through the origin (Table 3c) are similar to those obtained by Stewart [1988]. Figure 6 shows the diurnal trend of the ensemble mean measured and modeled # s- The modeled and measured data exhibited the same trend using either model. In both cases the model underpredicts in the morning. The AD1 model was selected to calculate r s in the evapotranspiration-interception m del because it has the lowest root-mean-square error (RMSE), and the slope of the line (m) is the closest to unity. 6. PERFORMANCE OF THE URBAN EVAPOTRANSPIRATION MODEL 6.1. Methods Used to Assess the Performance of the Model The model was initialized, with the surface dry, on Feb- mary 3 (YD 87/33) at 0000 LAT and terminated on June 28 (YD 87/179) at 0800 LAT, when field measurements ceased. During this 35'13-hour period, QE was modeled continu- ously, and there are 2944 hours of measured data available for comparison. The performance of the model is evaluated at a number of temporal scales: (1) hourly, each hour with both measured and modeled flux data, (2) daily, fluxes for all hours, when both measured and modeled ata are available on a day, are averaged and compared, because of missing data a "day" does not necessarily consist of 24 hours, (3) cumulative hourly, hours when both measured and modeled data are accumulated an  compared (because of missing data the values are less than the "true" total E for the time period (YD 87/33-179)), and (4) ensemble hourly, for each month and for the whole time period (YD 87/133-179) the values for a particular hour are averaged. The modeled surface water state was also compared with that from two wetness sensors. 6.2. "Base" Run of the Model For the "base" run of the model the values assigned to the various parameters were chosen tobe the most physically realistic (Table 3). They are not a best fit to the model. The statistics of the model's performance as compared to the measured Q E are presented on an hourly and daily basis in Table 4. The results for the model are better on an hourly than a daily basis primarily because of the large number of hourly data points (Figures 7a and 7b). In both figures the data re distributed around the one to one line, and the hourly RMSE of 27.7 W m -2 is of the same order as the error of the measurements. The cumulative E for the measure- ment period (Figure 8) shows that he model underpredicts up to about YD 87/65, overpredicts between YD 87/140-160, and matches for the remainder of the period. Figure 9 shows the time series of measured and modeled Qœ for two periods. In the early part of the year the model underpredicts during he middle ofthe day (Figure 9a). From about YD 87/75 the model mimics the measured data ex- tremely well under all conditions except during dewfall (Figure 9b). Figure 10provides a means tocompare the modeled surface state with that recorded by the wetness ensors. In 150 --- measured / --- base // _ // [..13 100- / .• - /? 5O 43 107 171 Time (D) Fig. 8. Cumulative m asured and modeled ("base" run) evapora- tion for Vancouver, British Columbia. each box of five traces the upper two plots show the water arriving at the surface; at the top is the rainfall received by all surfaces, and the second is the external water use (wu) which is only received by the irrigated grass portion of the (a) .... measured 20• •- ----- modelled 150' oo. I 't! . ,,. -so 4'8 ' so s2' s4 Day (b) ................... - , .... , , . ,.,, ,• ........... . . ., 200' .•, 150- 100' • 50- 0' -50 140 142 144 146 148 Day Fig. 9. Time series of hourly measured and modeled Q E ("base") for Vancouver, British Columbia. Note that the measured data are not continuous (a) YD 87/47-57 (b) YD 87/140-150. !750 GRIMMOND AND OKE' EVAPOTRANSPIRATION-INTERCEPTION MODEL (a) 3 P (mm) 0 WU 1 3 (2 (ram) whole 0 3 C (mm) 0 wet dry grass dry (b) 56 60 64 68 72 Day P (mm) WU wet dry grass dry 3 C (ram) whole li8 ' 1i2 126 130 134 Day Fig. 10. Time series of hourly precipitation, water use, surface water use, surface wetness, and modeled surface water state for Vancouver, British Columbia (see text for details) (a) YD 87/54-74 (b) YD 87/116-136. surface. External water use did not begin until YD 87/121. The third plot shows the status of the wetness sensors (i.e., whether wet or dry). Note that the sensor located on the grass (dashed curve) has a smaller range. The fourth trace is the modeled surface water status of the whole surface, and the fifth or lowest plot shows the modeled wetness status of the individual surface types. The six individual surfaces are not separately labeled, but the vegetated surfaces are those which hold more water on the surface, and the lesser values are the built and paved surfaces. Once irrigation begins (Figure 10b), the state of the irrigated grass surface (the dashed curve in the lowest plot) becomes distinguishable from the other surfaces. In general, the model does an extremely good job of indicating when the surface is dry. Visual observations indicate that the paved and built sur- faces dry more quickly than the vegetated surfaces, and as Figure !0 shows, the modeled paved and built surfaces mimic this pattern. This suggests that he drainage functions in the evaporation model simulate the surface water state appropriately. The only feature not well represented is dewfall which was visible on the grass ometimes but not simulated. 6.3. Sensitivity Analysis of the Model In order to assess the robustness of the model and its sensitivity to various input parameters and equations it performance was investigated when inputs were changed over a range. The effects of such changes are reported in comparison to both the measured results and the base run. Unless otherwise specified only one variable was changed from the base run (Table 3). Cumulative plots show the influence of the changes on the modeled evaporation (E) and allow seasonal effects to be gauged. 6.3.1. Time step. When this type of model is used in forests, time steps of 120 [Shuttlea,orth, 1988], 300 [Rutteret al., 1971], and 720s [Sellers and Lockwood, 1981] have been used. These rather short time steps allow tbr changes in the drainage rate, surface state, and evapotranspiration rate within an hour. The "base" time step used in this study is 300s, that is, 12 steps per hour. Decreasing the computational interval to 120 s does not improve the performance of the model and extending it to 720 s shows slight improvement. Changing to 3600 s causes a harther slight improvement in the statistics, but the surface dries too quickly so the surface water state shows a series of spikes rather than the more realistic gradual drying curves seen in Figure !0. We conclude that the model can be run at the longer 720-s time step and maintain a realistic surt•tce state, but extending to one hour results in less realistic predictions of surface water state. 6.3.2. Aerodynamic' resistance %. The base version0f the model uses the stability corrected equation for r a. Runs made using the neutral and the Thom and Oliver [1977] equations for r a indicate that the model is not very sensitive to which equation is used, and Qœ values are similar. Therefore it is not essential to have the input necessary to calculate the stability functions for the determination of r•. In many forest evaporation-interception studies the neutral method has been used alone and found to be satisfactory [Shuttleworth, 1988]. 6.3.3 Surface resistance rs. Two types of changes af- fecting rs were run: first, alterations in r s maximum and second, changes to the values of the parameters within the r s submodel (section 3.5). The model was run with five different methods to deal with large rs (or g s approaching zero) and negative Q*. In the base run, when g s approaches zero or becomes negative, rs is set to 9999 s m -•. In the other four cases, rs was et to 9999, 2500, 1250, or 750 s m -• when Q* was negative. The perforimance statistics how little difference between the base run and those with rsmax = 9999 s m-•. As the value of r smax is decreased, the performance d teriorates slightly. There is, a systematic influence on the cumulative flux, primarily due to the fact hat as r s max is reduced, nighttime values of Q• increase, and therefore accumulated Q•be- comes larger. However, the daytime performance is not improved. GRIMMOND AND OKE: EVAPOTRANSPIRATION-INTERCEPTION MODEL 1751 TABLE 5. Statistics of Model Performance for QE When the Storage Capacity (S) is Changed for Irrigated and Unirrigated Grass S Capacity (Irrigated, Unirrigated), for Grass, mm 0.5, 0.5 1.3,* 1.3 6.4, 6.4 7.6, 7.6 7.6, 10.2 Hourly Statistics Mean, W m-2 43.08 40.38 38.03 37.91 37.87 sd, W m -2 59.91 58.20 56.72 56.65 56.65 m 0.84 0.83 0.82 0.82 0.82 r 2 0.78 0.81 0.83 0.83 0.83 RMSE, W m -2 30.07 27.65 26.42 26.39 26.40 RMSEs, W m -2 10.58 10.72 11.69 11.75 11.76 RMSEu, W m -2 28.15 25.49 23.70 23.63 23.63 d 0.94 0.95 0.95 0.95 0.95 N&S? 0.77 0.81 0.82 0.82 0.82 Daily Statistics Mean, W m -2 44.66 41.87 39.44 39.31 39.25 sd, W m -2 26.09 23.87 22.45 22.39 22.42 m 0.98 0.93 0.98 0.89 0.89 r 2 0.67 0.71 0.75 0.75 0.75 RMSE, W m -2 15.60 13.00 11.58 11.53 11.52 RMSE s, W m -2 4.20 2.17 2.61 2.69 2.68 RMSEu, W m -2 15.02 12.82 11.28 11.21 11.21 d 0.89 0.91 0.93 0.93 0.93 N&S 0.49 0.64 0.72 0.072 0.72 References for S capacity values are as follows: 1.3, Zinke [1967]; 6.4, Brater [1968], Mechler and Riecken [1977], and Nouh [1986]; 7.6, Wright-McLaughlin Engineering [1969]; 10.2, Wright-McLaughlin Engineering [1969]. *Base value. ?Nash and Sutch.'ffe [19701. The r s parameters (G•-G6) were changed between the three sets of model parameters presented in Table 3 (AD2 is the base set). The base results are generally slightly better. Using AD1 and All parameters, the calculated daytime fluxes are slightly less than for the base case. The values obtained using AD1 are smaller than those using the All parameters. The AD1 and All versions of the runs remain closer to parallel with the measured cumulative E during the later months. The model responds toreasonably small changes in rs (as will be produced bychanging rsmodel parameters), but it is not so sensitive as to make the use of any of the models in section 5 unwise. 6.3.4 Storage capacity (S). Three types of changes were made to the model relating to the storage capacity: first, changing S for an individual surface while keeping the rest at the base value, second, changing all surfaces values, and third, changing S in conjunction with the drainage functions. The results of the latter are discussed in section 6.3.5. In this first est he value of S for each surface type was changed over a range of values determined from previous studies. It should be appreciated that the response of the model a so depends on the relative abundance of the surface type within the study area. Paved and building: When S was changed from 0.30 to 2.50 mm for paved and from 0.25 to 2.5 mm for buildings, there were only slight decreases in mean modeled values of Q• and its RMSE, RMSEs, and RMSEv. This is reassuring because it is difficult tomeasure and/or assign absolute values to this parameter [see Falk and Niemczynowicz, 1978]. Trees: Changing S for coniferous vegetation from 0.75 to 2.5 mm, the size of the winter and summer values of S for deciduous vegetation, and the timing of the transition be- tween winter and summer are virtually undetectable. How- ever, it should be noted that deciduous and coniferous vegetation occupy only a small area of the study site. Grass: S was varied over the large range from 0.5 to 10.2 mm that is found in the literature (Table 5). The irrigated and unirrigated grass are assigned the same value except for one case where the unirrigated is given 10.2 mm and the irrigated 7.6 mm. The model shows a distinct response to the size of S for grass because the changes are much larger than other surface types and because grass is the most extensive surface type in the study area. It is apparent from Figure 11 that the larger the value of S the smaller the modeled E, down to a limit of 6.4 mm above which S has no further influence. Once external water use begins (YD 87/121) on a portion of the grass, there is a distinct increase in the modeled values of E (i.e.,a steeper slope in Figure 11). The statistics show slightly contradictory effects of changing S at the hourly and daily scales. With increasing S the hourly statistics how a decreased slope m, and increased r 2, a decreased RMSE and RMSEt;, but an increased RMSEs. On a daily basis, however, increasing S is associated with a decrease in RMSE and in both the RMSEs and RMSEt;. The hourly ensemble data show that the largest difference from the base results occurs at 1800 LAT. In the winter the difference at this time is -<7 W m -2 for S = 0.5 mm and <20 1752 GRIMMOND AND OKE: EVAPOTRANSPIRATION-INTERCEPTION MODEL 150 -• • .• --- measured ......... 0.5 • 1.2t 6.4 ------ 7.6 ............. 7.6/10. 43 1{}7 Time (D) Fig. 11. Influence of the value assigned to the surface storage capacity of grass on cumulated E. W m -2 in June. With large values of S the modeled Q• is reduced by 2-3 W m -2 during the winter, and in June the maximum difference is approximately 6 W m -2' In the -2 morning and 15 W m in the late afternoon. it appears from this sensitivity analysis that the size of S for grass which influences the modeled Q v: is around the value of 1.3 rnm suggested by Zinke [1967]. In four runs all the values of S for the different surfaces were changed together (see Table 6 for details): first, all values of S were decreased to small values, second, all were assigned large values, third, all were set to 0.5 mm, and fourth, all were set to 1.0 mm. All of the runs underpredict Qœ until the onset of irrigation, as does the base run. After that the slope of the cumulative curve steepens considerably for all runs except for the case where all S values were large. After about YD 87/145 the slopes for all curves, except the base run, are similar to the measured values, suggesting they are performing well on a hour-to-hour basis. The comparison statistics suggest that when S is set to larger values the model performs better (Table 6). 6.3.5 Drainage fi•nctions (D). The drainage functions were varied over a range of coefficients and equation types (see Table 7 for details). Initially, the storage capacities for all surfaces were set to one value then changed to 0.5, 1.0, and 2.0 min. The largest cumulative values are obtained when S = 0.5 mm. With all of the models the beginning of external water use is associated with a steep rise in the slope of the cumulative curves. When S = 0.5 mm, all of the cumulative curves cross the measured E curve. When S = 1.0 mm and (7) is used, the cumulative E is less than the TABLE 6. Statistics of Model Performance for Qt' When the Storage Capacity (S) is Changed for All Surface Types S Capacity For All Surfaces, mm Small Base Large All All S pavement 0.30 0.48 2.50 0.5 1.0 S building 0.25 0.25 2.50 0.5 1.0 S coniferous 0.75 1.20 2.50 0.5 1.0 S deciduous* 0.20 0.30 0.40 0.495 0.995 S grass 0.50 1.30 7.60 0.5 1.0 Hourly Statistics Mean, W m -2 43.26 40.38 36.85 47.39 45.72 sd, W m -2 59.92 58.20 56.92 63.03 62.33 m 0.84 0.83 0.82 0.87 0.87 r 2 0.78 0.81 0.83 0.75 0.77 RMSE, W m -2 30.16 27.65 26.48 33.40 31.71 RMSEs, W m -2 10.66 10.72 !1.75 11.15 10.13 RMSE,,, W m -2 28.21 25.49 23.73 3!.49 30.05 d 0.94 0.95 0.95 0.93 0.93 N&S? 0.77 0.81 0.82 0.72 0.75 Daily Statistics Mean, W m-2 44.87 41.87 38.09 49.20 47.39 sd, W m-2 25.99 23.87 23.21 29.27 28.77 m 0.97 0.93 0.92 1.04 1.05 r 2 0.66 0.71 0.76 0.60 0.63 RMSE, W m -2 15.66 13.00 11.76 20.46 18.90 RMSEs, W m -2 4.42 2.17 2.87 8.76 6.97 RMSEu, W m -2 15.03 12.82 11.40 18.49 17.57 d 0.88 0.91 0.93 0.83 0.85 N&S 0.48 0.64 0.71 0.12 0.25 Deciduous values for summer, using the same headings as above, are 0.60, 0.80, 0.90, 0.5, and 1.0, respectively. *Deciduous value for winter ?Nash and Sutcliffe [1970]. GRIMMOND AND OKEi EVAPOTRANSPIRATION-INTERCEPTION MODEL 1753 TABLE 7. Statistics of Model Performance for QE When the Drainage Functions are Changed ands =0.5ram Drainage Functions Rutter et al. v* vl' (7), (5)õ (5)11 (7)$ (7)• [1971]ô Coefficient Do v v 10 0.013 0.018 32 12 0.0014 Coefficient b v v 3 1.71 1.76 1.5 1.5 5.25 Hourly Statistics Mean, W m -2 40.38 47.39 41.88 49.92 48.82 46.77 46.90 45.11 sd, W m-2 58.20 63.03 59.77 65.65 64.09 61.26 6!.32 60.67 m 0.83 0.87 0.83 0.88 0.87 0.84 0.84 0.84 r 2 0.81 0.75 0.77 0.72 0.73 0.75 0.75 0.76 RMSE, W m -2 27.65 33.40 30.26 36.78 35.32 32.83 32.88 32.13 RMSEs, W m -2 10.72 11.15 10.58 12.23 11.99 11.98 12.01 11.47 RMSEu, W m -2 25.49 31.49 28.35 34.69 33.22 30.57 30.61 30.02 d 0.95 0.93 0.94 0.91 0.92 0.93 0.93 0.93 N&S 0.81 0.72 0.77 0.66 0.69 0.73 0.73 0.74 Daily Statistics Mean, W m -2 41.87 49.20 43.41 51.94 50.71 48.31 48.44 46.75 sd, W m -2 23.87 29.27 26.77 28.96 28.49 28.56 28.54 27.95 m 0.93 1.04 1.01 0.97 0.99 1.03 1.04 1.02 r 2 0.71 0.60 0.68 0.54 0.57 0.62 0.62 0.64 RMSE, W m -2 13.00 20.46 15.37 22.72 21.35 19.24 19.30 17.99 RMSEs, W m -2 2.17 8.76 2.93 11.46 10.23 7.85 7.98 6.28 RMSE•, W m -2 12.82 18.49 15.09 19.61 18.74 !7.56 17.57 16.86 d 0.91 0.83 0.89 0.79 0.81 0.84 0.84 0.86 N&S 0.64 0.12 0.50 -0.09 0.04 0.22 0.22 0.32 Here, v stands for variable. *Base equations with base S capacities. ?Base equations but S = 0.5 mm for all surfaces. $Coefficients [Falk and Niemczynowicz, 1978]. õCoefficients [Calder and Wright, 1985]. [ICoefficients [Calder, !977]. ôCoefficients [Shuttleworth, 1988]. measured E, but by the end of the study period they converge, indicating overprediction in the latter period. The same quation underpredicts the cumulative E when S = 2.0 mm but maintains a more similar slope to that of the measured E. This indicates that care must be taken in assigning drain- age functions and coefficients. It appears that those selected in this study for the Sunset site simulate the drainage and surface water state well. 7. DISCUSSION AND CONCLUSIONS The objective of this study was not to model the evapo- transpiration f r a particular time and place but to develop a method to continuously model E over a range of meteoro- logical conditions in urbanized terrain. There has been no attempt to calibrate the model to the measured data. There- fore the good performance of the model, using physically realistic input values obtained from the literature, allows us to recommend this type of approach to modeling urban evapotranspiration a d surface water state in other cities. Although t e model does well on average it generally underpredicts in the early part of the year and performs better in the spring/early summer inVancouver. The mod- eled QE is always greater than zero because dewfall is not predicted even though it was observed on grass on some nights. This failure is one of the reasons why the model underpredicts in themorning and should be corrected. The measured Q E suggests that dewfall can be an important part of the urban energy and water balances in these months and should not be ignored [Grimmond and Oke, 1986]. If the model is to be applied to other cities, the changes required primarily relate to the assignment of parameter values. From the sensitivity analyses it can be concluded that if the model is to operate well, close attention needs to be paid to the drainage and storage capacities of the domi- nant surface types. Wetness sensors, located on a number of representative surfaces, appear to be a simple way of assess- ing the performance of the model, since measurements of Q E are not frequently available. If the model is applied in other hydroclimatological re- gimes, for example, in a tropical city, it would be advisable to conduct evaporation measurements so that the surface resistance submodel is appropriate. The robustness of the surface resistance submodel has not been assessed else- where. In the tropics and other locations where intense rainfall events occur, further investigation of the Massman [1983] drainage functions, which allow for the inclusion of the intensity of rainfall, may be worthwhile. Snow has been ignored, but its incorporation via an additional submodel should be reasonably easy with the appropriate melt equa- tions. The presence of frozen ground would require further adaptations ince the model would have to keep a record of water phases as well as depth within a store. We conclude that this physically based model is a prom- ising approach to model evapotranspiration and surface water state in urban areas. It has the advantage'of being 1754 GRIMMOND AND OKE; EVAPOTRANSPIRATION-INTERCEPTION MODEL dynamically responsive to both meteorological and seasonal conditions. The success of the one layer form of the model suggests that developments to include a multilayer version are worthwhile. The performance of the whole-canopy ap- proach to modeling the surface resistance is encouraging. The statistics indicate that it matches the results for forests, but the robustness of the parameters need to be tested in other urban areas. Acknowledgments. Special thanks are due to C. J. Souch, H. A. Cleugh, and H. P. Schmid, who helped in the field and for useful discussion. The field site was made available by B.C. Hydro. This work was supported by the Canada Department of Environment (Atmospheric Environment Service) and the Natural Sciences and Engineering Research Council of Canada through grants to T.O. S.G. was the recipient of a University of British Columbia Graduate Fellowship. REFERENCES Alley, W. M., D. R. Dawdy, and J. C. Schaake, Jr., Parametric- deterministic urban watershed model, d. Hy&aul. Div. Am. Sot'. Civ. Eng., 106(HY5), 679-690, 1980. Brater, E. F., Steps toward a better understanding of urban runoff processes, Water Resour. Res., 4, 335-347, 1968. Brutsaert, W., Evaporation Into the Atmosphere, 299 pp., D. Reidel, Hingham, Mass., I982. Calder, I. R., A model of transpiration and interception loss from a Spruce forest in Plynlimon, Central Wales, J. Hydrol., 33, 247- 265, 1977. Calder, I. R., and I. R. Wright, Gamma ray attenuation studies of interception from Sitka Spruce: Some evidence for a traditional transport mechanism, Water Resour. Res., 22,409-417, 1985. Chandler, T.J., Urban climatology and its relevance to urban design, WMO Tech. Note 149, 61 pp., World Meteorol. Organ., Geneva, 1976. Cleugh, H.A., and T. R. Oke, Suburban-rural energy balance comparisons in summer for Vancouver, B.C., Boundary Layer Meteorol., 36, 351-369, 1986. Davies, H., and T. Hollis, Measurements of rainfall-runoff volume relationships and water balance for roofs and roads, paper pre- sented at the Second International Conference on Urban Storm Drainage, Urbana, Illinois, 1981. Denis, J. E., Jr., D. M. Gray, and R. E. Welch, Algorithm 573 NL2SOL, an adaptive non-linear least squares algorithm, Assoc. Cornput. Mach. Trans. Math. Software, 7, 369-383, 1981. Dolman, A. J., J. B. Stewart, and J. D. Cooper, Predicting forest transpiration from climatological data, Agric. For. Meteorol., 42, 339-353, 1988. Dyer, A. J., A review of flux-profile relationships, Bo•tnda•n d Layer Meteorol., 7, 363-372, 1974. Falk, J., and J. Niemczynowicz, Characteristics of the above ground runoff in sewered catchments, in Urban Storm Drainage, edited by P. R. Helliwell, pp. 159-171, Pentech, London, 1978. Gash, J. H. C., An analytical model of rainfall interception by forests, Q. J. R. Meteorol. Soc., 105, 43-55, 1979. Gash, J. H., and A. J. Morton, An application of the Rutter model to the estimation of the interception loss from Thetford forest, J. Hydrol., 38, 49-58, 1978. Gash, J. H., W. J. Shuttleworth, C. R. Lloyd, J. C. Andre, J.P. Goutorbe, and J. Gelpe, Micrometeorological measurements in Les Landes forest during HAPEX-MOBILHY, Agric. For. Me- teorol., 46, 131-147, 1989. Gay, D. M., Remark on algorithm 573, Assoc. Cornput. Mach. Trans. Math. Software, 9, 139, 1983. Grimmond, C. S. B., The suburban water balance: Daily, monthly and annual results from Vancouver, B.C., M.Sc. thesis, 172 pp., Univ. of British Columbia, Vancouver, 1983. Grimmond, C. S. B., An evapotranspiration model for urban areas, Ph.D. thesis, 206 pp., Univ. of British Columbia, Vancouver, 1988. Grimmond, C. S. B., and T. R. Oke, Urban Water Balance 2: Results from asuburb of Vancouver, B.C., Water Resour. Res., 22, 1397-1403, 1986. Grimmond, C. S. B., H. A. Cleugh, and T. R. Oke, An objective hysteresis model for urban storage heat flux and its Comparison with other schemes, Artnos. Environ., in press, 1991. Grynning, S. E., A. A.M. Holtslag, J. S. Irwin, and B. Siverstscn, Applied ispersion modelling based on scaling meteorological scaling parameters, Atmos. lz'nviron., 21, 79-89, 1987. Hall, M. J., Urban Hydrology, 299 pp., Elsevier, New York, 1984. Halldin, S., H. Grip, and K. Perttu, Model for energy exchange0f a pine forest canopy, inComparison qfForest Water and Energy Exchange Models, edited by S. Halldin, pp. 59-75, International Society of Ecological Modelling, Copenhagen, 1979. Jarvis, P. G., The interpretation of the variations in leaf water potential nd stomatal conductance found in canopies in the field, Philos. Trans. R. Soc. London, Set'. B., 273, 593-610, 1976. Jarvis, P. G., G. B. James, and j. J. Landsberg, Coniferous forest, in Vegetation a d Atmospht, re vol. 2, Case Studies, edited by J. L. Monteith, pp. 171-240, Academic, San Diego, Calif., 1976. Kalanda, B. D., T. R. Oke, and D. I,. Spittlehouse, Suburban energy balance estimates tk..•r Vancouver, B.C., using Bowen Ratio energy-balance approach, J. Appl. Meteorol., 19, 791-802, 1980. Kramer, P. J., and T. T. Kozlowski, Physiology el' Woody Plants, p. 64, Academic, San Diego, Calif., 1979. Lacy, R. E., Climate and building, building research establishment report, 185 pp., Dep. of the Environ., London, 1977. Landsberg, H. E., The Urban Climate, 275 pp., Academic, San Diego, Calif., 1981. Lazaro, T. R., Urban !lydrology: A Multidiscipline Approach, 249 pp., Ann Arbor Science, Michigan, 1979. I.,indroth, A., Canopy conductance of coniferous forests related to climate, Water Rr, sottr. Res., 21,297-304, 1985. Massman, W. J., The derivation and validation of a new model f0r the interception of rainfitll by forests, Agric. Meteor&, 28, 261-286, 1983. McCaughey, J. H., D. W. Mullins, and M. Publicover, Comparative performance of two reversing Bowen Ratio measurement sys. tems, J. Atmox. Oceanic 7k, chm)l., 4, 724-730, 1987. Mcintosh, D. H., and A. S. Thom, Essentials q/' Meterology, 239 pp., Wykeham, London, 1972. Mechler, W. A., and D. G. Riecken, Review of methods for com- puting storm runoff from urban areas: History and state of the art, 57 pp., Greater Vancouver Sewerage and Drainage District, British Columbia, 1977. Monteith, J. L., Evaporation and environment, Syrup. Soc. Exp. Biol., 19, 205-224, 1965. Moore, C., Double-precision routines NL2SOL and NL2SN0, Tech. Note 214, 17 pp., Univ. of British Columbia Cornput. Centre, Vancouver, 1986. Mulder, J.P. M., Simulating interception loss using standard mete- orological data, in The Forest-Atmosphere Interaction, edited by B.A. Hutchinson and B. B. Hicks, pp. 177-196, D. Reidel, Hingham, Mass., 1985. Nash, J. E., and J. V. Sutcliffe, River flow forecasting through conceptual models, part I, A discussion ofprinciples, J.Hydrot., I0, 282-290, 1970. Nouh, M., Effect of model calibration on the least-cost design of storm-water drainage systems, in Urban Drainage Modeling, International Symposium on Comparison of Urban Drainage Models With Real Catchment Data, '86, Dubrovnik, Yugoslavia, edited by C. Maksimovic and M. Radojkovic, pp. 61-73, Perga- mon, New York, 1986. Oke, T. R., Methods in urban climatology, inApplied ClimatologY, Zurcher Geog. Schriften, vol. 14, edited by W. Kirchhofer, ^. Ohmura, nd W. Wanner, pp. 19-29, Eidgen6ssische Technishe Hochschule, Zurich, 1984. Oke, T. R., The urban energy balance, Prog. Phys. Geog., 12, 471-508, 1988. Oke, T. R., and H. A. Cleugh, Urban heat storage derived as energy budget residuals, Boundary Layer Meteorol., 39, 233-245, 1987. Penman, H. L., Natural evaporation from open water, bare soil and grass, Proc. R. Soc. London, Set. A., 193, 120-145, 1948. Pratt, C.J., J.J. Harrison, and J. R. W. Adams, Storm runoff simulation n runoff quality investigations, paper presented at he GRIMMOND AND OKE: EVAPOTRANSPIRATION-INTERCEPTION MODEL 1755 Third International Conference on Urban Storm Drainage, Gote- borg, Sweden, June 4-8, 1984. Randerson, D., Atmospheric s ience and power production, 850 pp., Office of Sci. and Tech. Inf., U.S. Dep. of Energy, Wash- ington, D.C., !984. RipIcy, E. A., and R. E. Redmann, Grassland, in Vegetation a d the Atmosphere, vol. 2, Case Studies, edited by J. L. Monteith, pp. 351-398, Academic, San Diego, Calif., 1975. Ross, S. L., and T. R. Oke, Tests of three urban energy balance models, Boundary Layer Meteorol., 44, 73-96, 1988. R0th, M., T. R. Oke, and D. G. Steyn, Velocity and temperature spectra and cospectra in an unstable suburban atmosphere, Boundary Layer Meteorol., 47, 309-320, 1989. Rutter, A. J., K. A. Kershaw, P. C. Robins, and A. J. Morton, A predictive model of rainfall interception i  forests, I Derivation of the model from observations in a plantation of Corsician pine, Agric. Meteorol., 9, 367-384, 1971. Rutter, A. J., A. J. Morton, and P. C. Robins, A predictive model of rainfall interception in forests, II, Generalisation of the model and comparison with observations in some coniferous and hard- wood stands, J. Appl. Ecol., 12,367-380, 1975. Schmid, H. P., and T. R. Oke, A model to estimate the source area contributing to surface layer turbulence at a point over a patchy surface, Q. J. R. Meteorol. Sot., 116(494), 965-988, 1990. Schmid, H. P., H. A. Clengh, C. S. B. Grimmond, and T. R. Oke, Spatial variability of energy fluxes in suburban terrain, Boundary Layer Meteorol., 54, 249-276, 1991. Sellers, P. J., and J. G. Lockwood, A computer simulation of the effects of differing crop types on the water balance of small catchments over long time periods, Q. d. R. Meteorol. Sot'., 107, 395-414, 1981. Sellers, P. J., Y. Mintz, and A. Dalchen, The design of the simple biosphere model (SiB) for use in general circulation models, J. Atmos. Sci., 43,505-531, 1986. Shuttleworth, W. J., Experimental evidence for failure of the Pen- man-Monteith equation in partially wet conditions, Boundary Layer Meteorol., !0, 91-94, 1976. Shuttleworth, W. J., A simplified one-dimensional theoretical de- scription of the vegetation-atmosphere interaction, Boundary Layer Meteorol., I4, 3-27, 1978. Shuttleworth, W. J., Evaporation, Rep. 56, 61 pp., Inst. of Hydrol., Wallingford, U K, 1979. Shuttleworth, W. J., Evaporation models in the global water bud- get, in Variations in the Global Water Budget, edited by A. Street-Perrott, M. Beran, and R. Radcliffe, pp. 147-171, D. Reidel, Hingham, Mass., 1983. Shuttleworth, W. J., Evaporation from Amazonian rainforest, Proc. R. Soc. London, Set. B, 323,321-346, 1988. Shuttleworth, W. J., Micrometeorology of temperate and tropical forests, Philos. Trans. R. Soc. London, Ser. B, 324, 199-224, 1989. Sman, H. T., C. C. D. F. van Ree, and M. Loxham, Water balance techniques applied to industrial complexes, paper presented atthe International Symposium on Hydrological Processes and Water Management i  Urban Areas, UNESCO, Duisburg, Federal Re- public of Germany, April 24-29, 1988. Stewart, J. B., Measurement a d prediction ofevaporation from forested and agricultural catchments, Agric. For. Manage., 8, 1-28, 1984. Stewart, J.B., Modelling surface conductance of pine forest, Agric. For. Meteorol., 43, 19-35, 1988. Stewart, J. B., and H. A. R. de Bruin, Preliminary studies of the dependence of surface conductance of Thefiord forest on envi- ronmental conditions, in The Forest-Atrnosphere Interaction, edited by B. A. Hutchinson and B. B. Hicks, pp. 91-104, D. Reidel, Hingham, Mass., 1984. Street, H. E., and H. Opik, The Physiology of Flowering Plants: Their Growth and Development, 3rd ed., 297 pp., Edward Arnold, Baltimore, Md., 1984. Thom, A. S., Momentum, mass and heat exchange of vegetation, Q. J. R. Meteorol. Soc., 99, 154-!70, 1972. Thom, A. S., and H. R. Oliver, On Penman's equation for estimat- ing regional evaporation, Q. J. R. Meteorol. Sot., 103,345-357, 1977. van den Ven, F. M. H., Urban hydrological cycle, paper presented at the International Symposium on Hydrological Processes and Water Management in Urban Areas, UNESCO, Duisburg, Fed- eral Republic of Germany, April 24-29, 1988. van Ulden, A.P., and A. A.M. Holtslag, Estimation of atmo- spheric boundary layer parameters for diffusion applications, j. Clim. Appl. Meteorol., 24, 1196-1207, 1985. Weiss, A., and D. L. Lukens, An electronic circuit for the detection of leaf wetness and a comparison of 2 sensors under field conditions, Plant Dis., 65, 41-43, 1981. WenzeI1, H. G., and M. L. Voorhees, Adaptation of ILLUDAS for continuous imulation, J. Hydraul. Div. Am. Soc. Civ. Eng., 106(HY 11), 1795-1812, 1980. Willmort, C. J., On the evaluation of model performance inphysical geography, in Spatial Statistics and Models, edited by G. L. Gaile and C. J. Willmott, pp. 443-460, D. Reidel, Hingham, Mass., 1984. Wright-McLaughlin Engineers Limited, Urban Storm Drainage Criteria Manual, 2 vols., Denver Regional Council of Govern- ments, Colo., 1969. Zinke, P., Forest interception studies in the U.S., in Forest Hydrol- ogy, edited by W. E., Sopper and H. W. Lull, pp. 137-161, Pergamon, New York, 1967. C. S. B. Grimmond, Climate and Meteorology Program, Depart- ment of Geography, Indiana University, Bloomington, IN 47405. T.R. Oke, Atmospheric Science Programme, Department of Geography, University of British Columbia, Vancouver, B.C., Canada V6T ! WS. (Received September 10, 1990; revised February 11, 1991; accepted February 15, 1991.)


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items