This paper was originally written in TeX. Since the advent of Netscape "plug-in" features that enable browsing TeX documents effortlessly, I provide the original TeX document in addition to the html version. For this html version, I have not spent a great deal of time cleaning up a lot of the "texizms," but equations are now inserted as small .gif images. The entire appendix is also included as a .gif image since it is difficult to interpret in TeX script.
Joint collaboration with the Cooperative Institute for Research in the Atmosphere, Fort Collins, CO
Boulder, CO
Horizontal Shape Matching (HSM) is a variational technique to obtain an analysis that merges the gradient structure from one data source (e.g., satellite) with a background field. This technique is modeled after the variational methods that fall into two general categories: strong constraint and weak constraint. HSM, which is essentially a specialized filter function, falls in the weak constraint category. Here a technique is demonstrated that effectively tunes the filter to the spatial resolution inherent in the data sources rather than on the absolute error of the data.
The HSM equation is developed for satellite gradient insertion including a cloud treatment. The filtering properties of the HSM equation are explored using Fourier analysis and an approximate numerical solution is shown to be viable for operational use. Applying this technique to a data set spanning more than a year proves the analysis has a positive effect through the statistical significance of a large sample of comparisons with RAOB data.
Horizontal Shape Matching (HSM) is an analysis technique used to incorporate structure into an analysis field based on variational methods (Sasaki 1958 and 1970, McGinley 1984). Such variational techniques have been used in purely diagnostic (Achtemeier et al. 1986) and prognostic applications (Chance 1986). At the Forecast Systems Laboratory (FSL) it is used diagnostically to insert satellite structure from Geostationary Operational Environmental Satellite (GOES)-7 data into a 10~km local water vapor analysis as part of the Local Analysis and Prediction System (LAPS). In this application, HSM solves a boundary value problem without time variation terms.
LAPS, a relatively new analysis system under development at FSL, addresses an anticipated need for analysis and modeling on the local scale that will arise in the Advanced Weather Interactive Processing System (AWIPS) era. LAPS uses both locally and nationally disseminated data, and its output can initialize local-scale models that address specific forecast problems. LAPS currently operates on a 10-km mesh approximately centered over Denver, Colorado. The square grid is 600~km on each side and the vertical coordinate is pressure with 50~hPa increments spanning surface to 100~hPa. LAPS analyzes all state variables in separate three-dimensional and surface domains along with special fields such as aircraft icing threat. Additional documentation on the LAPS system can be found in McGinley et al. (1991) and Birkenheuer (1992a).
FSL has operated LAPS with the HSM satellite gradient insertion technique for two years. It is the only part of LAPS that uses GOES Visible and Infrared Spin-Scan Radiometer (VISSR) Atmospheric Sounder (VAS) data. This paper develops the HSM technique, reviews the input data, and describes how the method applied to GOES-7 data provides measurable improvement to the water vapor analysis.
a. General Background
The advantage of adding only gradient information from a source such as satellite data is to avoid adding data bias. Although GOES satellite data have questionable absolute accuracy, we know that the features seen in moisture fields are realistic (Birkenheuer 1991). HSM uses that structural information without adversely perturbing field bias.
The fields used in this study to supply the gradient structure for our local scale moisture analysis are described in Birkenheuer 1991, with a slight modification (Birkenheuer 1992b). A physical retrieval algorithm adapted from the University of Wisconsin -- Madison (Hayden 1988) generates temperature and dew point temperature profiles using FSL's national-scale Mesoscale Analysis and Prediction System (MAPS) analyses or forecast fields (Benjamin 1991) as a first guess along with GOES-VAS infrared (IR) radiance data. The sounding retrieval domain covers a large part of the west-central United States. The retrieval algorithm produces a temperature and moisture profile every 60~km in clear regions. Integrating the moisture soundings over three independent layers spanning the vertical column provides vertical moisture structure. The set of integrated moisture profiles and their corresponding radiances form a sample from which the coefficients are derived through regression. A solution to the following equation is desired.
where A is an nx12 matrix containing n observations of brightness temperatures (12 channels) from the satellite data for a given period and s is an n-element vector containing the integrated precipitable water for a specified layer corresponding to each physical sounding retrieval (training set). These exist at each clear MAPS grid point. The coefficient vector c for a particular layer is solved by
This is solved separately for each layer using iterative refinement. Coefficients are computed uniquely for each period and applied only to data of that same time.
Following the solution of (2), the coefficients are applied to the 12 brightness temperatures for each field of view (FOV) over the LAPS domain. Four fields are computed in this manner. One for the case when c represents the total precipitable water (TPW) spanning the surface to 300~hPa. Three other fields are computed for separate atmospheric layers, layer precipitable water (LPW), by computing unique c vectors for each layer. In the LAPS water vapor algorithm, the TPW field is scaled to agree with ground-based, dual-channel passive microwave radiometers. Real-time total column water vapor measurements are available at two locations in the LAPS domain from ground-based passive microwave radiometers operated by NOAA's Environmental Technology Laboratory (Guiraud et al. 1979). These data help scale the water vapor fields by fitting analyzed integrated water to this real-time data source. After this step TPW is only used to assure consistency in the LPW field by forcing the integrated total to agree with the TPW field on a point-by-point basis. Cloudy areas are flagged and do not influence the final analysis (discussed in detail in section 3).
By using a relatively large geographic area for the coefficient determination, the problem of cloud obscuration is overcome because the dataset is large enough to provide ample clear scenes. Furthermore, the observations also represent a multistate area that extends over widely varying terrain. This ostensibly is a dataset that provides a valid representation of the air mass over the LAPS domain and is sensitive to topographical variability.
b. Selection of Layers
Layer selection results from two considerations: minimal error and meteorological application. The measure of error is the square root of the sum of squared differences (L2 norm) between the computed LPW using the regression coefficients and its ``raw data" counterpart from the retrieval package. The meteorological consideration attempts to partition the atmosphere into layers with some physical meaning. Three layers were subjectively selected after examining several plots, each containing monthly data that produced consistent error minima. The atmosphere is partitioned at 780 and 640~hPa defining three layers denoted as i where: i=1 corresponds to surface to 780 hPa (boundary layer over the high plains), i=2 from 780 to 640~hPa (low troposphere in the free atmosphere), and i=3 above 640~hPa. The surface pressure was arbitrarily assigned to be 840~hPa. Details on this partitioning can be found in Birkenheuer 1992a,b. Distributing the horizontal structure in each of these layers to the 21 levels comprising the LAPS vertical coordinate is discussed in section 7.
The general form of the HSM is derived by minimizing a functional of the form:
where $\phi$ is the desired function, $\phi^g$ is the field containing gradient information, and $\phi'$ is a background field. The terms $\alpha$, $\beta$, and $\gamma$ weight the response of the Laplacian (smoothing), gradient structure, and background respectively. The x and y subscripts denote partial derivatives. Note that the $\beta$ term can be a function of $x$ and $y$; because in our satellite application the $\beta$ term depends on cloud cover. Values for $\beta$ range from a maximum down to 0 in cloudy areas. One can envision this in terms of a constant $\beta$ modulated by a spatially varying clear fraction ($\psi$); ranging from 1.0 (clear) to 0.0 (cloudy). This approach assures that the analysis excludes gradient information in cloudy regions that obscure clear-air satellite radiances. In cloudy regions, the background field dominates the analyzed product. Furthermore, $\psi$ is directly available from the LAPS cloud analysis (Albers 1994).
The Euler-Lagrange equation (derivation shown in the Appendix) that results after applying variational calculus to (3), assuming constant boundary condidions, and including the clear function has the form
It should be understood that we choose to fix the boundary conditions for this solution. As a result, to be effective, the domian should be selected large enough so that the boundary approximation does not pose a problem.
One of the primary problems in applying HSM is the selection of the weighting coefficients ($\alpha$, $\beta$, and $\gamma$). The coefficient magnitudes should be inversely proportional to error. Lacking a good estimate of ``gradient error" and how it relates to absolute error in measurements or the background field can pose a problem in coefficient selection. Alternatives to absolute ``gradient error" are the relationships between gradients, spatial resolution of the analysis grid, and instrument resolution. The HSM equation is a very effective filter that blends the background field with the gradient information. By selecting coefficients appropriately, we can limit the background field's influence to its parent grid resolution while restricting the satellite gradient's influence to its intrinsic spatial resolution. The next step is to examine the HSM equation in terms of its filtering capability.
Our goal in applying HSM is to insert the gradient structure in a balanced fashion that offers the best blend of information from both the gradient and background fields. We simplify the problem by examining the clear case. Here $\psi$ is unity and its derivatives are zero. Equation 4 reduces to
The gradient response of (5) can be isolated by zeroing out the background field ($\phi' = 0$) and analyzing the response function. The filter response of this equation can be studied by applying it to waves in two dimensions. Let
where k and p are the wave numbers in the x and y dimensions, respectively. Substituting (6) and (7) into (3) and evaluating gives
after dividing through by $e^{i(kx+py)}$. The expression can be reduced to one-dimensional form by defining wavenumber $n$ such that $n^2 = k^2 + p^2$. Wavenumber (n) is related to wavelength by $n= {2\ \pi \over \lambda}$. The gradient response function, r, is the ratio of the analyzed wave amplitude A over the input amplitude B and equals
At this point we also recognize that the coefficients work in a relative way with one another. By assigning $\gamma=1$ we restrict coefficient determination to $\alpha$ and $\beta$.
To determine the peak response of the gradient term, we solve the max/min problem
and we specify that $\alpha > 0$ and $\beta > 0$, from which we conclude that the peak wavenumber is dependent only on $\alpha$.
We see the $\alpha$ term selects the peak frequency of gradient response and the $\beta$ term controls the magnitude of the gradient influence on the result.
However, the gradient response is not the only item to consider when selecting $\alpha$ and $\beta$. The background response is also affected by this choice. Specifically, (5) can be reassessed for an analytical background response in the same way the gradient response was derived producing,
The background high frequency ``cutoff" lowers in response to increasing $\beta$, given that $\alpha$ is constrained to enhance the peak gradient response (11). The $\alpha$ and $\beta$ term must be chosen so as not to obliterate any background effects while maintaining a desired level of gradient influence. However there is a caveat: If low-frequency gradient data are merged with a high-frequency background field, this technique probably will not work.
Coefficient selection takes into account both the gradient and background effects for the desired response. A desired response is one where the background has influence down to its information scale (estimated to be near 180~km wavelength for MAPS data), while allowing the gradient source to influence down to approximately 60~km, which is near the Nyquist wavelength for the spatially averaged water vapor fields derived from satellite data (Birkenheuer 1991). By selecting the gradient peak at 100~km a value of $\alpha$=64162~km$^{-4}$ is derived from (11). If we specify a background cutoff of 180~km we can compute a value of $\beta$=742~km$^{-2}$ using (12). With this selection of $\alpha$ and $\beta$ we observe a peak gradient response of 0.594 from (9).
plots the gradient response (9) using this set of coefficients. We see it has properties of a band-pass filter with zero response at $\lambda$=0 and asymptotically zero at $\lambda$=$\infty$. The slope of the response is steep in the 60~km wavelength region, dropping off rapidly for short wavelengths. The gradient response is greater than 0.5 between 175~km and 62.8~km. Hence, wavelengths shorter than desired will be attenuated while supplying favorable gradient and background information to the final result.
There are alternate ways of arriving at reasonable values for these coefficients. Again, one may start by assigning a peak gradient wavelength to obtain $\alpha$ using (11). Then using (9) directly one can plot $r$ as a function of $\beta$ observing both the limits of the response (region above the 0.5 response level) and the magnitude of the gradient response to settle on a reasonable value for $\beta$. However, one must keep in mind that the selected coefficients have an influence on the background field response.
a. The true solution
Fourier analysis can be used to study the influence of coefficients on response and the performance of the partial differential equation (PDE) integration. This is accomplished by using test fields with known wavelengths and amplitudes that represent both background and gradient fields. Input fields to (4) are two orthogonal standing waves with equal amplitude
(Figure 2a and 2b). Waves in one dimension (y) represent the background and those in the other dimension (x) depict a gradient signal. Cases for many combinations of the background and gradient wavenumbers can be tested. Furthermore, it directly indicates the quality of the numerical solution to the PDE.
Figure 3a shows the Fourier analysis of the result after combining the specific standing waves for the background and gradient in 2a and 2b through HSM (4) using the coefficients derived in section 5. The relative magnitudes (power spectrum) of the peaks (two for each wave) are indicative of the response of the HSM filter function. It is through Fourier analysis that we verify the numerical solution, since the computed response matches the analytical one. Slight differences are attributed to insufficient iterations in the PDE integration scheme and some minor dispersion of power (shown in Figure 3a as minor ``bumps" beneath each peak). The fact that we now have a tool to observe exactly what the numerical scheme is doing offers significant benefits. We can now check the performance of the numerical integrator, and also tune it. Approximately 2,000 iterations are needed to get the HSM equation to converge near the analytic solution. This computation is both time consuming and expensive. In practice, we need a quick solution because of meteorological data perishability; furthermore, we want to do this on inexpensive hardware.
b. An approximate solution
It was discovered through trial and error that by substantially increasing the value of the weighting coefficients in (4), a response very close to the true solution could be computed at a fraction of the iterations. Understand that this numerical calcuation if allowed to go beyond a specified number of iterations would converge to the ``wrong answer." These ad-hoc coefficients were developed by restricting the number of iterations and adjusting the $\alpha$ and $\beta$ values to best fit the response of the ``true" solution (Figure 1) through Fourier analysis. We will refer to this as the approximate solution.
The coefficients for the approximate solution were $\alpha$ = 3.943 x 10$^6$~km$^4$ and $\beta$ = 3.695 x 10$^4$~km$^2$. Table 1 shows the results when varying the number of iterations with the approximate solution.
The table shows that as the number of iterations increases from 50, the background cutoff changes to longer wavelengths; furthermore, there is a gradual increase in the gradient response, which is just slightly over 0.5 at 150 iterations. We also see that the background cutoff wavelength at 150 iterations is perhaps not as low as desired, but nevertheless offers a solution near the desired background cutoff of 200~km. The number of iterations was routinely fixed at 150 for the approximate solution.
Figure 3b plots the computed numerical gradient response (computed through Fourier power spectra) using the approximate method and a fixed background wavelength of 630~km and a series of gradient waves with wavenumbers ranging from 9.97x10$^{-3}$ up to 0.299~km$^{-1}$ (corresponding to wavelengths of 630 down to 21~km). Also plotted is the gradient response obtained from the ``true solution." The slightly lower response was deemed an acceptable trade-off for the more efficient numerical computation. It may be possible to improve the result of the approximate scheme further by adjusting the coefficients or the number of iterations. Of course, computing the true solution with $\alpha$ = 64162~km$^{-4}$ and $\beta$ = 742~km$^{-2}$ and 2000 or more iterations would always be preferred if possible.
Another consideration appears in Figure~3c, which plots the gradient response (derived from two-dimensional Fourier power spectra) of all combinations of background and gradient standing wavenumbers again using the approximate solution. It reveals a fairly uniform and satisfactory response over the whole range of standing wavenumbers.
To summarize, we can look to the analytic response functions along with given limitations in spatial resolution to determine suitable HSM coefficients. HSM can be applied by solving the equations exactly with those coefficients or by applying approximate methods that demonstrate a similar response. The approximate solution was used routinely and produced the results in this paper.
To distribute the structure from the three LPW fields into the 21 LAPS pressure levels, we first derive a representative specific humidity (SH) field from each LPW value, and then interpolate the structure in the vertical.
To merge LPW structure with SH data, LPW must be scaled to SH units. LPW in LAPS is expressed in units of length; however, since liquid water has a density of 1g cm-3, LPW also can be expressed in units of mass per unit area. A hypothetical SH field exists at the vertical midpoint of each LPW layer that can be computed from
where $LPW$ is precipitable water in units of g cm-2, $\Delta p$ is the pressure difference in the vertical over the layer of interest (hPa), g is the acceleration of gravity (cm s-2), and SH(g kg-1) is the layer's average specific humidity. We presume the structure in the LPW-generated SH layers can be linearly interpolated in the vertical (pressure space is linear with respect to mass) to portray the gradients at an analyzed SH level.
The $\beta$ term in the HSM equation can be divided into three parts, one for each SH level: $\beta_1$, $\beta_2$, and $\beta_3$. Their sum is constrained to equal $\beta$ as follows. For a specific LAPS level, we describe $\beta$ = $\sum \beta_i$ where $\beta_i \ =\ f_i\beta$ and
$\beta_i$ effectively distributes the structure from the three SH layers ($\phi^g_i$) at the level of interest. Table 2 lists the various fractional weights (fi) applied to the precipitable water layer as it relates to each LAPS pressure level. Though this treatment permits all three precipitable water layers to influence any LAPS pressure level, at this time only the two adjacent precipitable water layers apply to each LAPS pressure level; the fi values used produce a linear interpolation of structure in the vertical.
| LAPS pressure level (hPa) | fi=1 | fi=2 | fi=3 |
| 900 | 1.00 | 0.00 | 0.00 |
| 850 | 1.00 | 0.00 | 0.00 |
| 800 | 1.00 | 0.00 | 0.00 |
| 750 | 0.42 | 0.58 | 0.00 |
| 700 | 0.00 | 1.00 | 0.00 |
| 650 | 0.00 | 0.82 | 0.18 |
| 600 | 0.00 | 0.68 | 0.32 |
| 550 | 0.00 | 0.52 | 0.48 |
| 500 | 0.00 | 0.38 | 0.62 |
| 450 | 0.00 | 0.23 | 0.77 |
| 400 | 0.00 | 0.09 | 0.91 |
| 350 | 0.00 | 0.00 | 1.00 |
| 300 | 0.00 | 0.00 | 1.00 |
After partitioning the weighting coefficients and layer gradient fields, the ``operational" HSM equation valid at a particular LAPS level becomes
which can be solved numerically using over-relaxation techniques. If the $\psi$ term is zero, as obtained from the LAPS cloud analysis, no gradient structure will be added to the analysis because nothing would be gained by this step. Therefore, one of the first steps in the analysis was to integrate $\psi$ to evaluate the percentage of the scene that is clear. If less than 33\% of the scene is clear, it was discovered through trial and error that performing HSM could do more damage than good. All of the analyzed results contained in the paper ascribe to the criteria of $>$ 33\% clear.
Figure 4 shows the source of the gradient information derived from GOES-7 VAS for 1715 UTC 12 October 1994. Figure 4a shows the total column precipitable water and Figures 4b, 4c, and 4d show the precipitable water for the three levels described earlier. The three LPW fields share the same horizontal structure with an axis of higher moisture east of the Rocky Mountains extending from northeastern New Mexico across eastern Colorado into Nebraska.
Figure 5 shows stages of the LAPS specific humidity analysis, basically the field before and after HSM analysis for 1800 UTC 12 October. Figure 5a shows the 600~hPa background field for HSM that reflects inserted moisture where LAPS has analyzed clouds (Albers et al. 1994, 1995), possible boundary layer mixing and scaling to local radiometer measurements (Birkenheuer 1993). We note the presence of scattered clouds over the Rocky Mountains and south of the Nebraska panhandle. Cloud effects are easily detected since they appear as small-scale, tightly ``packed" gradients in the otherwise smooth moisture contours. In addition to the higher levels of moisture analyzed in the cloudy areas we see general high moisture regions in southern Wyoming, western Nebraska and northeast New Mexico.
Figure 5b shows the result of the HSM analysis in which $\alpha$ = 64162 km-4 and $\beta$ = 742~km-2 with 2000 iterations. Here we see that the Laplacian term in the HSM solution has smoothed the cloudy areas. In the southwestern portion of the domain we see that the cloudy area on the domain edge has not changed because of the fixed boundary conditions. We see that the level of moisture in southern Wyoming is near the background values, and finally note the very pronounced insertion of the gradient structure along the eastern Colorado high plains that is seen in the gradient fields. As Table 2 indicates, the structure at 600~hPa originates from the top two LPW fields (Figures 4c and 4d). We accept this as the ``optimal" result since it was solved with sufficient iterations and the coefficients ``tuned" for the spatial characteristics of the input data.
Figure 5c shows the HSM solution using the same background and gradient data as the 5b solution with the difference being $\alpha$ = 3.943 x 106~km4 and $\beta$ = 3.695 x 104 ~km2 and the total number of iterations is 150. Obviously this more efficient algorithm is still an approximation, and when we compare it with Figure 5b, we see the differences. Basically, values in the peaks and valleys differ up to 0.2~g/kg, with the predominant one being the analyzed high in southern Nebraska (1.7 vs. 1.5 g/kg). The two fields have similar shapes and both analyze a low-high couplet in southeastern Colorado, placing these features in the same location. Both analyses appear qualitatively to possess the same limits on high wavenumber noise.
This comparison shows the result that HSM offers and demonstrates that the 150-iteration solution approximates the more computationally expensive 2000-iteration solution adequately to be used operationally.
The primary purpose of this examination is to detect whether GOES water vapor gradients coupled into the LAPS moisture analysis through HSM have a measurable positive influence. To achieve this goal, statistics compare different sets of analyses from April 1993 to May 1994. The basic test statistic was the computed difference between the Denver Radiosonde Observation (RAOB) dewpoint measurement and the dewpoint derived from the LAPS specific humidity analysis [td (LAPS) - td (RAOB)] at all levels. We focus on two time periods referred to as t0 and t-1. t0 refers to the RAOB time, e.g., 1200 UTC and 0000 UTC, and t-1 refers to the hour preceding sounding time, e.g., 1100 UTC. The reason that one would expect this to be a ``fair" test is that both t0 and t-1 LAPS analyses do not contain RAOB data from that hour. Furthermore, soundings are typically launched in the t-1to t0 interval, often as early as t-1. From an analysis standpoint alone (where both t0 and t-1LAPS analyses were computed in the same way) one would expect that the t-1 would have an advantage, since the background field for LAPS one hour earlier would be based on a forecast field closer to initialization time. In our case we used MAPS to initialize LAPS; for example, the 1100 UTC LAPS analysis would use a MAPS 2-h forecast initialized at 0900 UTC and the 1200 UTC LAPS analysis would have used a MAPS 3-h forecast initialized at the same time for its background.
We subdivide these two datasets further by paying special attention to one set: the one in which the t0 LAPS analysis has undergone HSM satellite modification and the t-1 set has not. As a control, we also have at our disposal another set in which satellite data were not used at t0. In this analysis, we can make one further classification between clear and cloudy scenes. The LAPS cloud analysis gauged cloudy conditions near the Denver RAOB site. If clouds existed anywhere in the column at the LAPS grid point corresponding to the Denver RAOB site, the sample was classified cloudy.
* satellite modified at time t0 only
** no satellite modification at any time "control"
a. The clear situation
Table 3 shows the results obtained from the different sets of runs. Note first the figures from the ``control" group. Comparing the clear bias error for roughly 2000 observations, both the t0 and t-1 distributions had roughly the same bias of about 0.6~K, with standard deviations of about the same values as well, 4.1 to 4.2~K. This indicates that these distributions are the same. We also see a slight improvement in bias error and $\sigma$ at t-1, but the difference is subtle and not worth dwelling on except to speculate that it might be associated with the arguments presented earlier regarding the quality of the background field with respect to the MAPS forecast analysis time.
Looking now at the satellite modified t0 and t-1 values, we note that the t0 value receiving the satellite bias adjustment has an improved bias of 0.04~K versus the t-1 value of 0.27~K. This is a relative improvement of 83\%, but is it significant? To address this, two tests were performed on the distributions. First was a Student-t test that was applied in two different ways. One way assumed that the variances of each distribution were distinctly different. The Student-t statistic for this case was 2.82 representing a confidence level of 99.5\%. Assuming that the variances were the same, the results were nearly identical. Therefore, we may assume that, based on the Student-t test, these are separate distributions, truly showing a significant improvement.
The Kolmogorov-Smirnov test (Gibbons 1976) was also performed on the data. This approach simultaneously considers differences in variance and value. Here the probability of d (the computed statistic) being greater than observed was 0.019, and for this test, a strong significance level is indicated by 0.01 or smaller. Here this criterion is almost met. Coupled with the Student-t result, we can confidently state that there are two different distributions and the satellite modification shows a genuine positive result.
b. The cloudy situation
Again we begin by looking at the cloudy subsets in the ``control" group. Here we see that there was approximately one-half the number in the sample at t0 versus t-1. The bias statistics differ slightly, with the bias in t0 slightly worse. Both of the sample standard deviations are greater than their clear counterparts. Again, as in the clear case comparison, we see no striking differences between the two time periods.
Now looking at the satellite modified cloudy statistics, we first see that in approximately 10% of the cases the HSM step executed when the LAPS cloud analysis showed clouds over the Denver RAOB site. This is very different from the ``control" group for the following reason. In the generation of the satellite water vapor fields, there is a check for the cloud fraction in the scene. If this is greater than 33\%, the satellite water vapor fields are not produced. So we are not looking at extensive cloudy fields here. Furthermore, we are probably not looking at super clear cases either, but cases in between. For one sample of about 140 measurements, representing about 14 soundings (estimating 10 levels per sounding), HSM is performed when the whole field is deemed partly cloudy, and the satellite water vapor fields are generated while Denver remains cloud covered. The bias statistic in this situation is about 3~K or about 6 times worse than the bias at t-1. This by itself suggests that the HSM analysis does not handle the cloudy situation well. Both Student-t and Kolmogorov-Smirnov tests show these to be distinctly different distributions and the HSM analysis as applied for LAPS seems to do more harm than good in cloudy situations.
A host of reasons could contribute to this behavior of the analysis, therefore it is very difficult to sort out just what is happening here. The way LAPS treats clouds relative to the moisture analysis is complex. In a partly cloudy situation, it may be that the regression coefficients responsible for the satellite precipitable water fields are not as good, leading to erroneous gradients near clouds. Another possibility may be poor parameterization of RH as a function of cloud amount. Another possibility is that the analysis saturates the air in cloudy regions, and with a limited vertical resolution of 50~hPa a thin cloud may be analyzed as a much thicker entity. Many of these arguments also hold for the ``control" group cloudy subset as well; we must keep in mind that the control group is not as representative of the meteorological conditions here as it is in the clear situation. The major difference between the control and satellite modified subsets is that the control includes both partly cloudy and complete overcast situations; whereas the satellite modified group {\it only} includes partly cloudy situations. This latter fact points to other possibilities such as problems handling gradients in the vicinity of clouds or perhaps the exact location of clouds that might be caused by navigation errors that often occur in geostationary imagery. Satellite navigation error on the order of 5~km could be very detrimental to a mesoscale water vapor analysis in a partly cloudy situation. We also are assuming the that RAOB balloon goes straight up and does not drift horizontally. In a partly cloudy situation it is likely to drift from inside a cloud to a region outside the cloud. This action would also lead to a positive bias in LAPS as observed.
We conclude that the HSM gradient insertion contributes positively to the LAPS three-dimensional water vapor analysis. This is only true of situations in which we have clear skies analyzed over the verification area; HSM analysis in its current form does not perform well in cloudy situations. An alternate approach might be to improve the discrimination between clear and cloudy situations and skip HSM application for cloudy cases.
An effective method to determine HSM coefficients using analytical techniques was demonstrated. The corresponding response of the numerical solution was shown to be in agreement with the analytic theory. The statistical results show a positive effect when applying HSM.
An effective approximate numerical scheme using fewer integration steps was demonstrated to give similar results when compared to the ``ideal" solution. Fourier techniques are very effective when establishing the behavior of different sets of HSM coefficients.
It is interesting that the influence of HSM reduced the bias but did not significantly change the variance. Perhaps one way to interpret this evidence is to think of the HSM as pushing the analysis in the right direction perhaps without enough small-scale influence. Could the variance of dewpoint on a spatial scale of 65~km be the limiting factor? By reducing the spatial scale we could possibly sense an effect below this limit. Both an increase in the gradient influence and a broader (and more intense) gradient response could be achieved by accepting a lower gradient cutoff wavelength (i.e., less then 65~km). The current gradient cutoff is conservative even by GOES-7 standards; soon the algorithm will be modified to use GOES-8 data. Improved signal to noise will result in less spatial averaging. This, along with the better registration anticipated with GOES-8, should allow a reduction in the cutoff to perhaps 40~km or less.
Acknowledgments. Partial funding for this work was made through NESDIS as part of the GOES I-M Product Assurance Plan. The author is grateful for the advice in the peer and technical reviews both internal and external to FSL.
Achtemeier, G.L., H.T. Ochs III, S.Q. Kidder, and R.W. Scott, 1986: Evaluation of a multivariate variational assimilation of conventional and satellite data for the diagnosis of cyclone systems, Preprints, Second Conference on Satellite Meteorology/Remote Sensing and applications, Williamsburg, VA, American Meteorological Society, 76-81.
Albers S., J. McGinley and D. Birkenheuer, 1994: LAPS analyses of clouds and precipitation. Preprints, 6th Conf. on Mesoscale Processes, Portland, OR, Amer. Meteor. Soc., 158-161.
_______, 1995: The Local Analysis and Prediction System (LAPS): Analyses of clouds, precipitation, and temperature. submitted to Weather and Forecasting.
Benjamin, S.G., 1991: Short-range forecasts from a 3-h isentropic-sigma assimilation system using ACARS data. Fourth AMS International Conference on Aviation Weather Systems, 24-26 June 1991 Paris, France, 329-334.
Birkenheuer D., 1993: The meso-$\beta$-scale analysis of water vapor using LAPS. Preprints, Thirteenth Weather Analysis and Forecasting Conference, Vienna, VA, American Meteorology Society, Boston, 204-208.
________,1992a: The LAPS specific humidity analysis. NOAA Tech. Memorandum ERL FSL-1, 39 pp. [Available from NOAA, Forecast Systems Laboratory, 325 Broadway, Boulder, CO 80303].
________, 1992b: Three-dimensional water vapor analysis using VAS. Preprints, Sixth Conference on Satellite Meteorology and Oceanography, Atlanta, American Meteorological Society, 301-304.
________, 1991: An algorithm for operational water vapor analyses integrating GOES and dual-channel microwave radiometer data on the local scale. J. Appl. Meteor., 30, 834-843.
Chance, B.A., 1986: Application of satellite data to the variational analysis of the three dimensional wind field, NASA Contractor Report 4022, Illinois Univ., Urbana-Champaign, 96 pp.
Gibbons, J.D., 1976: Nonparametric Methods for Quantitative Analysis. Holt, Rinehart and Winston, 463 pp.
Guiraud, F.O., J.~Howard, and D.C.~Hogg, 1979: A dual-channel microwave radiometer designed for measurement of precipitable water vapor and liquid. IEEE Trans. Geosci. Electron., GE-17, 129-136.
Hayden, C.M., 1988: GOES-VAS simultaneous temperature-moisture retrieval algorithm. J. Appl. Meteor., 27, 705-733.
McGinley, J.A., S.C.~Albers, and P.A.~Stamus, 1991: Validation of a composite convective index as defined by a real-time local analysis system. Wea. and Forecasting, 6, 337-356.
________, 1984: Meteorological analysis using the calculus of variations (variational analysis). Rivista di Meteorologia Aeronautica, XLIV, 37-44.
Sasaki, Y., 1970: Numerical variational analysis with weak constraint and application to surface analysis of severe storm gust. Mon. Wea. Rev., 98, 899-910.
________, 1958: An objective analysis based upon the variational method. J. Meteor. Soc. Japan, 36, 77-88.
Figure 1. Response function showing the ratio of the analyzed gradient amplitude to its input using the operational coefficients. The function plotted is computed from $\alpha$~= 64160~km$^4$, $\beta$~= 742.5~km$^2$. Also shown is the 0.5 ``cutoff response level." The vertical line indicates a wavelength of 63~km, the high frequency cutoff used for satellite gradient information. back to text
Figure 2a. A hypothetical ``gradient" field with a wavelength of about 105~km. back to text
Figure 2b. A hypothetical ``background" field with a wavelength of approximately 157.5~km. back to text
Figure 3a. Two-dimensional plot of the Fourier ``power spectra" after combining 2a and 2b via HSM. The four peaks represent the resulting power at different wavelengths. There are 2 sets each with a pair of peaks representing a background and gradient. The tall peaks represent the gradient signal with a response level of approximately 0.503, and the lower peaks are the background response with a level of 0.343. back to text
Figure 3b. Plot of the analytic gradient response (denoted as the ``optimal" response) and the numerically obtained gradient response of the standing waves (plotted points) using the approximate method. back to text
Figure 3c. Two-dimensional plot of derived gradient response (from Fourier methods) as a function of standing gradient and background wavenumbers produced from the approximate method. The plot is oriented such that viewing from left to right is similar to Figure~3b. back to text
Figure 4. Analyses of precipitable water (cm) from GOES-VAS after Birkenheuer 1991 for 12 October 1994. Satellite derived total precipitable water (4a) is shown with layer precipitable water fields: surface to 780~hPa (4b), 780-640~hPa (4c), and above 640~hPa (4d). back to text
Figure 5. LAPS analyses of specific humidity (g/kg). The first frame (5a) shows the analysis before the HSM algorithm, essentially the background; 5b shows the HSM analyzed result with the coefficients developed in the text (true solution); 5c shows the HSM analyzed result using the approximate solution using fewer iterations. Note the ``packed" gradients around clouds. Cloud data force high levels of specific humidities up to saturation. Also note the invariance in the boundary in 5b and 5c, a result of fixed boundary condition back to text
