Open Collections

UBC Faculty Research and Publications

Mesoscale Analysis Method for Surface Potential Temperature in Mountainous and Coastal Terrain Deng, Xingxiu; Stull, Roland B. 2005-02-28

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata


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

Full Text

A Mesoscale Analysis Method for Surface Potential Temperature in Mountainous andCoastal TerrainXINGXIU DENG AND ROLAND STULLUniversity of British Columbia, Vancouver, British Columbia, Canada(Manuscript received 20 December 2003, in final form 19 July 2004)ABSTRACTA technique is developed to anisotropically spread surface observations in steep valleys. The goal is tocreate an improved objective analysis for the lowest, terrain-following numerical weather prediction (NWP)model level in mountainous terrain.The method is a mother–daughter (MD) approach, where the amount of information transferred fromone grid point (the mother) to all neighboring grid points (the daughters) depends on elevation differences.The daughters become mothers and further share information with their neighboring grid points. Thisiterative method allows information to follow valleys around ridges, while reducing spread over the ridgetop. The method is further refined to account for land–sea anisotropy.This approach is tested in the objective analyses of surface potential temperatures over the steep moun-tainous and coastal terrain of southwestern British Columbia, Canada. Analysis results are compared withother existing schemes using the Advanced Regional Prediction System Data Assimilation System (ADAS).It is found that the MD approach outperforms the other schemes over mountainous and coastal terrain.1. IntroductionObjective analysis transforms information from ran-domly spaced observing sites into data at regularlyspaced grid points (Krishnamurti and Bounoua 1996).It is often done by optimally combining observationsand a short-range NWP forecast, called a first guess(background). Analysis schemes, such as optimum in-terpolation (OI) (Daley 1991), successive corrections(Bratseth 1986), and 3D variational analysis (Larocheet al. 1999), require the specification of a covariancefunction to model the spatial correlations of back-ground errors. The background error covariance, usedto define the optimum weights for data spreading, isoften assumed to be isotropic.Isotropy, however, is a questionable assumption formesoscale analyses or for analyses in mountainous/coastal regions. Therefore, flow-dependent or anisotro-pic covariance models are receiving more attention,such as by introducing isentropic (Benjamin et al. 1991)or semigeostrophic (Desroziers 1997) coordinate trans-formations in analysis increments. Flow-dependentbackground error covariance can also be estimated ex-plicitly by Kalman filtering (KF), implicitly by 4Dvariational data assimilation (4DVAR) methods, or ap-proximated in the ensemble Kalman filter (EnKF)(Gauthier et al. 1993; Bouttier 1993; Thépaut et al.1996; Houtekamer and Mitchell 1998, 2001). But KF,4DVAR, and EnKF are computationally time consum-ing. Lanzinger and Steinacker (1990) incorporatedanisotropic correlation structures due to orography, toresolve sharp gradients across major mountain ridges(Alps, Pyrrenees). Miller and Benjamin (1992) used el-evation differences in correlation functions to accountfor terrain effects.Present computing power allows fine-grid NWPmodels to resolve small mesoscale flows within indi-vidual valleys. Such resolution is important, because insteep mountainous terrain the valleys contain most ofthe population centers, industries, and transportationroutes. To this end, the Mesoscale Compressible Com-munity (MC2) model (Benoit et al. 1997; Laprise et al.1997) has been running daily for several years at theUniversity of British Columbia (UBC). MC2 is a fullycompressible model using semi-implicit semi-Lagrangian techniques. Currently, the daily forecasts atUBC start with 108-km grid spacing, initialized fromthe Eta analysis, and then one-way nest down to 36, 12,4, and 2 km. The two finest grids can resolve differentweather in separate neighboring valleys. Accurate fine-resolution forecasts depend on accurate and represen-tative initial fields from which to start. This study fo-cuses on the surface data analysis at the lowest, terrain-following MC2 model level.Corresponding author address: Xingxiu Deng, Dept. of Earthand Ocean Sciences, University of British Columbia, 6339 StoresRd., Vancouver, BC V6T 1Z4, Canada.E-mail: xdeng@eos.ubc.caFEBRUARY 2005 DENG AND STULL 389© 2005 American Meteorological SocietyMWR2859The Advanced Regional Prediction System (ARPS)(Xue et al. 2000) Data Assimilation System (ADAS)(CAPS 1995; Brewster 1996) in ARPS5.0.0 Beta8 isused as an analysis tool. The ADAS employs the suc-cessive-correction method of Bratseth (1986), requiringless computational time than other advanced methods(i.e., KF and 4DVAR). This allows near-real-time,high-resolution analysis over most of British Columbia(BC), Canada, using moderately powerful computers.In the mountainous terrain of western NorthAmerica, many surface observation sites are located indeep valleys. Both the background error covarianceand representativeness of surface observations are af-fected by heterogeneous terrain. BC has a variety oflandscapes from seashores, straits, and fjords to moun-tains and valleys. Some valleys are narrow and longwith kinks and twists. The need for weather data atmany sites in a valley is often not matched by the avail-ability of weather stations. Typically, villages are strungalong valley floors, with ski resorts higher up the slopesand reservoirs in tributary valleys. Observations from asurface weather station at one village in the valley isoften not representative of conditions up- or downval-ley nor is it representative of conditions on nearbymountain slopes. Similarly, weather observations atridge top would not be representative of surfaceweather in the adjacent valley.This typical scenario stimulates the following ques-tion. At locations (in that one valley) distant from theobservation station, which of the following better rep-resents the initial air state: 1) the first guess (from aprevious high-resolution forecast) or 2) the observa-tions from the distant station? An experiment to ad-dress this question would be to densely instrument thevalley, and to determine the correlations. However, it islikely that the conclusions would apply only to that onevalley. It might also be risky to apply results from otherdensely instrumented valleys (e.g., the Mesoscale Al-pine Project) to the valleys in BC.Within the constraints of the present numerical ex-periments, we hypothesize that in a serpentine valley,the first guess (from a fine-resolution NWP model) isbetter than a distant observation from the same valley,in the absence of dense intravalley observations. Thishypothesis herein is referred to as the intravalley decor-relation assumption. An extension of this hypothesistreats two valleys separated by a high ridge and hereinis referred to as the intervalley decorrelation assump-tion. If there is an observation in only one valley, thenthe first guess in the noninstrumented valley representsthe air state better than the observation from the othervalley.We implement the intravalley and intervalley decor-relation hypotheses within the ADAS Bratseth schemefor data assimilation in complex terrain. The intravalleyassumption uses a Gaussian drop-off with distance fromthe observation, where distance is measured along thecircuitous path of the valley. The intervalley assump-tion uses an anisotropic term that reduces data spread-ing into terrain of differing elevation.A mother–daughter (MD) approach is proposed hereto calculate the anisotropic term and circuitous traveldistance (CTD). The anisotropic effect is introducedvia a “sharing factor,” namely, the fraction of informa-tion shared between an observation and analysis gridpoint (GP). Near coastlines, land–sea anisotropy is alsoincluded.Section 2 describes the ADAS, and modifications forthis study. Section 3 reviews basic concepts of atmo-spheric boundary layers (BLs) in mountainous terrain.The MD approach is presented in section 4. Numericalexperiments are given in section 5. In section 6, withcase examples we demonstrate characteristics and ad-vantages of the MD approach over two other schemes.Computational costs are also examined. Summary anddiscussions are in section 7, and a sensitivity analysis isin the appendix.2. The analysis toola. General descriptionThe ADAS was developed by the Center for Analy-sis and Prediction of Storms (CAPS, University ofOklahoma). It is a 3D mesoscale analysis system thatingests and analyzes meteorological data from differentobservational sources. The ADAS is efficient, easilyimplemented and attractive for high-resolution analysison model levels.The ADAS combines observations with a first guessusing the Bratseth (1986) scheme, which converges toOI. Like OI, the Bratseth scheme accounts for the rela-tive error between the observations and the first guess.Different from OI, this scheme iterates without the in-version of large matrices, and requires less computerresources.The ADAS has gained wide usage in research andoperational forecasting. Research applications for astrong cold-front passage (Ciliberti et al. 1999) and sen-sitivity experiments (Ciliberti et al. 2000) showed thatrealistic local and mesoscale structures can be analyzedwith the aid of local data. Operational applications in-cluded short-range weather forecasts over east-centralFlorida (Case et al. 2002), weather support for the 2002Winter Olympics (Horel et al. 2002), and near-real-timeapplications to complex terrain in the western UnitedStates (Lazarus et al. 2002). Lazarus et al. (2002) intro-duced into the ADAS a terrain factor for surface data.1The terrain factor affects only the analysis in the freeatmosphere.1This was expanded for 3D data (surface and sounding data) byKeith Brewster (CAPS 2004, personal communication) in a newerversion (ARPS5.1.0).390 MONTHLY WEATHER REVIEW VOLUME 133b. Modifications1) ADAS VERSUS MC2 GRIDSFor better analysis quality, analyses are performeddirectly on the MC2 grid using MC2 output as first-guess fields.2) QUALITY CONTROLFor surface observations, quality control includes cli-matological checks, temporal checks, horizontal consis-tency checks [Barnes (1964) technique], and observa-tion-to-background difference checks.For case studies with virtual observations (section6a), no quality control is needed. For case studies withreal observations (section 6b), the horizontal consis-tency checks and observation-to-background differencechecks are turned off for the following reasons. TheBarnes technique does not work well for quality controlof fine-resolution data in complex terrain, because anyone surface observation might not be represented wellby nearby observations at different elevations or on theopposite side of a mountain range. In steep terrain, theobservation-to-background difference might be largedue to a discrepancy between the actual station eleva-tion and the model GP elevation, rather than due tobad quality of the surface observation. We have re-duced the ADAS trilinear interpolation to bilinearwhen the first-guess fields are interpolated to the ob-servation locations. The ADAS trilinear interpolationalgorithm causes large interpolation errors by introduc-ing “free-atmosphere” background values (Lazarus etal. 2002); also, the algorithm causes large extrapolationerrors for the surface valley stations far below the low-est model level.Because the objective quality-control methods areproblematic in steep mountains, we instead do subjec-tive examinations on observation minus background re-siduals (not shown).3) ANISOTROPIC CORRELATION FUNCTIONThe ADAS Bratseth scheme works iteratively, up-dating the first guess using the difference between theobserved values and observational estimates derivedfrom the analysis (CAPS 1995; Lazarus et al. 2002). Theoptimum Bratseth weights rely on the correlation (H9267)ofthe background errors (and on the ratio of observationerror variance to background error variance). Over flatterrain, a Gaussian function [G(d)] is frequently used tosmoothly reduce the correlation with distance (d), caus-ing the impact of a single observation to be isotropic inthe grid space around the observation [H9267(d) H11008 G(d)].In mountainous terrain, we modify the Gaussiandrop-off to be a function of the CTD (s, defined insection 4) from the observation to an analysis point[H9267(s) H11008 G(s)]. As discussed in section 1, the correlationshould be weighted by an anisotropic term (S) that is afunction of the difference between elevations of theobservation (Zo) and the analysis point (Za). By sepa-rating horizontal and vertical effects into two factors, ananisotropic correlation becomesH9267H20849sH20850H11005 GH20849sH20850H11003 SH20849Zo, ZaH20850. H208491H20850An MD approach is proposed in section 4 to generateS and s for each observation station. This approach isbased on first-order BL characteristics in mountainousterrain, described next.3. Topographic anisotropy and boundary layersa. First-order boundary layer characteristics inmountainous terrainAtmospheric BLs always contain air of colder poten-tial temperature than the air higher in the free atmo-sphere (FA), because the average troposphere is stablystratified (Stull 1988, 2000). Between the BL and FA isa strongly stable layer that caps the BL. Trapped in theBL below this cap are pollutants, humidity, and heatreleased from the surface.In steep mountains, the ridge heights above valleyfloors are frequently the same order as BL heights. Toa zeroth order, the turbulent and advective communi-cation of physical, dynamic, and chemical states be-tween one valley and the other depends greatly onwhether the BL is shallower than the surroundingridges. For shallow BL situations, tracers emitted intoone valley BL are unlikely to reach neighboring valleys,due to blocking by mountains. For situations with BLsdeeper than ridges, BL air can mix between neighbor-ing valleys.To a first order, the BL top (zi) is often not level overcomplex terrain. Observations show a variety of behav-iors, ranging from BLs that follow the topography toBLs that seem relatively level (Lenschow et al. 1979;De Wekker et al. 1997; Kalthoff et al. 1998; Kossmannet al. 1998). Stull (1992) identified four archetypicalBL-top characteristics (Fig. 1). Gravitational forcestend to make a level BL top (Fig. 1c), much like anocean surface that is level over undersea mounts. Otherfactors including entrainment, advection, and frictiontend to make a terrain-following top (Fig. 1b). See Stull(1998) for a similarity-theory analysis of the competingfactors. We will focus on the level and terrain-followingcases (Fig. 2).BL parameterizations in most NWP models do notprovide information on levelness of the BL top withinfirst-guess fields for use during data assimilation. Norare observations of BL depth usually available at everysurface weather station. In an attempt to compensatefor the lack of information, we develop a backgrounderror correlation (H9267) that includes weighting factors forboth “level-top” and terrain-following BL effects.Namely, we define the sharing factor SF [S in (1)] dur-ing one iteration byFEBRUARY 2005 DENG AND STULL 391SodH20849H9263H11001 1H20850H11005 SomH20849H9263H20850H11003 WTFBLH20849H9263H11001 1, H9263H20850H11003 WLTBLH20849H9263H11001 1H20850, H208492H20850where subscript o represents the observation, subscriptsm and d represent a GP (mother) and its neighboringGP (daughter), respectively. Weights W represent theportions of the SF associated with terrain-followingBLs (TFBLs) or level-top BLs (LTBLs).b. Parameterization of terrain-following andlevel-top BL contributionsLet Z be elevation. Approximate the terrain-following BL effects on the SF, from one mother GPduring iteration H9263 to a neighboring daughter GP at it-eration H9263 H11001 1, byWTFBLH20849H9263H11001 1, H9263H20850H11005 1 H11002H20875| ZmH20849H9263H20850H11002 ZdH20849H9263H11001 1H20850|zref1H20876aH208493H20850for | Zm(H9263) H11002 Zd(H9263H11001 1) | H11021 zref1; otherwise WTFBLH11005 0.The terrain-following BL-depth parameter is zref1,while parameter a controls the analysis decorrelationrate.For illustration, suppose that zref1 H11005 1kmanda H11005 1.Let the weather station be collocated with a model GP(the first mother GP) at the base of linearly slopingterrain (Fig. 2a) such that each successive daughter GPis higher up the slope: z(H9263) H11005 0, 0.5, 1.0, 1.5 km . . . forH9263 H11005 0, 1, 2, 3....The iteration counter (H9263) coincideswith each successive daughter GP. The resulting se-quence of weights from (3) is WTFBLH11005 1.0, 0.5, 0.5, and0.5. These weights accumulate multiplicatively via theiteration given by (2). Considering only terrain-following effects, S(0) H11005 1, S(1) H11005 0.5, S(2) H11005 0.25, andS(3) H11005 0.125. The GP at Z H11005 1.5 km has nonzero SF, forthe case of a terrain-following 1-km-thick BL, whereturbulence can mix valley variables up the mountainslope.Similarly, we can define a level-top BL weight byWLTBLH20849H9263H11001 1H20850H11005 1 H11002H20875| ZoH11002 ZdH20849H9263H11001 1H20850|zref2H20876bH208494H20850for | ZoH11002 Zd(H9263 H11001 1) | H11021 zref2; otherwise WLTBLH11005 0.The level-top BL-depth parameter is zref2, while pa-rameter b controls the analysis decorrelation rate.For illustration with the same sloping GPs as in theprevious example, suppose that zref2 H11005 1kmandb H110051. The resulting sequence of weights from (4) is WLTBLH11005 1.0, 0.5, 0., and 0. These weights accumulate multi-plicatively via iteration as in (2). If we neglect the ter-rain-following effects, S(0) H11005 1, S(1) H11005 0.5, S(2) H11005 0.,and S(3) H11005 0. The GP at Z H11005 1.5 km has zero SF,because it is above the level BL top (Fig. 2b). At allFIG. 1. Archetypical categories of mixed-layer-top levelness (af-ter Stull 1992). (a) Hyper terrain following; (b) terrain following;(c) level; (d) contra-terrain following.FIG. 2. Shallow boundary layers in steep terrain, for (a) a ter-rain-following boundary layer and (b) a level-top boundary layerblocked by mountains. White dots represent surface grid points.The iteration counter H9263 coincides with each successive daughtergrid point. The lowest grid point at H9263 H11005 0 is collocated with thesurface weather station.392 MONTHLY WEATHER REVIEW VOLUME 133higher-elevation GPs, the SF is also zero. Those pointshigher on the mountain are assumed to be in a differentair mass than the valley station. It is worth noting thatzref2 is essentially the BL depth at the observation lo-cation. It will be shown later in the appendix that analy-sis results are far more sensitive to zref2 than to zref1.These illustrations, albeit somewhat primitive, repre-sent basic BL characteristics in a way that allows aniso-tropic correlations to be created, as shown next.4. The mother–daughter approacha. General conceptsWe first consider a surface observation that happensto lie directly on a model GP. The fraction of informa-tion shared between the observation and any analysispoint is called “sharing factor” (SF), which ranges be-tween 0 and 1. Let subscripts o, m, and d represent anobservation, a mother GP and a neighboring daughterGP, respectively. The iterative MD approach is ex-pressed bySodH20849H9263H11001 1H20850H11005 SomH20849H9263H20850H208771 H11002H20875| ZmH20849H9263H20850H11002 ZdH20849H9263H11001 1H20850|zref1H20876aH20878H11003H208771 H11002H20875| ZoH11002 ZdH20849H9263H11001 1H20850|zref2H20876bH20878if | ZmH20849H9263H20850H11002 ZdH20849H9263H11001 1H20850| H11021 zref1and | ZoH11002 ZdH20849H9263H11001 1H20850| H11021 zref2 H208495aH20850SodH20849H9263H11001 1H20850H11005 0if| ZmH20849H9263H20850H11002 ZdH20849H9263H11001 1H20850| H11350 zref1or | ZoH11002 ZdH20849H9263H11001 1H20850| H11350 zref2, H208495bH20850where H9263 is iteration counter; Z is elevation; and a, b,zref1, and zref2 are free parameters. The amount ofinformation that a mother transfers to each daughter isreduced by the elevation difference between them. Theiteration starts with Som(0) H11005 1 at the surface observa-tion location. As previously discussed, the first and sec-ond bracketed terms in Eq. (5a) contain the terrain-following [WTFBLin Eq. (3)] and level-top [WLTBLinEq. (4)] BL effects, respectively.The observation is treated as the first mother, whohas daughters at each of her eight neighboring GPs.This is the first iteration (H9263 H11005 1). These daughter SFsare calculated by applying (5) to every mother–daughter pair. The circuitous travel distance from theobservation to a neighboring daughter GP, after H9263 H11005 1,is equal to the straight-line horizontal distance if the SFat that daughter GP is not zero. For zero SF, the CTDis assigned a large value (s H11005 1.0e5 km), so that H9267(s)in(1) is zero due to the Gaussian drop-off. Namely, theobservation has no influence on that GP. Next, everydaughter becomes a mother and has her own eightdaughters around her, and so on (see below). Aroundany daughter might be several mothers who contributedifferent SFs to that daughter, but only the biggest SF(among all values from neighboring mothers and thevalue from the previous iteration) is kept. Thus, Eq. (5)has to be applied successively to every GP within agradually enlarged subdomain centered at the observa-tion location. This process is computationally expen-sive, as demonstrated in section 6b. From the seconditeration onward, the CTD from the observation to anydaughter GP is found by accumulating the incrementalpath segments from each successive mother to thedaughter GP. Usually there are many different pathsfrom the observation to a distant GP, but the pathsaved is the one with the largest SF at that daughter GP.Figure 3 illustrates the CTD (s) from the observation“O” to GP “D.”To illustrate calculations of the SFs and CTDs, con-sider a surface station that happens to lie directly on amodel GP at (x, y) H11005 (7, 5) (Figs. 4 and 5, indicated bya solid triangle). Horizontal grid spacing is 3 km. Ini-tially, the SF is 1.0 at the observation location, and zeroelsewhere (left-hand panel of Fig. 4a). The CTD is ini-tially 0 km at the observation point and 1.0e5 km oth-erwise (right-hand panel of Fig. 4a). The SFs and CTDsafter the first iteration are shown in Fig. 4b. The CTDfor a daughter GP remains 1.0e5 km if its SF is 0. Thus,the CTD (s) from the observation to the GP at (8, 5) is1.0e5 km, even though the straight-line horizontal dis-tance (d) is 3 km. A daughter GP with nonzero SF hasfinite CTD [e.g., s H11005 d H11015 4.2 km for GP (6, 6)]. In thesecond iteration, each daughter becomes a mother andhas her own eight daughters around her, yielding 24(8H1100116) daughters after H9263 H11005 2 (Fig. 4c). This process isrepeated until | Sod(H9263H11001 1) H11002 Sod(H9263) | H11021 0.01 for all GPs.The final SFs (Sod) and CTDs (s) are used in (1) todetermine the anisotropic correlation function [H9267(s)].Figure 5 shows the end results for this hypothetical sta-FIG. 3. Schematic of the straight-line horizontal distance (d) andthe circuitous travel distance (s H11005 d1 H11001 d2 H11001 d3) from the obser-vation location O to the grid point D, after the third iteration. Thegrid point D received the most contribution in its sharing factorfrom the grid point B, while the grid point B received the mostcontribution from the grid point A that received the most contri-bution from the observation O. This illustrates a situation such aswhen a mountain (indicated by shading) lies between O and D.FEBRUARY 2005 DENG AND STULL 393FIG. 4. (left column) Printed values of the sharing factors (zeroes are not printed except when they firstappear in an iteration; darker shading indicates higher elevations) and (right column) circuitous traveldistance (km) (large values of 1.0e5 are not printed). (a) iteration counter H9263 H11005 0; (b) H9263 H11005 1; (c) H9263 H11005 2.The observation is indicated by a solid triangle at (x, y) H11005 (7, 5). Grid spacing is 3 km. The parametersused are a H11005 2, b H11005 2, zref1 H11005 750 m, and zref2 H11005 750 m. Grid cells along the outside edge are cut offin this plot but are really the same size as all other grid cells.394 MONTHLY WEATHER REVIEW VOLUME 133tion, using a Gaussian drop-off with a standard devia-tion of 30 km.To examine the subtleties of the CTD, consider GP(6, 5) in Fig. 4. After H9263H11005 1, S H11005 0.31 and s H11005 3 km (Fig.4b). Obviously, the SF comes from the first mother (i.e.,the observation). After H9263 H11005 2, S H11005 0.38 and s H11015 7.2 km(Fig. 4c). Now the largest SF comes from a second-generation mother GP (6, 4). The CTD (s H11015 7.2 km)equals the sum of d H11015 4.2 km between the observationand the second-generation mother GP (6, 4), plus d H110053 km from GP (6, 4) to the target GP (6, 5).Parameters a and b in (5) control SF reduction by themother–daughter GP elevation difference. As a and bare reduced from 4 toward 1, the terrain effects on theanalysis increase. The optimum values of a, b, zref1,and zref2 can be obtained by analyzing a subset of theavailable observations, and verifying against the re-maining observations for a substantial period of dailyanalyses, assuming a dense observation network.For fixed GPs and observation locations within amodel domain, the SFs and CTDs for each observationstation can be determined. In the case of a surface ob-servation station that is not collocated with a model GP,the station is first approximated as being collocated atthe nearest neighboring model GP with minimum el-evation difference between them. The resulting SFs andCTDs are saved as a fixed file, which can be appliedunchanged for each day’s analysis to efficiently deter-mine the correlation in the ADAS Bratseth scheme. Inthe presence of temporally missing observations, a flag(i.e., H11002999.) can be assigned to the missing observa-tions. The SFs and CTDs for a newly added surfacestation can be calculated and appended to the fixed file.b. Refinement for coastal terrainTo consider land–sea contrasts, Hessler (1984) esti-mated two correlation functions from surface tempera-ture statistics near the Baltic Sea. Between pairs of seastations or pairs of land stations, he approximated cor-relations by exp(H110020.08d0.13) cos(0.4d). Between sea andland stations, he used exp(H110020.29d0.06)cos(0.2d).We introduce a similar approach by including an ad-ditional factor (i.e., the third term in brackets on theright-hand side):SodLSH20849H9263H11001 1H20850H11005 SomLSH20849H9263H20850H208771 H11002H20875| ZmH20849H9263H20850H11002 ZdH20849H9263H11001 1H20850|zref1H20876aH20878H11003H208771 H11002H20875| ZoH11002 ZdH20849H9263H11001 1H20850|zref2H20876bH20878H11003H208771 H11002| LSoH11002 LSdH20849H9263H11001 1H20850|KLSH20878H208496H20850for | Zm(H9263) H11002 Zd(H9263 H11001 1) | H11021 zref1 and | ZoH11002 Zd(H9263 H11001 1) |H11021 zref2; otherwise SLSod(H9263 H11001 1) H11005 0. The land–sea maskLSoat the observation location is equal to one for aland observation and zero for a water observation. Theland–sea mask LSdat a daughter GP is similar, butaveraged over the nine GPs immediately around andincluding the daughter GP. The average allows gradualtransition of the SFs across the coastlines in an attemptto slightly account for sea/land breeze effects. Param-eter KLScontrols the degree of the decorrelation be-FIG. 5. (top) Sharing factor isopleths (thin lines, after a Gaus-sian drop-off with standard deviation of 30 km) and (bottom) thecircuitous travel distance (km) for an idealized observation invalley A within a horizontal domain that spans 36 km H11003 36 kmwith a grid increment of 3 km. Terrain elevations (m) are shownby shading in the top plot. Darker shading indicates higher eleva-tions, with a maximum terrain-height difference of 2000 m. Theobservation is indicated by a solid triangle at (x, y) H11005 (7, 5). Theparameters used are a H11005 2, b H11005 2, zref1 H11005 750 m, and zref2 H11005 750 m.FEBRUARY 2005 DENG AND STULL 395tween land stations and analysis points over water, orbetween water stations and analysis points over land;KLSwas subjectively set to 1.0 for the case studies insection 6b, but could be greater than 1.0 for strongsea-breeze events.Ideally, the location of the land–sea mask shouldshift as the sea/land breezes evolve every day; however,we neglect this shift. This allows a computationally ef-ficient constant land–sea mask that still captures first-order effects.5. Numerical experimentsFour assimilation experiments (Table 1) are per-formed for comparison purposes. GAUSS, the nativescheme in the ADAS (ARPS5.0.0 Beta8), includes ter-rain effects with a Gaussian function of elevation dif-ference (H9004z) between the GP and observation station,and horizontal Gaussian decay to define the region ofinfluence. In experiment TERR_DIFF, the Miller andBenjamin (1992) method is added to the ADAS, butslightly modified from the original correlation function2for potential temperature used in Miller and Benjamin(1992).In experiment MD, the SFs from Eq. (5) are used toreplace the Gaussian elevation term in GAUSS. MDfurther differs from GAUSS and TERR_DIFF by usingthe CTD, rather than the straight-line horizontal dis-tance, in the horizontal part of the correlation function.This effectively limits the region of influence in com-plex terrain and restricts valley observation informationto follow the valleys around the ridges, while reducingspread over the ridge top. Similarly, ridge-top informa-tion follows ridges, with little spread into valleys.MD_LSMG is the same as MD except that land–seaanisotropy is included in the SFs (section 4b).Unless stated specifically, the horizontal correlationlength scale (Rh) is 90 km for the first and second Brat-seth passes, 60 km for the third pass, and 30 km for thefourth and fifth passes.6. Case studiesThe utility of the MD approach is demonstrated be-low for surface potential-temperature analyses in tworegions. Virtual and real observations are used formountainous and coastal terrains, respectively, insouthwestern BC.a. Virtual observationsWe first use “fraternal-twin” experiments to evaluatethe MD approach in mountainous terrain. In fraternal-twin experiments, the “truth” and “analysis” modelsare similar, but not identical (Arnold and Dey 1986).The truth model (MC2 at 2-km grid spacing) is firstintegrated 12 h to generate a reference truth atmo-sphere. To generate a different first guess, the analysismodel (MC2 at 3 km with radiation turned off to pur-posely introduce errors) is integrated 12 h with thesame boundary conditions, but using perturbed initialconditions. Amplitudes of random perturbations fortemperature (u wind, v wind, logarithm of pressure per-turbation, specific humidity) are 5.0°C (5.0 kt, 5.0 kt,0.0005, 0.005 kg kgH110021). A dynamic initialization algo-rithm built into MC2 is activated for the 3-km run toremove spurious gravity waves excited by the addedperturbations.The fraternal-twin experiments are performed for7–8 March 2003, characterized by an Arctic outbreak,when cold shallow air masses from northern Canadaswept into the valleys in BC.The MC2 runs consist of five self-nested grids withgrid spacings of 108, 36, 12, 4, and 2 or 3 km (see the lastthree grids in Fig. 6). The Eta analysis and forecasts at90.7-km grid spacing from the National Centers for En-vironmental Prediction (NCEP), valid at 0000 UTC 7March 2003, are used as the initial and boundary con-ditions for the coarsest grid. The MC2 at 2 km is startedat 1200 UTC 7 March from MC2 4-km output and in-tegrated 12 h to generate the truth atmosphere, fromwhich virtual surface observations are extracted. Allvirtual observations are assumed to be perfect, withzero observation error. The MC2 at 3 km is also startedat 1200 UTC 7 March, but from the randomly per-turbed output of the 4-km run, and integrated 12 h to2exp(H11002| d |2/2r2)/(1 H11001 Kz· |H9004z|2), r H11005 300 km, KzH11005 7 H1100310H110027mH110022.TABLE 1. Experiments performed, where Rhand Rzare corre-lation length scales in the horizontal and vertical, respectively,which define the regions of influence; dijand H9004zijare straight-linehorizontal distance and elevation difference between an analysisgrid point and the observation station, respectively; Kzis a coef-ficient; sijis the circuitous travel distance from the observation toan analysis grid point determined in the mother–daughter pro-gram; Sodis the sharing factor at an analysis grid point, whichrepresents how much the analysis grid point shares the observa-tion information; SLSodis the same as Sodexcept that the former isfrom the refined mother–daughter approach. In this study, Rhis90 km for the first and second iterations, 60 km for the thirditeration, and 30 km for the fourth and fifth iterations, during theADAS assimilation cycle.Experiment Correlation function (H9267)GAUSS expH20873H1100212| dij|2Rh2H20874· expH20873H1100212|H9004zij|2Rz2H20874TERR_DIFF expH20873H1100212| dij|2Rh2H20874H11408H208491 H11001 Kz· |H9004zij|2H20850MD expH20873H1100212sij2Rh2H20874· SodMD_LSMG expH20873H1100212sij2Rh2H20874· SodLS396 MONTHLY WEATHER REVIEW VOLUME 133provide a first guess. The analyses are performed at thelowest terrain-following model level of the MC2 3 kmat 0000 UTC 8 March.Figure 7 shows six virtual surface stations. Stations o1(elevation 920 m) and o2 (elevation 536 m), located indifferent valleys, are used for analysis. Station o1 is inthe Elaho River Valley, and o2 is in the Lillooet RiverValley. Four stations (b1, b2, b3, and b4) in the LillooetRiver Valley are used for verification. These stationsare selected to evaluate the impact of observation o1 onthe neighboring Lillooet River Valley. Elevations forstations b1, b2, b3, and b4 are 894, 881, 765, and 830 m,respectively.1) IMPACT OF A SINGLE VALLEYSURFACE OBSERVATIONFor objective analysis, the interpolated first guess isusually subtracted from the observation to give an “ob-servation increment,” which is analyzed to produce“analysis increments” (AIs) at the GPs. The final analy-sis at each GP is then the first guess plus the analysisincrement. The AIs produced by o1 show how data arespread to surrounding GPs.Figure 8 shows the impact of o1 on the GP analysis.The observation increment is H110023.35 K. The Bratsethanalysis is confined to a single pass for convergence inthe presence of a single perfect observation. All experi-ments show that the observation mostly influences thevalley GPs.GAUSS and TERR_DIFF both yield large AIs invalleys disconnected from the Elaho River Valley, eventhough both valleys may not necessarily share the sameair mass. The Lillooet River Valley (containing b1, Fig.7) is such an example. In MD, the large AIs are in theElaho River Valley (containing o1), with much smallervalues in the Lillooet River Valley.2) IMPACT OF TWO SURFACE OBSERVATIONS FROMDIFFERENT VALLEYSThe analysis at one GP using observations from ad-jacent but disconnected valleys, where each valley con-tains a different meteorological regime, is a problem incomplex terrain. To examine this, consider two obser-vations in different valleys (stations o1 and o2 in Fig. 7).The observation at o1 is 3.35 K colder than the firstguess, while the observation at o2 is 3.36 K warmer.Four stations (b1, b2, b3, and b4) in the Lillooet RiverValley are used to verify the analyses.First, only the observation at o2 is used. To evaluatethe three schemes, the GP analyses are interpolated tothe four verifying stations (not collocated with any GP)by cubic polynomial interpolation. The root-mean-square error (rmse) between the observations andanalyses at those sites are calculated. For comparisonpurposes, the parameters for GAUSS (RzH11005 500 m),TERR_DIFF (KzH11005 7.0 H11003 10H110026mH110022), and MD (a H11005 2,b H11005 2, zref1 H11005 750 m, zref2 H11005 750 m) are chosen so thatthe three schemes produce similar rmses. The values ofthese parameters are kept fixed for all subsequent ex-periments.Rmses for GAUSS, TERR_DIFF, and MD are0.6179, 0.6468, and 0.6559 K, respectively. Rmses forthese three schemes are slightly greater than that be-tween the observations and the first guess, which isFIG. 6. MC2 grid domains for 4, 3, and 2 km. MC2 4-km outputprovides nesting files to initialize 2- and 3-km runs. MC2 2 km isrun to provide a reference atmosphere to generate virtual surfaceobservations for analysis and verification. MC2 3 km is run toprovide first-guess fields for analysis.FIG. 7. Virtual surface observation stations, indicated by solidtriangles, used in the analysis and verification. The stations arepositioned at the truth-model terrain height in the Coast Moun-tains north of Vancouver. Station o1 (50.326°N, 123.578°W) isnear the mouth of the Elaho River. Station o2 (50.3344°N,122.767°W) is near the town of Pemberton. Terrain elevations (m)are from the truth model. Darker shading indicates higher eleva-tions, with a maximum elevation difference of 2055 m in thisfigure.FEBRUARY 2005 DENG AND STULL 3970.5395 K. The analyses and the first guess versus theobservations for the three schemes are shown in theleft-hand panel of Fig. 9. It is evident that the threeschemes produce similar analyses.The observations from both sites (o1 and o2) arethen used to produce an analysis. While the same set ofobservations at b1, b2, b3, and b4 are used for verifi-cation. The additional observation at o1 has a negativeimpact on both the GAUSS and TERR_DIFF analyses(see the right-hand panel of Fig. 9). Rmse increases to2.1464 K for GAUSS and 2.4380 K for TERR_DIFF.However, in MD, the additional observation at o1 has aminor impact on the analysis at the verifying stations, asdesired because the two valleys with observations arenot strongly coupled. The slight influence from theadded observation reduces the rmse to 0.4807 K. Thermse with two observations is reduced about 11% fromthe first guess.3) IMPACT OF A SINGLE HIGH-ELEVATIONSURFACE OBSERVATIONIn complex terrain, most surface observation sites arelocated in densely populated valleys. There are, how-ever, some surface stations located at high elevations inthe mountains, such as at ski areas and in subalpineforest ecosystems. This subsection examines the analy-sis impact of high-elevation observations.Figure 10 shows the AIs produced by o3 (elevation1954 m). The observation increment is H110021.56 K. Allschemes produce minimal spread into valley GPs, withMD producing zero spread into the valleys.For GAUSS and TERR_DIFF, the observation in-crement is spread to all of the surrounding high moun-tains, including those separated by deep valleys. In MD,the observation corrects the first guess only for thoseGPs on the same mountain as o3. The GAUSS andTERR_DIFF analyses in Fig. 10 look more realisticthan MD for two reasons. During shallow cold-air pool-ing events, the high-mountain tops are all penetratinginto the free atmosphere, where they would experiencethe same weather. Second, during deep BL events withBL top above the mountain top, vigorous turbulencewould also cause mixing across the valleys. In futurework, we plan to modify MD to allow such spreading.b. Real observationsThe MD approach is also tested with real observa-tions in the coastal terrain of the Georgia Basin, BC.Hourly data are combined from several agencies, in-cluding the BC Ministry of Transportation and High-ways (MOTH), BC Ministry of Water Land and AirProtection (WLAP), Environment Canada (EC), andBC Hydro (HYDR). The higher density of observa-tions over the Georgia Basin (Fig. 11) allows us to with-hold some for verification (Table 2). The land–seamask is based on the model-resolved coastlines. Theratio of observation-error variance to background-errorvariance is 0.08 as in Miller and Benjamin (1992). Thelow observation-to-background error implies that theobservations are heavily weighted.A3–4 February 2003 case, characterized by a strongridge over the Georgia Basin, is used to examine ifweak land–sea thermal contrasts can be properly ana-lyzed. Figure 12 shows the evolution of potential tem-peratures for several surface stations on land and waterFIG. 8. Analysis increments (isopleths) for potential tempera-ture at the lowest model level in response to a single surfacepotential-temperature observation at station o1 indicated by asolid triangle. Contour interval is 0.3 K. Darker shading indicateshigher elevations, with a maximum elevation difference of 1940 m,which is less than that in Fig. 7 because of the smoothed terrainused in the NWP model. (top) GAUSS; (middle) TERR_DIFF;(bottom) MD.398 MONTHLY WEATHER REVIEW VOLUME 133during 0000–1200 UTC 4 February 2003. In late after-noon at 0000 UTC (1600 PST), all land stations exceptVancouver INTL ARPT were slightly warmer than thewater stations. After a transition period from 0100 to0300 UTC, all land stations were colder than the waterstations. Land–sea thermal contrast is the most promi-nent in early morning at 1200 UTC (0400 PST).The Eta analysis and forecasts from NCEP, valid at0000 UTC 3 February 2003, are used to drive the coars-est grid (108-km grid spacing), which in turn drivesgrids of 36, 12, and 4 or 3 km. MC2 at 3 km providesfirst-guess fields and is started at 1200 UTC 3 Februaryfrom MC2 4-km output. Analyses are performed everyhour from 0000 to 1200 UTC 4 February, by blendinghourly surface observations with the first guess at thelowest terrain-following model level valid at the sameFIG. 9. Analysis, indicated by a circle, and first guess, indicated by a plus, vs observation. (left column)Observation at o2 only is used in the analysis. (right column) Observations at both o1 and o2 are usedin the analysis. (top) GAUSS; (middle) TERR_DIFF; (bottom) MD. The four points (circles or pluses)correspond to the four verifying stations b1–b4.FEBRUARY 2005 DENG AND STULL 399time. The hourly analyses are not incorporated into theforecast cycle; thus, each analysis is independent of pastobservations. Eight stations are used for analysis and 11others for verification (Table 2), but the available num-ber of reporting stations varies with the analysis time(Table 3).Figure 13 shows time series of the first guess (FSTG)or the analyzed potential temperature and observationsfor one water station (Entrance Island) and one landstation (White Rock), separately. The 3-km run pro-duces a good first guess for White Rock from 0300 to1200 UTC 4 February, but underpredicts the potentialtemperatures for Entrance Island. The analyses fromall experiments show more agreement with the obser-vations than the first guess. Comparatively, MD andMD_LSMG both produce better analyses than GAUSSand TERR_DIFF. MD_LSMG is superior in maintain-ing the thermal contrast between Entrance Island andWhite Rock. Figure 14 further demonstrates advan-tages of MD_LSMG, based on separate averages fromthree land and three water stations. Averaging all landstations, including the mountain stations, makes theland–water contrast very small. However, advantagesof MD_LSMG can be seen from separate averages ofall water stations and all land stations except the moun-tain stations (not shown).The verification statistics (Table 4) show that allschemes improve FSTG as measured by bias, mean ab-solute error (mae), and normalized root-mean-squareerror (nrmse). TERR_DIFF is slightly better thanGAUSS, with smaller mae and nrmse. MD is betterthan GAUSS and TERR_DIFF, and MD_LSMG is thebest. Percentage improvement of rmse over FSTG isbelow 40% for GAUSS and TERR_DIFF but above50% for MD and MD_LSMG. Compared to MD,MD_LSMG gains 5% more improvement.To compare computational efficiency, each experi-ment is performed on one processor of a High-Performance Computing Linux Super-Cluster, with1-GHz Pentium III CPU. The results from a surface-only analysis at 0100 UTC 4 February 2003 are pre-sented in Table 5. For the purpose of computationalcomparison, all 169 available temperature observationsare used within a domain having 257 H11003 205 H11005 52 685GPs, with 3-km grid spacing. The execution time forGAUSS is slightly less than that for TERR_DIFF. Theexecution time for the MD approach is split into twoparts: 1) the calculation of the SFs (which needs to bedone only once for each surface station within an analy-sis domain, as mentioned in section 4a) and 2) the ap-plication of the SFs during the analysis (which is re-peated for each analysis). Part 1 of both MD andMD_LSMG is computationally more expensive com-pared to GAUSS and TERR_DIFF, while part 2 isnearly the same as TERR_DIFF. MD_LSMG is com-putationally cheaper than MD, because adding land–water contrasts limits the region of influence and hencespeeds convergence. The timing for part 1 is done with-out limiting the observation’s influence region. Thus,the computational costs are affected by the domain size.To eliminate the domain-size effect, the SFs are cal-culated within an influence region (a square here) ofthe observation. The radius of the circle inscribed inthe square was set to 300 km, which is slightly greaterthan the radius of influence (ROI; 273.14 km) calcu-lated in the ADAS based on Rhfor the first pass. Nowthe execution time for calculating the SFs is reduced36.5% for MD and 33% for MD_LSMG [see item (*) inTable 5].FIG. 10. Analysis increments (isopleths) for the potential tem-perature at the lowest (terrain-following) model level in responseto a single high-elevation observation at o3, indicated by a whitesolid triangle on the mountain top. Contour interval is 0.2 K.Darker shading indicates higher elevations (m). (top) GAUSS;(middle) TERR_DIFF; (bottom) MD.400 MONTHLY WEATHER REVIEW VOLUME 133FIG. 11. The northern half of the Georgia Basin area. The surface stations are shown bypluses, below which are station IDs. Station information is in Table 2. The bottom plot iszoomed to the top center of the top plot. Terrain elevations (m), indicated by shading, arefrom unsmoothed terrain data resolved by the NWP model at 3-km grid spacing. Darkershading represents higher elevations, with a maximum elevation difference of (top) 1374 mand (bottom) 1200 m.FEBRUARY 2005 DENG AND STULL 4017. Summary and discussionAn MD approach and its shoreline refinement areproposed for surface-data analysis on terrain-followingsurfaces in mountainous/coastal terrain. The resultingSFs and CTDs are then incorporated into the ADAS.For the cases examined here, the new scheme is com-pared with the original scheme (GAUSS) in the ADAS(ARPS5.0.0 Beta8), and the scheme (TERR_DIFF) de-veloped by Miller and Benjamin (1992). It is found thatthe MD approach can better account for both terrain-height and valley differences in the analysis of obser-vation increments over mountainous regions. The re-fined MD approach can further account for land–seacontrasts and is recommended for use in mountainous,coastal terrain such as southwestern BC.The one-time, initial computational cost of the MDapproach is greater than that of either GAUSS orTERR_DIFF because of the need to calculate the SFs.However, for any fixed set of the free parameters, theresulting SFs can be saved as a fixed file and can beapplied unchanged for each day’s analysis.The results from the real-observation case might befurther improved if the land/sea-breeze frontal locationis used to define the edges of the land mask, rather thanthe actual coastline. However, determining the positionof a sea/land-breeze front increases the complexity andrequires that the SFs be recomputed during everyanalysis. The real-data case shows encouraging resultsFIG. 12. Observed potential temperature (derived from ob-served temperature) for stations over water (lines with opencircles) and over land (lines with filled circles) during 0000–1200UTC 4 Feb 2003. Stations over water are Entrance Island (620,solid line), Pam Rocks (615, dashed line), and Sand Heads CS (671,dotted line). Stations over land are West Vancouver (673, solidline), Vancouver INTL ARPT (597, dashed line), and White Rock(674, dotted line). Breaks in the graph denote periods for whichno observations were reported. 0000 UTC corresponds to 1600PST (late afternoon), and 1200 UTC is 0400 PST (early morning).TABLE 3. Number of reporting stations for analysis andverification at each analysis time on 4 Feb 2003.No. of stationsTime (UTC) Analysis Verification0000 8 110100 8 110200 8 110300 8 110400 4 70500 7 110600 6 110700 7 110800 7 110900 7 111000 7 111100 7 111200 7 11TABLE 2. Surface stations for analysis and verification. Station IDs are from the Emergency Weather Net database. For theland–sea (LS) mask, 1 stands for land and 0 corresponds to ocean or water.Full name Station ID Lat (N) Lon (W) Elevation (m) Agency LS PurposesPort Mellon 38 49.5228° 123.4822° 3 WLAP 1MT Strachan 132 49.4167° 123.2000° 1420 MOTH 1Gold CK 417 49.4472° 122.4750° 794 HYDR 1Abbotsford ARPT 604 49.0333° 122.3667° 58 EC 1 AnalysisSaturna Island 621 48.7833° 123.0500° 24 EC 0Pitt Meadows 641 49.2000° 122.6833° 5EC 1Point Atkinson CS 654 49.3333° 123.2667° 35 EC 0Nanaimo ARPT 582 49.0500° 123.8667° 28 EC 1Crofton 59 48.8803° 123.6458° 20 WLAP 0Harmac Pacific 61 49.1353° 123.8475° 23 WLAP 0Deeks Peak 117 49.5333° 123.2167° 1280 MOTH 1MT Strachan Precip 138 49.4044° 123.1911° 1220 MOTH 1Vancouver INTL ARPT 597 49.1833° 123.1833° 2ECPam Rocks 615 49.5000° 123.3000° 10 EC 0 VerificationDiscovery Islands CS 617 49.4167° 123.2333° 15 EC 1Entrance Island 620 49.2167° 123.8000° 3EC 0Sand Heads CS 671 49.1000° 123.3000° 15 EC 0West Vancouver 673 49.3333° 123.1833° 178 EC 1White Rock 674 49.0167° 122.7667° 13 EC 1402 MONTHLY WEATHER REVIEW VOLUME 133for the MD approach, even with land–water contrastsfixed at the coastlines.To use the MD approach in mountainous/coastal ter-rain, fine-resolution NWP models (i.e., H9004x H11021 5 km) arerecommended, for two reasons. First, the SF is a func-tion of the elevations of analysis GPs; thus, there is aneed to resolve all important ridges and valleys. Sec-ond, the quality of the analyses in data-sparse ridges orvalleys depends strongly on the quality of the first guessfrom the NWP models.The MD technique is developed under the intraval-ley and intervalley decorrelation assumptions. None-theless, one can imagine scenarios where the observa-tion in one valley might happen to be well correlatedwith the air state in a neighboring valley. During a fair-weather event with strong solar heating, the BL couldbe deeper than the mountain-ridge height, allowing tur-bulence to mix air between the two valleys. Or duringstrong synoptic forcing, high winds and intense turbu-lence could inject similar air into both valleys and/orhomogenize the air in both. For these situations, run-ning the forecast model after initialization will mix theair in neighboring valleys due to resolved advection andparameterized turbulence, even if the initializationspreads the observations only within one valley.Sensitivity tests on the MD free parameters are per-formed for the real-data case (appendix). A subset ofavailable observations are analyzed using various val-ues of the free parameters, and the analyses are verifiedagainst the remaining observations. Nrmses are foundFIG. 13. Time series of observed potential temperatureand the [first guess (FSTG) or analyzed] potential tem-perature from four assimilation experiments for two sta-tions. Observations at Entrance Island over water: dashedline with open triangles. Observations at White Rock overland: dashed line with solid triangles. Analyses at EntranceIsland: solid line with open circles. Analyses at WhiteRock: solid line with filled circles. Analyses are betterwhen the solid lines are closer to the respective dashedlines.FEBRUARY 2005 DENG AND STULL 403to be far more sensitive to zref2 than to zref1, and thenrmse sensitivity to zref2 shows a daily cycle, corre-sponding to the characteristics of the BL. This impliesthat one could automate operationally detecting decor-related weather versus well-mixed events by using ob-served BL depth (if available) as zref2 during dailyanalyses. Further experiments will be required to inves-tigate the sensitivities for other times and locales. Theoptimum values of the free parameters might be depen-dent on weather, season, and topography. The resultsshown here suggest that the free parameters can bedetermined via sensitivity tests over a substantial pe-riod of daily analyses.For GAUSS and TERR_DIFF, a shorter Rhcould beused to successfully reduce the intervalley spreadingFIG. 14. Time series of averaged potential temperatureof three stations over water and land, respectively. Thethree stations over water, used for the average, are En-trance Island, Pam Rocks, and Sand Heads CS. West Van-couver, Vancouver INTL ARPT, and White Rock are thethree stations used for the average over land. Observationsand analyses over water are denoted by the dashed linewith open triangles and the solid line with open circles,respectively. Observations and analyses over land are rep-resented by dashed line with solid triangles, and solid linewith filled circles, respectively. FSTG stands for the firstguess. Closer correspondence between the solid anddashed lines indicate a better analysis.TABLE 4. Verification of analyzed potential temperatures interms of bias, mean absolute error (mae), and normalized root-mean-square error (nrmse) for all reporting verification stationsduring the analysis period from 0000 to 1200 UTC 4 Feb 2003. Nequals number of stations times the number of observation times.Nrmse is root-mean-square error (rmse) for each experiment nor-malized by the rmse of the first guess (FSTG). Percent improve-ment of the rmse is with respect to FSTG.Experiment NBias(K)Mae(K)Nrmse(K)PercentimprovementFSTG 139 H110022.6832 2.7624 1.0 0GAUSS 139 H110020.5625 1.7565 0.6610 34TERR_DIFF 139 H110020.6401 1.6482 0.6317 37MD 139 0.2344 1.3069 0.4943 51MD_LSMG 139 H110020.0051 1.1956 0.4433 56404 MONTHLY WEATHER REVIEW VOLUME 133from a disconnected neighboring valley (appendix).However, this strategy also reduces the potentiallygood intravalley spreading. In addition, a universal Rhmight be hard to find in highly variable terrain. TheMD approach suppresses intervalley spreading whileallowing intravalley spreading over a longer distance.This approach automatically accounts for the elevationsof surface stations and analysis GPs.The ultimate goal is to provide better input fields forfine-resolution NWP models. To produce improvedanalyses to drive the MC2, we plan to combine themesoscale surface analyses from the MD approach withthe coarse-resolution Eta analyses from NCEP, and tocompare verification statistics with MC2 control runsdriven by unmodified Eta initial conditions. We furtherplan to modify MD to spread free-atmosphere moun-tain-top observations differently than BL valley obser-vations.Acknowledgments. We thank Joshua Hacker fromNCAR for his suggestions, and Stephane Chamberlandof RPN for MC2 user support. We also thank HenrykModzelewski, George Hicks, Yongmei Zhou, JohnSpagnol, and other colleagues at UBC. Thanks also tothree reviewers, including Steven Lazarus, John Horeland Keith Brewster. Grant support came from the Ca-nadian Natural Science and Engineering ResearchCouncil, the BC Forest Investment Account, Environ-ment Canada (EC), and the Canadian Foundation forClimate and Atmospheric Science. Observation datafrom the Emergency Weather Net were provided byEC, BC Ministry of Transportation and Highways, BCMinistry of Water Land and Air Protection, and BCHydro. Geophysical Disaster CFD Center computerswere used, funded by the Canadian Foundation for In-novation, the BC Knowledge Development Fund, andUBC.APPENDIXSensitivity TestsThis appendix provides a sensitivity analysis ofMD_LSMG to different values of the free parametersfor the case from section 6b. The analysis and verifica-tion observations remain the same. This single case isinsufficient to determine the optimum values. How-ever, it does indicate some sensitivity characteristics.First, we vary zref1 and zref2 from 250 to 3000 mevery 250 m, while holding a H11005 2 and b H11005 2. Recall thatzref1 and zref2 are the maximum effective BL heightsthrough which surface observations are felt.Figure A1 shows nrmses versus zref1 and zref2 forTABLE 5. Execution time (s) for each analysis experiment onone processor of a High-Performance Computing Linux Super-Cluster, with 1-GHz Pentium III CPU. The results are obtainedfrom a surface-only analysis at 0100 UTC 4 Feb 2003, when all 169observations are used for analysis within a domain of 257 H11003 205grid points, with 3-km grid spacing. The parameters used in themother-daughter approach are a H11005 2, b H11005 2, zref1 H11005 750 m andzref2 H11005 750 m. For the two MD approaches, item 1 is the one-timesetup cost without limiting search radius, item 2 is the executiontime for its application to this test domain, and item * is also theone-time setup cost but with limiting search radius to 300 km.Experiment Execution time (s)GAUSS 22.99TERR_DIFF 28.12MD (1) 3902.03(2) 29.96(*) 2476.80MD_LSMG (1) 2833.89(2) 29.93(*) 1897.60FIG. A1. Nrmse of analyses vs zref1 and zref2 when a H11005 2 andb H11005 2. Contour interval is 0.025 K. (top) Late afternoon (0000UTC 4 Feb 2003); (bottom) early morning (1200 UTC 4 Feb2003). Smaller nrmse is better.FEBRUARY 2005 DENG AND STULL 4050000 UTC (1600 PST) and 1200 UTC (0400 PST) 4February 2003. At both 0000 and 1200 UTC, nrmses areless than 1.0 for any values of zref1 and zref2 within therange 250–3000 m. Thus, the analyses agree with theobservations better than the first guess. Nrmse is notsensitive to zref1 when zref1 H11022 750 m, or zref2 H11021 500 m.This is because the elevation differences between amother GP and its immediate daughter GPs are usuallysmall. Nrmse is far more sensitive to zref2 than to zref1,as zref2 defines the maximum effective height abovethe observation, and because the elevation differencesbetween the observation and the analysis GPs vary con-siderably.Figure A1 shows different features at 0000 and 1200UTC. In late afternoon (0000 UTC), nrmse shows smallsensitivity to zref1 H11022 750 m and zref2 H11022 750 m. In earlymorning (1200 UTC), lower nrmse is confined to a nar-row band from zref1 H11005 500 to 3000 m and from zref2 H110051000 to 1250 m. Similar sensitivity tests from 0100 to1100 UTC indicate a transition at 0200 UTC (notshown). The behaviors correspond to the characteris-tics of the evolving BL.Second, we vary a and b from 0.5 to 4.0 every 0.5 for1200 UTC 4 February 2003. For these experiments,zref1 H11005 2250 m and zref2 H11005 1000 m, which were the bestvalues from the previous tests. As shown in Fig. A2, aH11005 2 and b H11005 2 (or a H11005 1 and b H11005 2.5) yield lowernrmses.Using the CTD in the MD approach effectively limitsthe region of influence in complex terrain (section 5).GAUSS and TERR_DIFF ROI sensitivity tests areperformed by varying horizontal correlation lengthscale Rh(used to calculate ROI in the ADAS) from 100to 10 km every 10 km for the first and second Bratsethpasses. For each test, the third Bratseth pass reduces Rhby 70%, while the fourth and fifth passes reduce Rhby70% from the third pass.For the virtual-observation case (section 6a), rmsedecreases rapidly with decreasing Rhfor both GAUSSand TERR_DIFF (Fig. A3), as the correction from o1becomes negligible in the Lillooet River Valley. WhenRhH11005 20 km, a minimal rmse is found to be 0.3385 K forGAUSS and 0.4105 K for TERR_DIFF. The rmse withRhH11005 20 km is reduced about 37% for GAUSS and 24%for TERR_DIFF from the first guess. While such smallvalues of Rhsuccessfully reduce intervalley spreading,they do so at the expense of reducing the potentiallygood intravalley spreading. The MD approach avoidsthis problem.FIG. A2. Nrmse of analyses (1200 UTC 4 Feb 2003) vs a and bwhen zref1 H11005 2250 m and zref2 H11005 1000 m, which were the bestvalues for 1200 UTC from Fig. A1. Contour interval is 0.025 K.Smaller nrmse is better.FIG. A3. Rmse of analyses (0000 UTC 8 Mar 2003) vs Rhusedin the first and second Bratseth passes. The third Bratseth pass hasan Rhthat is reduced by 70% from the second pass, while thefourth and fifth Bratseth passes have a Rhthat is reduced by 70%from the third pass. (top) GAUSS; (bottom) TERR_DIFF.406 MONTHLY WEATHER REVIEW VOLUME 133Similar tests are performed for the real-observationcase (section 6b). Figure A4 shows nrmses versusRhused in the first and second Bratseth passes for1200 UTC 4 February 2003. For both GAUSS andTERR_DIFF, nrmses are less than 1.0. Nrmse showssmall sensitivity to Rhwhen RhH11022 40 km but relativelylarge sensitivity when RhH11021 40 km. A minimal nrmse isfound at RhH11005 20 km for GAUSS and at RhH11005 10 km forTERR_DIFF. Identical tests are performed for 0000UTC (Fig. A5). Different from the 1200 UTC results,RhH11005 100 or 90 km produces lower nrmses for bothGAUSS and TERR_DIFF.The ROI sensitivity for the real-observation case ismuch smaller than that for the virtual-observation case.This is not surprising, considering that only two obser-vations in the two adjacent valleys were used for thevirtual-observation case. These results indicate that aprescribed Rhdoes not work everywhere in highly vari-able terrain. An Rhthat depends on the elevations ofsurface stations and analysis GPs might help to improvethe analyses from GAUSS and TERR_DIFF. In com-plex terrain, an effective Rhmight not only depend onthe average station separation, but also on the geo-graphic regions, seasons, and weather.REFERENCESArnold, C. P., and C. H. Dey, 1986: Observing systems simulationexperiments: Past, present, and future. Bull. Amer. Meteor.Soc., 67, 687–695.Barnes, S. L., 1964: A technique for maximizing details in numeri-cal weather map analysis. J. Appl. Meteor., 3, 396–409.FIG. A4. Nrmse of analyses (1200 UTC 4 Feb 2003) vs Rhusedin the first and second Bratseth passes. The third Bratseth pass hasan Rhthat is reduced by 70% from the second pass, while thefourth and fifth Bratseth passes have a Rhthat is reduced by 70%from the third pass. (top) GAUSS; (bottom) TERR_DIFF.Smaller nrmse is better.FIG. A5. Same as Fig. A4, but for 0000 UTC 4 Feb 2003. (top)GAUSS; (bottom) TERR_DIFF. Smaller nrmse is better.FEBRUARY 2005 DENG AND STULL 407Benjamin, S. G., K. Brewster, R. Brummer, B. Jewett, T. W.Schlatter, T. L. Smith, and P. A. Stamus, 1991: An isentropicthree-hourly data assimilation system using ACARS aircraftobservations. Mon. Wea. Rev., 119, 888–906.Benoit, R., M. Desgagne, P. Pellerin, S. Pellerin, Y. Chartier, andS. Desjardins, 1997: The Canadian MC2: A semi-Lagrangian,semi-implicit wide band atmospheric model suited for fine-scale process studies and simulation. Mon. Wea. Rev., 125,2382–2415.Bouttier, F., 1993: The dynamics of error covariances in a baro-tropic model. Tellus, 45A, 408–423.Bratseth, A. M., 1986: Statistical interpolation by means of suc-cessive corrections. Tellus, 38A, 439–447.Brewster, K., 1996: Application of a Bratseth analysis schemeincluding Doppler radar data. Preprints, 15th Conf. onWeather Analysis and Forecasting, Norfolk, VA, Amer. Me-teor. Soc., 92–95.CAPS, 1995: ARPS version 4.0 user’s guide, Supplement 1: 3-Danalysis with ARPS data assimilation system: ADAS version2.3. Center for Analysis and Prediction of Storms, Universityof Oklahoma. [Available online at]Case, J. L., J. Manobianco, T. D. Oram, T. Garner, P. F. Blottman,and S. M. Spratt, 2002: Local data integration over east-central Florida using the ARPS Data Analysis System. Wea.Forecasting, 17, 3–26.Ciliberti, C. M., J. D. Horel, and S. M. Lazarus, 2000: Sensitivityexperiments with a high resolution data assimilation scheme.Preprints, Ninth Conf. on Mountain Meteorology, Aspen,CO, Amer. Meteor. Soc., 413–416.——, J. D. Horel, and S. M. Lazarus, 1999: An analysis of a coldfrontal passage over complex terrain in northwest Utah. Pre-prints, Eighth Conf. on Mesoscale Processes, Boulder, CO,Amer. Meteor. Soc., 459–462.Daley, R., 1991: Atmospheric Data Analysis. Cambridge Univer-sity Press, 457 pp.Desroziers, G., 1997: A coordinate change for data assimilation inspherical geometry of frontal structures. Mon. Wea. Rev., 125,3030–3038.De Wekker, S. F. J., M. Kossmann, and F. Fiedler, 1997: Obser-vations of daytime mixed layer depths over mountainous ter-rain during the TRACT field campaign. Preprints, 12thSymp. on Boundary Layers and Turbulence, Vancouver, BC,Canada, Amer. Meteor. Soc., 498–499.Gauthier, P., P. Courtier, and P. Moll, 1993: Assimilation of simu-lated wind lidar data with a Kalman filter. Mon. Wea. Rev.,121, 1803–1820.Hessler, G., 1984: Experiments with statistical objective analysistechniques for representing a coastal surface temperaturefield. Bound.-Layer Meteor., 28, 375–389.Horel, J., T. Potter, L. Dunn, W. J. Steenburgh, M. Eubank, M.Splitt, and D. J. Onton, 2002: Weather support for the 2002Winter Olympic and Paralympic Games. Bull. Amer. Meteor.Soc., 83, 227–240.Houtekamer, P. L., and H. L. Mitchell, 2001: A sequential en-semble Kalman filter for atmospheric data assimilation. Mon.Wea. Rev., 129, 123–137.——, and ——, 1998: Data assimilation using an ensemble Kal-man filter technique. Mon. Wea. Rev., 126, 796–811.Kalthoff, N., H.-J. Binder, M. Kossmann, R. Vogtlin, U. Cors-meier, F. Fiedler, and H. Schlager, 1998: Temporal evolutionand spatial variation of the boundary layer over complexterrain. Atmos. Environ., 32, 1179–1194.Kossmann, M., R. Vogtlin, U. Corsmeier, B. Vogel, F. Fiedler,H.-J. Binder, N. Kalthoff, and F. Beyrich, 1998: Aspects ofthe convective boundary layer structure over complex ter-rain. Atmos. Environ., 32, 1323–1348.Krishnamurti, T. N., and L. Bounoua, 1996: An Introduction toNumerical Weather Prediction Techniques. CRC Press, 293pp.Lanzinger, A., and R. Steinacker, 1990: A fine mesh analysisscheme designed for mountainous terrain. Meteor. Atmos.Phys., 43, 213–219.Laprise, R., D. Caya, G. Bergeron, and M. Giguère, 1997: Theformulation of the André Robert MC2 (mesoscale compress-ible community) model. Numerical Methods in Atmosphericand Oceanic Modelling, The André J. Robert Memorial Vol-ume, Atmosphere—Ocean, 195–220.Laroche, S., P. Gauthier, J. St-James, and J. Morneau, 1999:Implementation of a 3D variational data assimilation systemat the Canadian Meteorological Center. Part II: The regionalanalysis. Atmos.–Ocean, 37, 281–307.Lazarus, S. M., C. M. Ciliberti, J. D. Horel, and K. Brewster, 2002:Near-real-time applications of a mesoscale analysis system tocomplex terrain. Wea. Forecasting, 17, 971–1000.Lenschow, D. H., B. B. Stankov, and L. Mahrt, 1979: The rapidmorning boundary-layer transition. J. Atmos. Sci., 36, 2108–2124.Miller, P. A., and S. G. Benjamin, 1992: A system for the hourlyassimilation of surface observations in mountainous and flatterrain. Mon. Wea. Rev., 120, 2342–2359.Stull, R. B., 1988: An Introduction to Boundary Layer Meteorol-ogy. Kluwer Academic, 666 pp.——, 1992: A theory for mixed-layer-top levelness over irregulartopography. Preprints, 10th Symp. on Turbulence and Diffu-sion, Portland, OR, Amer. Meteor. Soc., 92–94.——, 1998: Studies of mountain climates, and a theory for mixed-layer-top levelness over complex topography. GeophysicalDisaster Computational Fluid Dynamics Center Tech. Rep.WFRT-98-001, University of British Columbia, Vancouver,Canada.——, 2000: Meteorology for Scientists and Engineers. 2d ed.Brooks/Cole, 502 pp.Thépaut, J.-N., P. Courtier, G. Belaud, and G. Lemaitre, 1996:Dynamical structure functions in a four-dimensional varia-tional assimilation: A case study. Quart. J. Roy. Meteor. Soc.,122, 535–561.Xue, M., K. K. Droegemeier, and V. Wong, 2000: The advancedregional prediction system (ARPS)—A multi-scale nonhy-drostatic atmospheric simulation and prediction model. PartI: Model dynamics and verification. Meteor. Atmos. Phys., 75,161–193.408 MONTHLY WEATHER REVIEW VOLUME 133


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