RL926 MICROWAVE RADIOMETRY OF SNOWCOVERED GRASSLANDS FOR THE ESTIMATION OF LAND-ATMOSPHERE ENERGY AND MOISTURE FLUXES John F. Galantowicz September, 1995 RL-926 = RL-926

MICROWAVE RADIOMETRY OF SNOW-COVERED GRASSLANDS FOR THE ESTIMATION OF LAND-ATMOSPHERE ENERGY AND MOISTURE FLUXES by John F. Galantowicz A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy (Electrical Engineering and Atmospheric, Oceanic and Space Sciences) in The University of Michigan 1995 Doctoral Committee: Professor Anthony W. England, Chair Professor Sushil K. Atreya Professor Henry N. Pollack Assistant Professor Kamal Sarabandi Professor Fawwaz T. Ulaby

..... 1.1-1 i!............ I..........

( John F. Galantowicz 1995 All Rights Reserved

To Sara 11

ACKNOWLEDGMENTS I thank the committee members for their insightful comments and suggestions incorporated into this thesis and for their generous support for my work. This work owes a great deal to the guidance of my advisor, Prof. Tony England, who contributed many of the best ideas that tie my thesis together and has always been willing to let curiosity lead. I am grateful to the many individuals who volunteered their services and expertise in ways that have enhanced this work. I thank Dr. Dave Meyer of the USG(S EROS Data Center for his hospitality and support of my experimental work in Sioux Falls, South Dakota, and the staff of EROS who made the REBEX-1 field study possible. Thanks goes to the staff of the Marshall Space Flight Center's Distributed Active Archive Center in Huntsville, Alabama, for providing SSM/I data. I would especially like to thank Mr. Ron Hartikka for his assistance and ideas in the laboratory. I also thank the students who helped me as a part of their own studies: Ms. Cynthia Ballard, Mr. Howard Chang, Mr. Stan Modjeski, and Mr. Eric Beaubien who assisted in the production of SSM/I images appearing in the thesis. The intellectual stimulation, camaraderie, and hard work of my current and former fellow graduate students have both made my work less difficult and lifted my spirits for six years. I thank Dr. Richard Austin, Mr. Paul Dahl, Ms. Jasmeet Judge, Dr. John Kendra, Mr. Ed Kim, Dr. Brian Kormanyos, Dr. Joe Landry, Mr. Yuei-an Liou, 111

Dr. Adib Nashashibi, Dr. Leland Pierce, Mr. Paul Siqueira, Mr. Jim Stiles, and Dr. Brian Zuerndorfer for their technical advice and assistance. I also thank them and the many other Radiation Lab students who have been great friends and willing collaborators during this project. For their financial support, I thank the Office of Naval Research, which funded the first three years of my post-graduate work, and also the National Aeronautics and Space Administration and the National Science Foundation. I also thank the Rackham School of Graduate Studies for providing funds for travel to technical conferences. Finally, I would like to thank my wife, Sara, who, through her faith and encouragement, has been my fundamental inspiration. iv

TABLE OF CONTENTS DEDICATION................................................. ii ACKNOWLEDGMENTS........................................ iii LIST OF TABLES............................................. viii LIST OF FIGURES................................x LIST OF APPENDICES....................................... xvi CHAPTERS 1 INTRODUCTION....................................... 1 1.1 Linking snowpack, atmosphere, and radiometry............. 1 1.2 Approach and premise of the thesis............................. 3 1.3 Questions addressed by this dissertation..................... 6 1.4 Format of the dissertation.................... 6 2 BACKGROUND....................................... 8 2.1 Introduction........................... 8 2.2 Great Plains snowpacks...................................... 8 2.3 Microwave radiometry of snowpacks............................ 9 3 SIMULATING THE DYNAMIC SNOWPACK............. 12 3.1 Introduction............................. 12 3.2 Overview of Snowflow................................... 13 3.3 Heat flow in soil........................................ 16 3.3.1 Liquid water content in sub-freezing soil.............. 18 3.3.2 Enthalpy change in freezing soil........................ 19 3.3.3 Soil thermal conductivity.............................. 19 3.4 Snowpack formation and heat and moisture fluxes...... 24 3.4.1 Heat and vapor transfer within the snowpack..... 25 3.4.2 Snow-soil interactions........................27 3.4.2.1 Handling an insulating layer of grass.... 28 v

3.4.3 Shortwave attenuation and absorption............... 28 3.4.4 Snow-atmosphere interactions.................... 29 3.4.4.1 Radiative energy exchange................ 30 3.4.4.2 Latent and sensible heat exchange........... 31 3.4.4.3 Heat from precipitation..................33 3.4.5 Liquid water seepage in the snowpack............... 33 3.4.6 Compaction and mass balance of snowpack layers.. 36 3.4.7 Metamorphism of snow grains............. 38 3.5 D iscussion............................ 39 4 MODELING RADIOBRIGHTNESS OF THE SIMULATED SNOWPACK..................................... 40 4.1 Introduction........................... 40 4.2 Conventional radiative transfer theory................ 40 4.3 Solution for the simulated snowpack.............. 44 4.3.1 The scattering phase matrix and Mie scattering... 48 4.3.1.1 Attenuation coefficients: Empirical scattering reduction...................... 52 4.3.2 The method of invariant imbedding........... 55 4.3.3 Air-snow and snow-soil boundaries.................60 4.3.4 Dielectric properties of water, ice, snow, and soil... 61 4.4 D iscussion............................ 64 5 GROUND BASED RADIOBRIGHTNESS OBSERVATIONS IN THE NORTHERN GREAT PLAINS: THE FIRST RADIOBRIGHTNESS ENERGY BALANCE EXPERIMENT. 65 5.1 Introduction........................... 65 5.2 A pparatus............................ 66 5.2.1 Micrometeorological Subsystem and system integration 66 5.2.2 Design of the microwave radiometers and the Tower Mounted Radiometer System.............. 70 5.2.2.1 Microwave radiometer calibration...... 73 5.3 Installation............................ 78 5.4 Experim ent log.......................... 84 5.5 Correcting radiobrightness errors................ 85 5.5.1 Removing out of range radiobrightness values..... 86 5.5.2 Sky reflector positioning errors............. 87 5.5.3 Revision of the day 403 37 GHz radiometer calibration 88 5.5.4 Estimation of actual sky radiobrightnesses...... 90 5.5.5 Calibration of the 85 GHz radiometer for days 347-403 94 5.5.6 An alternative calibration parameterization...... 97 5.5.7 Sensitivity of radiobrightness to antenna efficiency assum ptions........................ 99 vi

5.6 D iscussion............................ 99 6 MICROWAVE RADIOMETRY FROM SPACE: THE SPECIAL SENSOR MICROWAVE/IMAGER................. 100 6.1 Introduction............................. 100 6.2 Data from the SSM/I............................... 101 6.3 Compensating for atmospheric attenuation and emission... 103 6.3.1 Using the rawinsonde atmospheric profiles............104 6.3.2 Compensation without a priori information.........109 6.4 Comparison of SSM/I and REBEX-1 radiobrightnesses ~. 114 7 COMPARING MODELS TO OBSERVATIONS.............124 7.1 Introduction............................ 124 7.2 Model initialization and inputs........................ 124 7.3 Evaluation of the snowpack SVAT...................... 127 7.4 Radiobrightness comparisons.................. 133 7.5 Comments on snowpack structure, grain size, and wetness.. 149 8 CONCLUSIONS, CONTRIBUTIONS, AND RECOMMENDATIONS...................................158 8.1 Conclusions............................. 158 8.2 Contributions............................ 162 8.3 Recommendations.......................... 163 APPENDICES..................................... 166 BIBLIOGRAPHY.................................. 216 vii

LIST OF TABLES Table 4.1 Mie-Rayleigh scattering parameter comparison. Qsca is in units of m2. 49 5.1 Micrometeorological instruments and parameters. N/A indicates data not available................................................. 68 5.2 Microwave radiometer specifications. N/A indicates data not available. See Appendix B for radiometric resolution calculations............... 70 5.3 Estimated loss parameters.......................... 76 5.4 Calibration parameters used during REBEX-1 from day indicated to day of next calibration. Some 37 and 85 GHz parameters were later modified. See Sections 5.5.3 and 5.5.5................ 78 5.5 Summary of REBEX-1 site coverage conditions and hardware problems. 86 5.6 Corrected calibration parameters for the period from days 403 to 471. 90 5.7 Reflector efficiencies from (5.18), the number of clear sky profiles used, and the standard deviation of 7rrfl. Also shown are the percentage of TSKY values from (5.19) that were less than zero........... 93 5.8 Revised 85 GHz calibration parameters for days 347-403............ 97 5.9 Alternative calibration parameters................................ 98 6.1 SSM/I sensor specifications [59]. EFOV on earth is in km.......... 101 6.2 Average, e, and standard deviation, r7, of difference between SSM/I (TA and TTER) and REBEX-1 brightnesses. All brightnesses are h-pol. 107 viii

6.3 Comparison of TTER from the rawinsonde and dry standard atmosphere methods showing average, e, and standard deviation, oe, of the difference e = TTER(rawin.)-TTER(dry atm.). Also given are these statistics for the magnitude of the correction, that is e = TTER- TA. All brightnesses are h-pol.............................. 112 6.4 Statistics of difference between REBEX-1 and SSM/I terrain brightnesses in January. All brightnesses are h-pol. DR is the January dynamic range of the REBEX-1 brightnesses.............. 119 7.1 Comparisons of modeled snowpack variables to REBEX-1 measurements. A = model - observation, a is standard deviation, and DR is the dynamic range of the measured data over the test period....... 128 7.2 Comparison of modeled brightnesses (K) to h-pol. observations from REBEX-1 and h-pol. and v-pol. observations from SSM/I. A = model - observation, a is standard deviation, and DR is the dynamic range of the measured data over the test period...... 136 A.1 Snowflow parameters set by the user......................... 169 A.2 Snowflow empirical parameters and their sources................ 169 A.3 Other physical parameters used in Snowflow..................... 170 B.1 Calibration check results. TECCO is the temperature of the Eccosorb target at the time the TAP were measured. There is insufficient data to calculate 85 GHz radiometric resolution........................ 172 B.2 Soil moisture sampling data. Masses are in grams and include tare mass. Drying was first done at 70~C then at 105~C. Gravimetric moisture content is tabulated at both temperatures and volumetric moisture content is tabulated only at 105~C....................................... 175 B.3 REBEX-1 day to calendar date conversion chart................. 202 C.1 SSM/I operational data for the F08 and Fll satellite platforms [60].. 210 C.2 Parameters for the Ml and Mh EASE-Grid projections.......... 214 ix

LIST OF FIGURES Figure 3.1 Schematic of the Snowflow simulated snowpack.....................13 3.2 Schematic of liquid water seepage through the simulated snowpack (after [38]).................................34 4.1 Conventional radiative transfer......................41 4.2 Upwelling and downwelling brightness in a snowpack layer............ 44 4.3 Scattering geometry for an isolated particle.....................50 4.4 Particle scattering geometry in the snowpack coordinate system.... 51 4.5 (a) Normalized scattering coefficients and (b) the independent scattering reduction factor............................ 54 4.6 Brightness balance for the kth thin snowpack layer........... 56 4.7 Components of the layer scattering source term, T..............57 5.1 Interior of the trailer on site sheltering data acquisition and device control electronics and the Macintosh computer running the FluxMon HyperCard stack which controlled the experiment.................. 67 5.2 The TMRS-1 radiometer housing............................71 5.3 The 85 GHz radiometer. The 19 and 37 GHz radiometers have layouts which are comparable component by component............ 72 5.4 Microwave radiometer block diagram................... 74 5.5 View of the REBEX-1 site from the east........................79 x

5.6 View of the radiometer housing with its back cover removed. Inserted into the housing are (from the far side) the 85 GHz radiometer, the IR radiometer and video camera module, the center connector box, the 37 GHz radiometer, and the 19 GHz radiometer...... 80 5.7 Plan view of the REBEX-1 site............................ 81 5.8 Graduated PVC pipe with alternating 1.27 cm (0.5 in) black and white stripes used for gaging snow depth from video stills of the REBEX-1 site........................................................ 82 5.9 Insertion of the soil temperature probes. At the time of this picture, I had already inserted the probes in the side of the trench and refilled it, burying the 64, 32, and 16 cm probe cables. Cables leading to the 2, 4, and 8 cm deep probes protrude from the side of the trench.... 83 6.1 Comparisons of SSM/I antenna temperatures and REBEX-1 terrain brightnesses................................................. 105 6.1 (Continued from previous page.)................................. 106 6.2 Variation with time of the difference between TTER from SSM/I and TTER from REBEX-1............................ 108 6.3 Brightness temperatures upwelling from the atmosphere as a function of surface water vapor density...................... 111 6.4 Simulated antenna temperatures as a function of surface water vapor density given terrain brightnesses of 200, 235, and 270 K.......... 111 6.5 Comparisons between SSM/I terrain brightness estimates, TTER, made with a dry standard atmosphere (TA - dry atm.) and a rawinsonde atmosphere (TA - rawin. atm.). All brightnesses are h-pol.. 113 6.6 19 GHz terrain radiobrightness from SSM/I and REBEX-1.......... 115 6.7 37 GHz terrain radiobrightness from SSM/I and REBEX-1.......... 116 6.8 85 GHz terrain radiobrightness from SSM/I and REBEX-1........ 117 6.9 Difference between REBEX-1 and SSM/I terrain brightnesses..... 118 6.10 19 GHz SSM/I terrain brightness at vertical and horizontal polarization and the v-h difference...................................... 121 6.11 37 GHz SSM/I terrain brightness at vertical and horizontal polarization and the v-h difference...................................... 122 xi

6.12 85 GHz SSM/I terrain brightness at vertical and horizontal polarization and the v-h difference............................. 123 7.1 Schematic of the Snowflow and Esnow inputs and products.......... 125 7.2 Modeled and observed snow depths with wet model soil........... 128 7.3 Modeled and observed snowpack surface temperature with wet model so il..................................... 129 7.4 Modeled and observed soil temperature at 2 cm depth with wet model so il..................................... 13 1 7.5 Modeled and observed soil heat flux at 2 cm depth with wet model soil. Positive values indicate heat flow into the soil................132 7.6 Modeled and observed soil temperature at 2 cm depth with dry model soil..................................... 134 7.7 Modeled unfrozen soil water content for model soils with initial moisture volume fractions of 0.43 (top) and 0.20 (bottom) before freezing. 135 7.8 Modeled and observed 19 GHz h-pol. terrain brightness with wet model soil and ground-based observations....................... 137 7.9 Modeled and observed 37 GHz h-pol. terrain brightness with wet model soil and ground-based observations.................... 138 7.10 Modeled and observed 85 GHz h-pol. terrain brightness with wet model soil and ground-based observations.................... 139 7.11 Modeled and observed v-pol. terrain brightness at 19, 37, and 85 GHz with wet model soil and space-based observations from SSM/I..... 140 7.12 Modeled h-pol. terrain brightness at 19, 37, and 85 GHz with both wet and dry m odel soil............................. 141 7.13 Plots of modeled 19 GHz h-pol. terrain brightness, modeled soil surface unfrozen water content, observed 19 GHz h-pol. terrain brightness, observed 2 cm soil temperature, and modeled and observed snowpack depths.................................... 143 7.14 Modeled h-pol. terrain brightnesses with dielectric layering artificially enhanced..................................... 146 Xll

7.15 SSM/I 19 GHz v- and h-pol. brightnesses (top graph and overlays) compared to model results with wet soil, dry soil, and wet soil plus enhanced dielectric layering....................... 147 7.16 SSM/I 19 GHz v- and h-pol. brightnesses and their difference at Langdon, North Dakota................................. 148 7.17 Modeled and observed h-pol. 37 GHz terrain brightnesses and modeled total snowpack liquid water content. The arrows indicate modeled partial snowmelt events that correspond to observed brightness jumps. 150 7.18 Snow liquid water content profiles during a cycle of partial melt and refreeze................................. 151 7.19 Snow grain diameter profiles from the Snowflow model....... 152 7.20 SSM/I 19 GHz h-pol. images from February 4, 1992 at 23:24 UTC (top) and February 8, 1992 at 12:44 UTC. See text for description... 155 7.21 SSM/I 37 GHz h-pol. images from February 4, 1992 at 23:24 UTC (top) and February 8, 1992 at 12:44 UTC. See text for description... 156 7.22 SSM/I 85 GHz h-pol. images from February 4, 1992 at 23:24 UTC (top) and February 8, 1992 at 12:44 UTC. See text for description... 157 B.1 Summary of terrain apparent radiobrightnesses at 19, 37, and 85 GHz. 177 B.2 Summary of reflector-measured sky radiobrightnesses at 19, 37, and 85 GHz and thermal infrared sky temperature................... 178 B.3 Summary of reflector-measured sky radiobrightness and estimated sky radiobrightness at 19 GHz........................ 179 B.4 Summary of reflector-measured sky radiobrightness and estimated sky radiobrightness at 37 GHz........................ 180 B.5 Summary of reflector-measured sky radiobrightness and estimated sky radiobrightness at 85 GHz......................... 181 B.6 Summary of REBEX-1 global and net radiation, 10 m wind speed, subsurface soil temperatures, and soil heat flux at 2 cm depth..... 182 B.7 Summary of REBEX-1 relative humidity, rainfall per experiment cycle, thermal infrared surface temperature, air temperature, and soil temperature at 2 cm depth........................ 183 xiii

B.8 Summary of terrain apparent radiobrightness spectral gradients and TPR reference load gain factors...................... 184 B.9 Summary of REBEX-1 snow depths estimated from video stills of the site and from National Weather Service reports at Sioux Falls, SD... 185 B.10 October terrain (surface) apparent radiobrightnesses and reflector-measured sky radiobrightnesses at 19 and 37 GHz, net and global radiation, and wind speed at 10 m...............................187 B.11 October subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature........................ 188 B. 12 November terrain (surface) apparent radiobrightnesses and reflectormeasured sky radiobrightnesses at 19 and 37 GHz, net and global radiation, and wind speed at 10 m..................... 189 B.13 November subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature......................... 190 B.14 December terrain (surface) apparent radiobrightnesses and reflectormeasured sky radiobrightnesses at 19, 37, and 85 GHz, net and global radiation, and wind speed at 10 m........................ 191 B.15 December subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature......................... 192 B.16 January terrain (surface) apparent radiobrightnesses and reflector-measured sky radiobrightnesses at 19, 37, and 85 GHz, net and global radiation. and wind speed at 10 m........................ 193 B.17 January subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature....................... 194 B.18 February terrain (surface) apparent radiobrightnesses and reflectormeasured sky radiobrightnesses at 19, 37, and 85 GHz, net and global radiation, and wind speed at 10 m................... 195 B.19 February subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature....................... 196 xiv

B.20 March terrain (surface) apparent radiobrightnesses and reflector-measured sky radiobrightnesses at 19, 37, and 85 GHz, net and global radiation, and wind speed at 10 m.......................... 197 B.21 March subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature.................................... 198 B.22 April terrain (surface) apparent radiobrightnesses and reflector-measured sky radiobrightnesses at 19, 37, and 85 GHz, net and global radiation, and wind speed at 10 m.......................... 199 B.23 April subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature.............................. 200 C.1 SSM/I scan geometry.......................... 205 C.2 A section of the SSM/I swath...........................208 C.3 Antenna patterns: (a) 19 GHz V-pol. effective antenna pattern, (b) 37 GHz V-pol. interpolated to 19 GHz pattern at center of scan, (c) 37 GHz V-pol. interpolated to 19 GHz pattern between edge and center of scan, (d) 85 GHz V-pol. interpolated to 19 GHz pattern at center of scan, (e) 85 GHz V-pol. interpolated to 19 GHz pattern between edge and center of scan. The contours are at 3, 6, 12, and 24 dB......... 212 C.3 (Continued from previous page.)......................... 213 C.4 SSM/I 85 GHz h-pol. images from February 2, 1992 at 12:21 UTC demonstrating resampling to low resolution (top) and high resolution (bottom ) EFOV................................ 215 XV

LIST OF APPENDICES APPENDIX A SNOWFLOW PARAMETERIZATIONS............... 166 B DATA FROM REBEX-1....................... 171 C RESAMPLING SSM/I RADIOBRIGHTNESSES TO A COMMON GRID.............................. 203 xvi

CHAPTER 1 INTRODUCTION The brightness of microwave radiation emitted by snow varies measurably with the moisture, structure, and substrate conditions of the snowpack. This feature has been widely studied and exploited in mapping seasonal, alpine, and polar snowfields. This thesis employs the sensitivity of snowpack microwave emissions in a different way: as a tool in the determination of energy and moisture fluxes between snow-covered land and the atmosphere. 1.1 Linking snowpack, atmosphere, and radiometry The cornerstone and major contribution of the thesis is REBEX-1, the first Radiobrightness Energy Balance Experiment. REBEX-1 combined continuous measurements of terrain brightness at three microwave frequencies (19, 37, and 85 GHz) with simultaneous monitoring of micrometeorology. The microwave radiometers simulated the channels and observation angle of the space-borne Special Sensor Microwave/Imager (SSM/I) with the added advantage of continuity at a single locale. The experimental site near Sioux Falls, South Dakota, is typical of the northern Great Plains grasslands in climate and vegetation cover. The experiment lasted from October, 1992 through April, 1993 and spanned vegetation senescence, snowpack 1

2 formation and evolution, and spring thaw. REBEX-1 was the first experiment linking the radiobrightness of terrain to local weather over this length of time. Terrain radiobrightness (also called apparent radiobrightness) is the the combined intensity of microwave radiation emitted by and reflected off of the ground and ground-cover. The source of microwave emission is thermal radiation in the emitting medium transmitted through the surface and into the air. Weather conditions are linked to emission through their control of the moisture content and thermodynamic state of the terrain. Because the microwave dielectric properties of water are strong functions of frequency, temperature, and phase (liquid or solid), radiobrightness is a sensitive indicator of moisture variation. But the relationship is complicated because moisture content and temperature may vary throughout the column of ground cover and soil that constitute the distributed source of emitted radiation. The radiobrightness of snow-covered terrain is an extreme example of the distributed source phenomenon. The dielectric loss factor of ice is around three orders of magnitude smaller than that of water at SSM/I frequencies, and the low thermal conductivity of snow means that the temperature change from soil to air through the snowpack may be more than 30 K. Consequently, the effective temperature of the emitting source may have little relationship to the snowpack surface temperature. In addition, the inhomogeneity of the snowpack dielectric constant scatters and reflects radiation on its way to and from the snowpack surface, modifying both the emitted intensity and the reflectivity of the terrain. And when water is present in the snowpack at low volume fractions, it changes the snowpack into a dense, absorptive cloud with emissivity approaching unity. Formation and evolution of a snowpack are integrally tied to local meteorology, and

3 the link between snowpack radiobrightness and weather is similarly strong. Besides the obvious fact that snowpacks form from precipitation, the evolution of snowpack structure is a function of internal temperature and temperature gradients that are directly related to snowpack-atmosphere energy fluxes. And snowpack characteristicsfor example, low thermal conductivity, high shortwave and low longwave reflectivity, and high moisture availability in turn affect the mechanisms of land-atmosphere energy and moisture transfer. In seasonal snowpacks, melt and ablation (rapid vaporization) are the conclusion of this process and signal the transition to a new land-atmosphere transfer regime. 1.2 Approach and premise of the thesis This thesis presents two numerical models for use in analyzing REBEX-1 snowseason data. The first, called Snowflow, is a model of snowpack evolution and atmospheric interactions, and the second, called Esnow, calculates the radiobrightness of the Snowflow-modeled snowpack. Snowflow is among a class of numerical models for soil-vegetation-atmosphere transfer (SVAT) calculations. SVAT models are land-atmosphere interaction parameterizations that, when linked to an atmospheric model, provide a boundary condition for the vertical fluxes of latent and sensible heat at the bottom of the atmosphere. Snowflow differs from most SVATs in the level of detail used in modeling the surface medium. Atmospheric models are computationally intensive and must trade off the need for computation speed with vertical and horizontal spatial resolution. For example, the spatial resolution of Pennsylvania State University/National Center for Atmospheric Research Mesoscale Model version 4 (MM4) is 60 km, and GCMs (General Circulation Models) characterize the earth's surface in cell sizes ranging from

4 2.5~-10 in latitude and longitude [1]. At these levels of resolution, many processes in the atmosphere are reduced to sub-grid scale parameterizations-cloud formations, for example-and the high sub-grid heterogeneity of the earth's surface is best modeled by a parameterization scaled to fit the grid-sized needs of the atmospheric model. Simple SVAT parameterizations employ both the macroscopic reduction of ecotomes into basic classes and the microscopic reduction of soil and soil cover classes into abstract parameter sets [2] [3]. A parameterized description is not designed to match the temperature regime or physical structure of even an idealized terrain. Instead, terms of the parameterization correspond to the fluxes required by the atmospheric model and the parameterized information available from it. As a basic example, the bucket model of soil has a threshold water capacity-for example, 15 cm-beyond which runoff is generated [4]. Evaporation may be calculated as a function of the fullness of the bucket (soil moisture) and the potential evaporation. Although the model is based on physical arguments, the parameterization says little about conditions in the soil itself. In contrast to a simple SVAT, Snowflow must produce a comprehensive simulation of near-surface conditions matching the level of detail needed to estimate terrain radiobrightness. As discussed above, the distributed nature of the thermal emission source in a snowpack requires an emission model with detailed temperature and structural information, while the fluxes of moisture and heat between the snow surface and the atmosphere can be parameterized using only snowpack surface conditions. Yet the motivation for SVAT model simplicity-computational speed-remains if the SVATlinked emission model is to be used as the boundary condition for an atmospheric model. The independent variable of a SVAT-and the atmospheric model to which it is linked-is time, and the length of time modeled may be months and the tem

5 poral resolution may be just a few minutes. Although emission calculations are not necessary at all SVAT model times, the SVAT-linked radiobrightness model needs to evaluate emission two to six times per day to adequately capture the dynamics of land-atmosphere interactions. Consequently, for a combined SVAT-emission model to be useful as an experimental test-bed, its evaluation speed must be comparable to that of the SVAT model alone. A well-tested numerical analog to global or mesoscale meteorology would be a unique indoor laboratory in which climate predictions and sensitivities could be tested. But verification of an atmospheric model is hindered by the lack of long-term global scale meteorological data against which the simulated conditions could be compared. Similarly, any atmospheric model must be initialized with large scale data that accurately describes the state of the real system. In some cases remote sensing measurements have been used in both the initialization and verification stages of model development. For example, Dickinson et al. [5] initialized vegetation parameters for the Biosphere-Atmosphere Transfer Scheme (BATS) using both satellite remote sensing data and conventional maps when linking BATS to MM4. The BATS/MM4 model was tested against historical remote sensing data including seasonal snowpack depth. The premise of this thesis is that remote sensing may be used more intensively in mesoscale climate model verification by exploiting the satellite-measured radiobrightness records of large-scale seasonal snowpacks and comparing them to landatmosphere interaction model dynamics. The Great Plains seasonal snowpack is well suited to this use because it both responds to meteorological conditions through structural and thermal changes and is a direct product of meteorology. The Snowflow-Esnow model described here presents a simplified but not simple attempt to verify the concept of a SVAT-linked radiometric emission model. The REBEIX-1 data

6 set serves as a surrogate atmospheric model while providing frequent measurements of terrain radiobrightness for comparison to the models. 1.3 Questions addressed by this dissertation This thesis addresses the following questions: * How does emission from snowpack-covered terrain respond to long-term atmospheric conditions? * What characteristics of the Great Plains snowpack and its substrate soil must be modeled to simulate dynamic emission accurately? * Under what conditions do microwave radiometric measurements from space correlate to those made from ground-based instruments? * What SVAT processes have strong enough radiobrightness signatures that radiometric measurements may be used to monitor them? * How can measurements of microwave radiobrightness be used in conjunction with SVAT-linked emission model predictions to improve the SVAT model simulation? 1.4 Format of the dissertation The dissertation is organized based on the chronology of SVAT-linked emission model development and testing. Chapter 2 first covers background material regarding snowpack radiometry and land-atmosphere transfer. Then Chapter 3 describes Snowflow, the snowpack thermal simulation model. Snowflow's critical parameterizations, the governing equations, and their solutions are developed. Appendix A

7 includes several additional subordinate parameterizations and the model's fixed input parameter set. Chapter 4 describes the model of snowpack microwave brightness called Esnow. The chapter explains how Esnow uses the detailed snowpack simulation from Snowflow to calculate snowpack emission and reflection. The experimental work of the thesis is presented in Chapters 5 and 6. Chapter 5 explains the REBEX-1 apparatus, methodology, and post-experiment data processing. Graphical presentation of REBEX-1 data is included as Appendix B. Chapter 6 describes data from SSM/I, comparing observations from the satellite to those from REBEX-1. Processing methods for the SSM/I data are included in Appendix C. Chapter 7 brings together theory and observation. The chapter discusses the accuracy of the Snowflow snowpack simulation when it is driven with atmospheric data from REBEX-1. Esnow emission estimates are compared to both the ground-based and satellite-acquired radiobrightnesses over a 55 day period. Chapter 8 summarizes the conclusions, implications of the results, contributions of the thesis, and recommendations for future work.

CHAPTER 2 BACKGROUND 2.1 Introduction This chapter briefly covers two background topics: (a) the role of Great Plains snowpacks in the climate system and (b) research on the microwave radiometric properties of snow. The discussion is meant to place this thesis in the context of previous work, lay the groundwork for the analysis of modelled and observed radiobrightness in Chapters 6 and 7, and motivate the application of land-atmosphere transfer modeling in microwave remote sensing. 2.2 Great Plains snowpacks Mid-continental regions are climatically more sensitive to the land-atmosphere boundary condition because the influence of oceanic weather systems decreases away from coastal regions. In particular, the Great Plains of Canada and the United States straddle the boundary between the wet east, which receives moisture from the Gulf of Mexico, and the dry west, with less than 50 cm of precipitation per year [1]. Wintertime snow usually covers the northern Great Plains and is a factor in springtime water availability. The presence of a snowpack delays springtime warming through heat absorption during melting, and in winter snow dramatically changes the thermal balance-primarily through its high albedo. In addition, the smooth 8

9 snowpack alters the aerodynamic roughness by burying vegetation, decreasing the efficiency of energy transfer. The net climatic effect of the seasonal snowpack is of interest in climate modeling because of the feedback response the snowpack may have under climate change scenarios [4]. Snow has a high albedo and its presence will decrease the absorption of incident solar radiation. Snow's thermal infrared emissivity is also high, leading to faster cooling of the surface during cold, clear-sky conditions and more efficient warming when the sky is radiometrically warm. Consequently, the climatic effect of snow cover may be either net cooling or warming: higher albedo and infrared emissivity means an increase in radiative cooling under clear-sky conditions but high IR absorptivity reduces radiative cooling under cloudy skies. Depending on cloud conditions-an uncertainty in climate change models current GCMs indicate that the net seasonal effect of snowcover on regional climate may be weakly positive to negative net cooling [6]. 2.3 Microwave radiometry of snowpacks Observations of snowpack radiobrightnesses have been conducted with groundbased, airborne, and satellite radiometers. Early ground based studies showed that although radiobrightness is linked to the hydrologically important snowpack parameters of depth, water equivalent, and wetness other factors must be taken into account [7] [8]. The fundamental conclusion of early studies and the basis for further work has been that radiobrightness is reduced dramatically by the presence of a snowpack over bare ground [7][8][9][10]. Saturation of this trend has usually been shown to occur near the water equivalent depth of 25 cm (that is. the depth of water which would be measured if the snowpack were completely melted). Under some circumstances

10 radiobrightness has been shown to increase with water equivalent above 25 cm in highly evolved snowpacks [11]. In a study in which consecutive layers of a snowpack were removed, it was found that the presence of a 30 cm depth hoar layer contributed to the bulk of radiobrightness reduction in snowpacks that were from 64 to 83 cm thick [12]. Dry snow is primarily a mixture of air and ice which has low dielectric loss at microwave frequencies (1-100 GHz) compared to soil or water. Consequently, emission depth in dry snow is large-from 10 to 100 wavelengths [13]. The contribution to radiobrightness of the medium underlying the snowpack is therefore significant except with very thick or wet snow layers [7][14][10]. Since emitted radiation originates at significant depths within the snowpack, emissivity cannot be measured directly due to the difficulty of defining the thermometric temperature of the source. For thin snowpacks, the substrate temperature has been used [12] as well as the average temperature of the snowpack [10] to define effective emissivity values so that apparent snow radiobrightness could be corrected for the variable contribution of reflected downwelling sky brightness. Volume scatter darkening effects in the snowpack are a consequence of large penetration depths. Since snow grain sizes range typically from 0.1-5.0 mm and may grow up to more than 3 cm in depth hoar [15], scattering effects are increasingly large for frequencies over 10 GHz. The presence of water in the snowpack surface layer reduces emission depth and volume scattering effects [16][17][9]. The resultant absorbing "cloud" has a high emissivity especially at the higher frequencies and its radiobrightness is easily distinguishable from that of dry snow. Diurnal melt/freeze cycles can be recognized by the contrast between appropriately timed radiobrightness observations [11].

11 Experimental studies of microwave emission from snowpacks and the availability of data from spaceborne radiometers have led to the practical use of satellite radiobrightnesses in monitoring snowpack parameters (for example, [18]). There have been considerable obstacles to corroborating radiometrically derived snow parameters from satellite data including low measurement resolutions, geolocation errors, slope effects, vegetation effects, and atmospheric interference [19] [20] [21] [22] [23 1. Several studies of the correlation of radiobrightness or spectral gradient to snow depth or snow water equivalent have found that good statistical agreement can be found but only when seasonal or regional variation is removed [24] [14] [22]. The large contrast in radiobrightness between wet and dry snowcovers has been observed from satellite platforms as well [25] [22] [26]. The mapping of hemispherical snowcover by radiometry has been achieved based primarily on the negative spectral gradients (usually a function of TB(37 GHz)-TB(19 GHz)) characteristic of dry snowpacks [27] [28][29][30]. Data from the SSM/I have only recently become available but there has been some investigation of the use of the higher frequency 85.5 GHz channel for retrieval of snow parameters [23] and classification of snowcover [28][23]. Several studies have noted the need for detailed ground truth measurements of not only snow water equivalent, depth, moisture, and underlying soil conditions but also crystal size distribution, stratification, and temperature structure [31][11][32]. Techniques have been developed for retrieving detailed descriptions of snow grain properties [33] but have not been widely used in microwave studies because they are difficult to master and time consuming [34]. Researchers have suggested the use of snowpack emission models which track grain size as a function of the meteorological conditions which control metamorphism, but none have been published to date [11][32].

CHAPTER 3 SIMULATING THE DYNAMIC SNOWPACK 3.1 Introduction This chapter presents a numerical model of snowpack development, temperature, moisture movement, and atmospheric interactions called Snowflow. Snowflow is a SVAT (soil-vegetation-atmosphere transfer) model whose snowpack simulation includes information necessary for radiobrightness calculations-that is, far more information than is needed for energy and mass transfer alone. Although Snowflow is a new model, I have borrowed many of its formulations from other approaches. The solution method and most of the parameterizations are taken from Anderson [37]. Formulations for snowpack water seepage, compression, and snow-grain growth are primarily from the more comprehensive SNTHERM.89 model by Jordan [38], who also draws on the earlier Anderson work. Snowflow differs from these models primarily in its treatments of (a) soil thermal conductivity and unfrozen water content (section 3.3), (b) latent and sensible heat exchange (section 3.4.4.2), and (c) an insulating grass layer (section 3.4.2.1). Where necessary, Snowflow's parameterizations are designed to simulate conditions specific to the grasslands of the northern Great Plains. Consequently, several sections in this chapter refer to the REBEX-1 experiment site, described in more detail in Chapter 5. 12

13 Qsw Qiw Qemit Qprcp Qair Qlatent Air N s................... Snowpack layers Grass 0 layer Soil layers Ng ----------------------------------------------------------------- Qg=O Figure 3.1: Schematic of the Snowflow simulated snowpack. 3.2 Overview of Snowflow Snowflow's model components fall into four broad categories: snowpack energy, snowpack structure, soil temperature, and atmospheric interactions. Snowpack energy (enthalpy) change includes temperature changes and-if the snow is at the freezing point-phase changes. Snowpack structure includes moisture movement by gravity, density changes, snow grain growth, and snowpack development from precipitation. Snowflow can run without a snowpack in which case only soil temperature and atmospheric interactions are active while the model checks for new snowfall. Figure 3.1 shows the conceptual elements of Snowflow. The surface flux boundary condition is Snowflow's link to an atmosphere that drives the simulation and determines its particular locality and time. In Snowflow's current configuration the atmosphere is independent and may be modeled on historical data or be the result of

14 an atmospheric model. With minor modifications, Snowflow could also interact with an atmospheric model as it operates, providing a dynamic boundary condition for the model atmosphere. The following energy boundary condition constrains the total enthalpy change of the snow and soil, AQ, in a time interval, At: Qnet = Qsw + Qlw - Qemis + Qair + Qlatent + Qprcp + Qg = AQ (3.1) where Qnet is the net heat added to the snow/soil system per unit area, Qsw is heat from solar (shortwave) radiation, Qlw is from atmospheric (longwave) radiation, Qair is from sensible heat transfer with the air, Qiatent is from latent heat transfer with the air, Qprcp is from heat exchanged with precipitation, and Qg is from deep in the ground. Snowflow defines a thermally active region that extends down to the point where the temperature gradient falls below a given criterion-effectively assuming that Qg = 0. By conservation of energy, the change in stored heat per unit area is given by the sum of enthalpy changes in all the layers to the depth where Qg = 0: N, Ng AQ = A.HA + E AHgi (3.2) i=o i=0 where AH. is the enthalpy change per unit area of the ith layer of medium x, and x is either s (snow) or g (ground). Alternatively, a lower boundary could be set at a depth at which Qg could be estimated from historical data or a multi-year model. Snowflow divides both the soil and snowpack into discrete layers and each layer interacts only with its immediate neighbors. For the ith layer of snow or soil we can define an equalization function, Exi, based on one-dimensional heat transfer: Exi(Ux,i+l UX ) - A x, - QHi (3.3) where Qxi is heat added to the layer, uZ is the unknown variable describing a layer's thermodynamic state, x is s for snow or g for soil (ground). The vector of solutions,

15 I!I' - - U, satisfies Ez,i(U) = 0 (3.4) for all layers. In soil, the thermodynamic variable is always temperature, T, but in snow it may be either T or liquid water content, W (expressed as a depth). In the later case, the snowpack temperature is known and is equal to the freezing point of water, To. The top and bottom layers of the snow-soil system are special cases of (3.3). For the top layer (either snow or soil), Snowflow uses the net energy balance (3.1) to constrain the thermal state: E,top(U) -Qx,net(x,top) - AQ(u) (3.5) where u is the vector of the unknown state variables for all the layers. The bottom soil layer is always constrained by the condition Qg = 0 which implies that Tg,Ng+l = Tg,Ng. Consequently, for the bottom soil layer: Eg,bot = Eg,bot(Ui, Ui-) - AHIbot - Qg,bot (3.6) Snowflow solves the set of soil and snow energy equations simultaneously using the Newton-Raphson iteration technique. First, (3.3), (3.5), and (3.6) are expressed by a set of first order multivariate Taylor series:1 Ei(u) = Ei(Ua) + E (Ua)u +.. (3.7) where E.T is the transpose of the vector whose jth element is: E( ) dE (u) (3.8) m s u=Ux 'The medium subscript, x, is now implicit.

16 Ua is an estimated solution for the unknown vector and Au = u - Ua. If U is the exact solution then Ei(U) = 0 and: E (Ua)AU = -E(Ua), (3.9) where the (i, j) element of the matrix E is the derivative E'j. Because all but the top -_! layer equation are dependent on only three unknowns, E (Ua) is a singly bounded band-diagonal matrix. At a particular time step, t + At, Snowflow solves for Ut+'t by taking Ua = Ut initially, solving (3.9) for AU by decomposition and forward- and back-substitution, and determining a new estimate of the unknown vector Ut+t = Ut + AU. The process is repeated until all of the elements of |AU| are less than a threshold value. When the unknown for any snowpack layer changes from T to W, the process is automatically repeated until the solution for the new unknown converges. The remainder of this chapter details the parameterizations that Snowflow uses to specify (3.3), (3.5), and (3.6) for snow, soil, the snow-soil boundary layers, and the boundary with the atmosphere. Appendix A includes additional formulas for physical variables and tabulated values for Snowflow parameters. 3.3 Heat flow in soil Snowflow models soil using the one-dimensional heat flow equation [39]: -(T)_ = ~d (Kisol(T z) (t) (3.10) where T(z,t) is temperature at time t and depth below the surface z, HV(T) is enthalpy per unit volume (J/m3), and IK,,soi(T, z) is soil thermal conductivity (W/mK). We want to derive a heat transfer equation in the form of (3.3). Expanding the depth derivative, we have: aHV(T, z) aKso01 &T 02T at = z az + K~ (z2.1 V(3.11)

17 We then apply (a) an implicitly formulated discrete time approximation using the average heat exchange at times t and t + At, and (b) a finite difference approximation for a layer of thickness dg,. Then the equalization function (3.3) for soil (ground) is: Eg,i = AH. - Qg,i soil it T ^r^ t a 2 \ t =AHA 0.5 At dgi [+ (T ~ 'i aZ ) Z z K+t+A t T +( so il (T t+ & z J +s l ) &z2 /i 1 (3.12) Snowflow uses Egi for the intermediate soil layers in (3.9). The top soil layer is unique because of its interaction with the snowpack and the addition of heat from shortwave radiation. Section 3.4.2 deals with these questions after a description of snowpack heat flow has been formalized. Snowflow approximates the depth derivatives in (3.12) and elsewhere using the finite difference approach: (v 0.5 + i+ - (3.13) k\z Z - Z_-1 Z+1l - Zi ( 2 (.0 _ V- V (3.14) az2 Zi+l - Z- Zi+l - Zi Zi - Zi-l where V is the function of depth whose derivatives are to be found and zi is the midpoint depth of layer i. With this expansion, the derivatives of Egi with respect to the vector of soil unknowns Tg are easily found and applied in (3.9). For brevity, this set of derivatives will not be written out here. The following subsections describe Snowflow's parameterizations for soil enthalpy change and conductivity. In a freezing soil these variables depend primarily on the residual unfrozen water content, so we first turn to the characterization of this phenomenon.

18 3.3.1 Liquid water content in sub-freezing soil Experimental observations have demonstrated that soil water freezes over a broad range of temperatures below the freezing point of free water. Water in soil forms weak chemical bonds which inhibit the formation ice crystals. Some soil-water is bound more closely to solid particles and, in very moist soils, a portion of the soil water will behave like free water. The proportion of water bound in any particular state is dependent more on the soil type-and primarily its specific surface area than on the total water content. That is, given the soil type and adequate water availability, the amount of unfrozen water at some temperature below freezing is independent of the total water content. Unfrozen water content can be estimated from the empirical formula: xw if Tsoil > Tfpd, XU Pbu (3.15) pWb (To-T )OIwu otherwise, where xw is the volume fraction of water when the soil is above freezing, awu and i3wu (which is negative) are parameters of the soil type, Pb is the soil bulk (in situ) density, and Pw is the intrinsic density of water. Tfpd is the freezing point depression temperature, the lowest temperature below freezing at which xu = xw: / P 1/\ WQ Tfpd =To - ( Pw (3.16) CtwuPb I Experimental values of awu and,,3 have been tabulated for a wide variety of mineral soils but there is scant data on soils with high organic contents. Soils with a large amount of organic material are characterized by a slower decrease in the unfrozen water content with decreasing temperature than the more extensively investigated mineral soils [40]. Snowflow uses awu and AWU values for kaolinite soil from [41] because it also has a shallow freezing curve.

19 3.3.2 Enthalpy change in freezing soil The enthalpy change of a soil layer, A H. in (3.12) is the integral from Tt to Tt+/t of the apparent heat capacity of the soil layer per unit area, Ci: Tt+At AH. = T C.(T)dT (3.17) gt- t g, where C T) Pw(CpwXu{ + Cp'i(Xw - Xu)) + pbCpav - If u < T fpd C',(T) - (3.18),dgi(PwCpwXw + PbCpav) T > Tfpd. When T < Tfpd, the first term in parenthesis is the combined specific heat capacity of the ice and water content of the soil, Cpav is the average heat capacity of the dry soil constituents, and the last term is due to freezing/melting, where If is the latent heat of fusion and /3, < O. If both integral limits are below the freezing point depression temperature, then integrating (3.17) yields: Hg,i [(PbCpa + pwCpxw)(Tt+t - Tt) H' - ~,~d [(pbCp% + pCP )T+ T ) (Cp - cpi)pbawu (wu+ _ + l 0 3w+a w - Iwu -1 wu 1 2 - ) (3.19) where 01 =T- Tt and 02 = T- T t+/t. In general, (3.17) is integrated piecewise across the discontinuity at T = Tfpd. 3.3.3 Soil thermal conductivity Freezing soil is a multi-phase, multi-constituent matrix whose components have intrinsic thermal conductivities spanning three orders of magnitude. Appendix A lists intrinsic conductivities for the minerals, ice, water, and air which constitute the Snowflow model soil matrix. To calculate the matrix conductivity, Snowflow applies the de Vries method with modifications to include ice in the matrix [42][43] [44].

20 The general form of the de Vries conductivity mixing formulation is: Isoil = x00 + knxnAn (3.20) x0 + E knxn where Xn and A, are the volume fraction and intrinsic conductivity, respectively, of the nth constituent. The 0th constituent is generally the one that is most nearly continuous-for example, air in dry soil or water in very wet soil. The weighting factor, k=, relates the microscopic temperature gradient in the nth constituent with that in the background: k n= (az) (3.21) (O9T/9z)o' Assuming that the granules of constituent n are ellipsoidal and randomly oriented, de Vries gives the weighting factors as: kn=- E [1+ ( )](3.22) x=a,b,c L \OAJ where, for example, 1 ro du ga 3 J-abc (a2 + U)3/2(b2 + u)1/2(c2 + )/2 (3.23) g, depends only on the ratio of the ellipsoid axes (a, b, c), and ga + gb + gc = 1. For spherical particles, ga = gb = 9c = 1/3. Use of the de Vries conductivity model (3.20) requires the choice of an appropriate background medium and determination of the go for the remaining constituents. The de Vries model is fundamentally empirical and has been shown to work well when the g, are tuned for a particular soil type. For example, de Vries offers a correction factor for the dry soil case when using ga = gb = 0.125 for a mineral soil with dry air filling the voids: KIdry = 1.25 (oidvoid + knxnAn) (3.24) Xvoid + Ekxn 32

21 For a saturated soil, the voids are filled with water and: I Xsat XW + E k-xn7An (3.25) To find Ksoii over a range of moisture conditions, de Vries interpolates between the saturated and dry cases. The following discussion uses this method and extends it to include ice as well as air, water, and solids. In moist, freezing soil, there are three possible background media: air, water, and ice. But, as discussed in Section 3.3.1, water bound to the skeleton of the soil matrixthe soil grains-exists in soils to temperatures well below freezing. Consequently, we can take water to be the background medium in moist soils whether below freezing or not, being careful to check that in the limiting case of completely frozen soil the soil conductivity reduces to that of an air-ice-solids matrix. Beginning with limiting cases, the weighting factors for two types of solid soil constituents (to be specified later) in a background of constituent 0 are given by: 3= (1+ (A - I1) gp + (A ( g )) (3.26) 3= i(+ ( -)aP + I2ga+ Ip) ) (327) where ga,p is the principal size factor of the soil particles and the background medium is either air or water. To find the weighting function for air in a water background, kairXwU as a function of xw, consider a soil that is very nearly saturated with small air inclusions that are approximately spherical in shape. In this case, ga,ajr(xzL = xv) = 1/3 is the shape factor and the corresponding weighting coefficient is ka,w(^v), where X, is the volume fraction of voids. On the other extreme. in a very dry soil with air nearly filling the voids, water will coat the soil particles and the weighting factor for

22 water in an air background is: _ 1 _2 1 _ _ _ lim kw,a =- +. X1i mk 3 (1 + A i, ) gap 1 + (Ar ) ( - 2gap) (3.28) XA -ir A -ai (3.28) Equation (3.21) gives us a way to exchange the background medium with the nth constituent. Switching the background in (3.28) from air to water yields: lim ka = 1 xw-o ' lim=-_o kwa 1* 2 1 3 V + (Aai - )ga,air(O) 1+ ( r- 1)(1 - 2ga,air(0)) (329) We can then solve a quadratic equation for ga,air(xw = 0). Aair varies with temperature when the air contains water vapor which in turn is a function of the saturation state of the soil. In Snowflow, the conductivity of the soil air is divided in to dry, Aair,d, and vapor, \air,v, regimes according to: J Aaird Xf d < Xw < Xv Aair = (3.30) Aar,v(T) 0 < zX < Xfld where xfld is the so called field capacity of soil moisture and is given Appendix A. We can linearly interpolate between the extreme values of ga,air(Xw), yielding an expression for the size factor of air in a water background for the complete range of w: ga,air(Xv)XV - (ga,air (X) - ga,air,v(0)) Xf d < Xw < Xv Xv ga,air(Xw) = ga,air,d(0) + -(ga,air(xfld) - ga,air,d(O)) 0 < Xu < fld (3-31) Xv where ga,air(xv) and ga,azr,v(0) are solutions of (3.29) with moist air and ga,air,d(O) is the solution of (3.29) with dry air. For ice particles in a water background, Snowflow assumes that the shape factors are all 1/3 and uses this value in (3.22) to calculate ki,. Then the soil conductivity

23 as a function of temperature and total soil moisture content is: XA + xlk1Al + X2k2A2 + XairkairAair + Xik\X_ - - -- -- - -- _ - --- - --- - - - - -,X Wj ^_~ X ads Xu + x1kl + X2k2 + Xairkair + xiki K5soil(Tx~)= = KIdry - -W (dry - Ksoil(T. Xads)) X < Xads332) X ads where Xads is the volume fraction of adsorbed (bound) water and xi and Ai are the volume fraction and intrinsic conductivity of ice. Below this value, Ksoil is interpolated down to the empirically determined dry soil conductivity given by (3.24). Dry soil particles are of either mineral or organic composition. Of the mineral components, quartz content is most critical because of its high intrinsic conductivity (8.16 W/mK). de Vries characterizes soil as a combination of quartz and a composite of other minerals with a mean conductivity of 2.93 W/mK. The conductivity of organic matter is an order of magnitude smaller but highly variable-de Vries gives a value of 0.25 W/mK. Because scant empirical evidence is available on the matrix conductivity of organic soils, Snowflow simply assumes a quartz plus mineral mixture. The bulk (dry) density of this mixture is given by: pb = (1 - Xv)(WqPq + WmPm) (3.33) where wq is the weight fraction of quartz, pq is the intrinsic density of quartz, Wm is the weight fraction of other matter in the dry soil, and Pm is the mean intrinsic density of that matter. A typical bulk density of 1520 km/m3 was chosen, with the corresponding weight fractions found from (3.33) using xv from the REBEIX-1 test site. These parameters are listed in Appendix A.

24 3.4 Snowpack formation and heat and moisture fluxes Snowflow handles the snowpack processes of heat transfer, compaction, liquid water seepage, grain growth, and interaction with the atmosphere in separate but interdependent modules. Following discussion of simulated snowpack formation. this section examines each of these processes individually. Snowflow evaluates each precipitation event and flags it as snow if the wet bulb temperature of the air, T7b, is below freezing. Snowflow determines the characteristics of a new snow layer based solely on its initial temperature, Ts, which is set to T7b. The initial density of new snow is (modified from Anderson [37]): 75 + 1.7(Ts - 258.16)1.5 if Ts > 258.16 K Ps =75 if Ts < 258.16 K. (334) The initial grain size (a diameter in m) is a function of density [38]: gs = 0.16E-3 + l.lE-13p, (3.35) and the thickness of the new layer is given by: P ds = - (3.36) Ps where P is the precipitation in kg/m2. Other parameters discussed in the following sections are initially set to 0, including the liquid water content. Snowflow adds snowfall in consecutive time steps to the same model layer unless the layer has exceeded a maximum thickness. When additional snow is added, the characteristics for the layer given above are recalculated as a weighted average of the new snow values and the values for the existing layer. When a layer exceeds the maximum thickness, it is divided into a layer of optimal thickness dsopt and a layer of thickness ds - dsopt.

25 3.4.1 Heat and vapor transfer within the snowpack The Snowflow simulation treats the problem of heat transfer in the snowpack in a manner similar to that for soil. In addition to conductive heat transfer through snow grains, energy is transferred through shortwave radiation, water vapor, and longwave radiation. In practice, the transfer of heat by conduction through snow grains, conduction through the air in voids, and longwave radiation cannot be experimentally separated. Consequently, these mechanisms are combined into a single flux characterized by an effective conductivity, Keff in a manner similar to (3.10) yielding the one-dimensional heat transfer equation: HVs (T) IKeff,OT T+ t = z az +Keff az2 dHSVK is the differential enthalpy change per unit volume of snow due to conduction alone. Equation (3.37) only applies to dry snow. If the snow is wet its temperature must be uniformly To and dHVK = 0. Consequently, we can substitute dHVK = cipsdTsK in (3.37), where ci is the heat capacity of ice and dT,^- is the differential snow temperature change due to conduction alone. The substitution yields: 9T,_K Keff &T a2T ipsat a9z =z + Keff 2 (3.38) The effective conductivity as a function of snow density is given by [38]: IKef = Xair,d + (7.75E-5 + 1.105E-6ps)(A - air,d). (3.39) Since Kieff is a function of Ps only, it will be convenient to apply the chain rule to its derivatives, which yields: aT, 9KRff ap, 9T 92T cips = psff + Keff z2(3.40) is at -ap, az az az2 Snowflow uses Anderson's derivation of heat transfer by the processes of vapor diffusion and vapor transfer by temperature gradients [37]. Water vapor diffusion

26 is a function of the vapor density, f(T) = vs,-which is at its saturated value in a snowpack-and the snowpack effective diffusion coefficient, De(T), given in Appendix A. Defining 9f/aT = f', the snow temperature change due to vapor diffusion, dTs,D, is given by: Ts,D 02f ODe Of Of cipS = lsDe +l z z t (3.41) 2 f O(TnT)Of 9f = ICoD~TTM + lICoD -~ l z at 2 z z t = lsCDeT D T 2_ + lsCDe sD enDD1 -I = sCDeTnD f a2T + ( a + IsCDenDT dT\2 fI Of fz2 Kaz j \9z) at Combining (3.40) and (3.42) and noting that lsf' << Cips. we have the complete equation for internal heat transfer in the snowpack: _Ts_, aT,iK aTs,D cip st = CpS s + CiPs, (3.42) \2T 4[Ker ap] aT = [Kf + sCDef'T ]2 + [ I - z az2 [ap, z az +CDels [fInDT- + T D f"] (OT)2 In addition to internal heat transfer there are two sources of external heat: (a) the heat from shortwave radiation penetrating the snowpack and (b) the heat absorbed/released by melting/freezing. The enthalpy and temperature change due to these external sources is: aTsE aFs -- p aw (3.43) CZ Sat az w at where Fs,, is the flux of shortwave radiation at depth z in the snowpack. Combining (3.43) and (3.43) and converting to discrete time, we have the finite difference

27 approximation for the equalization function in snowpack layer i: E,i = AHi - Qs,i = dspsci(Tt+t - Tt) + If p(Wt+ - Wt) - 0.5At(Fti + Ft+,) -O.5.Atds [(Keff + IsCDefTD j I 2 [a fap 1 &T az2 ap, az z + +Des l +T dz2 ap,s 9z z + CDel.sflnDTlnD-l + TnDf] aT)t+t] -z (3.44) Snowflow uses Es,i in a manner comparable to Egi, (3.12). 3.4.2 Snow-soil interactions Snowflow parameterizes the heat exchange between the snowpack and soil using the steady-state heat flow equation: FlsoFlgo(Tso - Tgo) FloFlSO(Ts -T t-f g Q O.5At I( +0F0(T0-T0) V ( FlsoFlgo(Tso - Tgo) + QsFlodso/2 + Flod, /2 + Flgodo/2 + Flod o/2/ (34 where dso/2 and dgo/2 are the midpoints of the bottom snow and top soil layers, respectively, and F1 is a composite conductivity factor given in Appendix A. The steady-state approximation is applicable here because the insulating effect of the snowpack is likely to minimize the temperature gradient and, consequently, temperature fluctuations at the bottom of the snowpack. Snowflow also uses the steady-state formulation for heat flow between the top and second soil layers and between the bottom and next-to-bottom snow layers. The intermediate layer formulations are not usable for these two interfaces because it is not possible to calculate second derivatives there using (3.14) if either (3.12) or (3.44) are used as the equalization functions.

28 3.4.2.1 Handling an insulating layer of grass Snowflow treats the grass layer between snow and soil as a massless, thermally resistive thin layer which effectively reduces the snow-soil conductivity. We can generalize (3.45) as a problem of heat flow between two media with uniform temperature through a virtual slab: Qs,g= A(o - T0) (3.46) Rsg where Rs,g = L~/IK is the thermal resistance of the virtual slab, IK is its conductivity, and L, is its thickness. By equating (3.46) and (3.45) and assuming that do = dgo = 8mm, then for typical values of Fls and Flg the thermal resistance of the virtual slab is R,,g = 0.044 K/W. Now we can insert the grass layer as a parallel slab with thermal resistance Rgrass such that: Qs=,grass At (Tso T9o) (3.47) h s7,g + Rgrass A rough approximation of the conductivity of the grass layer is as a grass-air mixture with 1% grass by volume. Taking conductivities of Xorganc = 0.25 W/mK and Aair = 0.025 W/mK, the conductivity of the mixture is 0.027 W/mK. Then a grass layer about 3 cm thick has a thermal resistance of 1.06 K/W, or about 24 times greater than the typical value for Rs,g. Snowflow implements this approximate correction as: (Tsn - QtI) - s' Qs,g grass = At (T- T) - Q,g (3.48) 25R,9 25 3.4.3 Shortwave attenuation and absorption The shortwave heat source in (3.44), Fs,, is calculated in Snowflow by [37]: F, = Q (e-( - e- V(Z Z,-)) (3.49) SW(eT. z$

29 where z'r is the height of the top of the snowpack and z' is the height of the top of layer i. v is a single band attenuation coefficient given by: / \ 1/2Ps (z =o ) 1 P (3.50) \gs) pi PI where vo is an experimentally determined constant, vi is the attenuation coefficient for solid ice, and pi is the intrinsic density of ice. Further research since the publication of Anderson's thesis has shown that the attenuation rate in the snowpack is highly frequency dependent, with near-infrared wavelengths absorbed within the top 10 cm of the surface and shorter wavelengths penetrating to great depths. Brandt and Warren [45] suggest that a 10-band model would adequately describe shortwave snowpack attenuation. Any approach for attenuation in pure snow will require some modification for the relatively shallow snowpacks and short grasses of the Northern Great Plains. Since research concerning these conditions is unavailable, Snowflow uses the simpler, single band model as a starting point. 3.4.4 Snow-atmosphere interactions Fluxes between the Snowflow snowpack surface and the atmosphere constrain the Snowflow snowpack solution through the energy balance equation for the top snowpack layer (3.5). The atmosphere provides radiative fluxes from the sun and sky that add heat independent of snowpack temperature. The remaining fluxes depend on the surface temperature and the flux values are determine as part of the solution. The following sections describe Snowflow's parameterizations of atmospheric fluxes.

30 3.4.4.1 Radiative energy exchange Radiative energy exchange between the snowflow snowpack and the atmosphere has three components: absorbed longwave, absorbed shortwave, and emitted longwave radiation. In Snowflow, the absorbed shortwave heat per unit area over a time interval At is: Q ( 1-A) Q (3.51) where A is the albedo (shortwave reflectivity) of the snowpack and Qsw,inc is the incident shortwave radiation from the atmosphere (approximately wavelengths 0.3-3 Jum). Similarly, the absorbed longwave heat is: Qlw = (1 - el)Qwinc (3.52) where ei, is the longwave emissivity of the snowpack and Qluwinc is the incident longwave radiation from the atmosphere (3-100 tm wavelength). Longwave radiation emitted by the snowpack goes as the fourth power of the top snowpack layer temperature: Qeit -= t el a T (3.53) where a is the Stefan-Boltzmann constant. Qlw,inc is not available from REBEX-1, so Snowflow estimates it from air temperature and relative humidity [39]: Qlwin = AtLrT4ir (0.61 + 4.33E-3P(Tir)R) (3.54) Estimating Ql, inc by this method neglects the contribution of clouds to the longwave flux. (Clouds always increase the flux.) Although it is impossible to address the question of clouds precisely, Snowflow avoids very low values of Qlwinc by setting an artificial lower limit of 180 W/m2.

31 3.4.4.2 Latent and sensible heat exchange The only near-surface atmospheric conditions available from REBEX-1 are air temperature and humidity at a single reference height above the ground and wind speed at a second reference height. This strictly limits the methods available to estimate the fluxes of sensible heat and moisture, many of which are based on detailed field measurements of wind, temperature, and/or humidity profiles. Snowflow uses the energy balance (or Bowen ratio) method for calculating these fluxes [46] 47]. Over a time step At, sensible heat transferred, Qair, and potential latent heat transferred, Qpot are given by: Q = At ( cppairk - (3.55) [ln (zd) - /m3] [Iln (zd) - + ] Q =At ~ Ilkpairk2U3(q2 q- ) (3.q6) [In z3) - 3] [In (z) - q2 + 4q] where ka is the von Karman constant (0.4), cp is the specific heat of air, Pa is the air density, u3 is the wind speed at height Z3 above the surface, T2 is temperature at height Z2, T1 is temperature at height z1, d is the so-called zero-displacement height, zo is the roughness length of the surface, and q2 and qi are specific humidities. 1k is either lV if the surface is at or above the freezing point and Is otherwise. The II- where x may be m for momentum, q for moisture, or h for temperature-are functions of the atmospheric stability and are discussed further below. The energy balance method says little about interaction of the surface with the atmosphere except through the roughness length parameter. Snowflow uses the energy balance method heat transfer equations with modifications for estimating the fluxes from Snowflow's surface temperature, Ts. First, the height z1 is set at the top of the snowpack and T\ = T,. Also, the specific humidity at the surface of the snowpack is assumed to be at the saturated value for Ts so Qlatent = Qpot

32 The shear functions, 6z, and profile functions,, are found using an iterative method [47]. Initially assuming neutral stability (pm = 1,' m = 0), the Richardson number, Ri is calculated as: n ( Pz- ) 2 Ri = Rb Zo m (3.57) Rb= g22(dT+ ) (3.58) Tu2 ( dz ) where z = (z1 + z2)/2, Rb is the bulk Richardson number, u is the wind speed at z, T is the temperature at z, Yd is the dry adiabatic lapse rate (9.8E-3 K/m), and g is the acceleration of gravity. Stability is a function of Ri: Ri > 0 implies stability, in which case: z Ri L 1 - 5Ri' m- =1 + 5, L 2 m = -5 -Z3 Vm3 - -5 —. Oq2 = -5Z2 gl = -5 ZL' nditions: L = Ri, x(Z)' - = -ln((1 + x(2)2)(1 + x(z))2 8)2 ) For stable co (3.59) (3.60) (3.61) (3.62) (3.63) (3.64) (3.65) (3.66) (3.67) (3.68) (3.69) (3.70) I7 - 2 tan-1 x(Z) + - 2 tan-1 x(z3) + 2 2 =q3 = 21n( 1 + - 16 ql —=21n [ + 1-16 -

33 where / z 0.25 x(z) = (1 - 16L. (3.71) 3.4.4.3 Heat from precipitation When Snowflow determines that precipitation falling on the snowpack is unfrozen (see the beginning of section 3.4 for this criterion), the water temperature is artificially set to To. The water then interacts with the snowpack through the seepage parameterization in section 3.4.5. In nature, when rain water comes in contact with a snowpack, either some portion of precipitation is frozen or some portion of the ice in the snowpack is melted during the event. Snowflow's approach omits the ability of rainwater to melt ice in the snowpack. The omission does not affect results in this thesis because there was little rainfall on the snowpack during the test period described in Chapter 7. 3.4.5 Liquid water seepage in the snowpack Liquid water seeps through the simulated snowpack under the pull of gravity when the liquid water content of a layer exceeds its holding capacity. Snowflow's parameterization of seepage comes from Jordan [38] with some modifications. In the numerical simulation, the solution of the heat transport problem (see section 3.2) precedes and is independent of the seepage solution. When the heat transport problem is solved, the initial effective saturation, se, of each layer is given by:2 St Sr e 1 - Sr (3.72) 12Superscripts in this section indicate the snowpack layer. 2 Superscripts in this section indicate the snowpack layer.

34 i=Ns l a y e r s * - - - - -.................................................................. i= 0 ''''"':'=0'" '""""""' '' Figure 3.2: Schematic of liquid water seepage through the simulated snowpack (after [38]). and Si= 7l (3.73) pwci ~~= 1- I (3.74) Pi i = - (3.75) - WVPw (3.76) where s is saturation, Sr is the irreducible water saturation constant, 4 is porosity, 7Yi is bulk ice density, and 7/ is the bulk density of liquid water. If the si of a set of adjacent layers is greater than zero, then there is downward flow in the region. Figure 3.2 shows the geometry of a snowpack with two active flow regions. Snowflow finds the post-flow effective saturation of each layer in a flow zone by solving a set of linear equations of the form: As'1+ A+A's = (7 (3.77)

35 where the coefficients are, for the top (nth) layer in the flow zone, A' = 0, (3.78) A' = c7i + 0.5Akims, (3.79) B" = -O.5(U+1 - Uli) + c7is, (3.80) and for the interior layers of the flow zone, A'= -0.5Aki+lm 1, (3.81) Ai' = c7i + 0.5Akmi, (3.82) B' = -0.5(U+1 - U[) + c7is + 0.5Aki+lb+ - 0.50.5Akcib. (3.83) The parameters Ak and c7 can be calculated independent of position of the layer in the flow zone: Aki= pg,max (3.84) c7' = - ) (3.85) where I:M a = 0.077gs2e- 0078. (3.86) pi is the dynamic viscosity of water (1.78E-3 Ns/m2) and gs is the average snow grain size (diameter). Away from the edges of the flow region (i' # 0, i' $ n'), the constants b, and ms are given by: m' = 3(s,)2 (3.87) b = -2(s)3. (3.88) At the flow boundaries, bs = 0 and ms is found by first estimating the effective saturation, Se,est, by solving the cubic equation: 3) c7 c7zs - + 0.5U 3( + 0.5Ak- --,e 0 (3.89) 05k 0.5Aki

36 Then m, = se,est The liquid water flow rates, UI, is the final product of the flow calculations. Its value in the equations above is estimated from its value in the previous time step. After (3.77) is solved for the effective saturation of the layers in the flow zone, then Ul is given by: U = -Aki(s)i. (3.90) U/ is always negative and represents the flow of water out of the bottom of the snowpack layer. Water that flows out of the flow zone into a dry layer will freeze until an equilibrium is reached. Appendix A describes this calculation. 3.4.6 Compaction and mass balance of snowpack layers The compaction rate of a snowpack layer, CR, is a function of the structure of the layer and the overburden of snow above [38]: 1 ad, d, dsat C — = 2.778E-6c3c4eo0-04(To-T) + e-5(T-T)-6 (3.91) 77o where Ps is the overburden pressure, r/o is the viscosity coefficient (at T = To and p, = 0), C5 and ce are empirical constants, and c3 and c4 are given by: C3 = C4=1 if 7 = 0 and 7y < 150 kg/m3 C3 = e-0.046(yi-150) if 7y < 150 kg/m3 (3.92) C4 = 2 if yT > 0. Water and vapor movement in the snowpack change the mass of the snowpack layers. Mass changes coupled with compression alter the layer's density. Section 3.4.5 provides the solution for fluid the fluid flux, U1. The net fluid flow for an interior layer

37 i is Unet = U - U/. For the top layer (i = Ns) the net flow is: ~' -lP U - U1N. (3.93) At Note that positive Ui,net means net flow out of the layer. The vapor flux from snowpack layer i to layer i - 1 below is given by: 2CD T nr -n Ip/ — l " - 1 2CeTL ifT L(Ti - T1) (3.94) dsf[-_ + ds,i-lfiI The net flux for interior layers is Uet = Ui-1 + Ui+1 where Ui1- = -Ul-i. The net flux of the top layer is: r _ Q- atent if TN' > To, At I v,net = (3.95) uNf - Q-atent otherwise. V At is Then the change in the total mass of the layer is given by: A(psd)i = (sds)t+ - (psds) = -At(Ul,net + U,net). (3.96) The compaction rate then gives the depth of the new layer: (psdst+At dt, (1- At CRN ) N,) for the top layer, = RN ds)t+At + UN, (Psa t _psd sN, v,net (3.97) dt (1 - At CRi) otherwise. And the new snow density is: (pd,\t+At UNz ( ) + Utnet for the top layer, $sN, Ps,i = (3.98) pd+ otherwise. n n s,2 With this parameterization, the net vapor flux of the top layer changes its thickness and not its density.

38 3.4.7 Metamorphism of snow grains The size of grains is important in the determination of the structural stability of a snowpack for avalanche forecasting and the support of weight. In the energy and moisture transfer equations, grain size affects the permeability of the snowpack to fluid flow (section 3.4.5) and the shortwave extinction coefficient (section 3.4.4.1). The grain size is also of critical importance in calculating the brightness of microwave radiation emitted by the snowpack, discussed further in Chapter 4. Snowflow uses (3.35) from Anderson [37] to initialize snow grain size and follows the method of Jordan [38] for modeling growth of grains. Snow grain growth is an area of continued research effort and these semi-empirical parameterizations are interim solutions. Snowflow calculates the change in grain size, gst+ t- gst, in dry snow as a function of the water vapor flux, Uy: gs = gs + At gl (3.99) gst where Cgsl is an empirical constant and Uv = CDeT a f I (3.100) In wet snow, grain growth is accelerated as smaller grains-whose equilibrium fusion temperature is smaller than that of larger grains-melt and feed the growth of larger grains: gst + At 2(xI + 0.05) 0 < xi < 0.09 gst gst+At= (3.101) gst + At g2 (0.14) xI > 0.09 gs where Cgs2 is an empirical constant and xi is the volume fraction of liquid water.

39 3.5 Discussion Chapter 7 describes the results of Snowflow's snowpack simulation when it is driven by the REBEX- 1 atmosphere. That chapter discusses the short-comings and sensitivities of the snow and soil models and the way that REBEX-1 data are linked to the Snowflow parameterizations described in this chapter.

CHAPTER 4 MODELING RADIOBRIGHTNESS OF THE SIMULATED SNOWPACK 4.1 Introduction This chapter describes Esnow, a snowpack emission model designed to take full advantage of the detailed data available from Snowflow's snowpack simulation described in Chapter 3. The Snowflow model provides a snowpack stratified in up to 50 layers with temperature, density, grain size, and wetness information in each layer. Beneath the snowpack is a soil half-space with temperature and moisture information. The objective of the Esnow model is to estimate the brightness of microwave emission and the reflection characteristics of the terrain over a broad range of snowpack conditions with reasonable computation speed. 4.2 Conventional radiative transfer theory Conventional radiative transfer (CRT) theory is based on a heuristic development of the problem of intensity transmission through attenuating media. In a medium containing scatterers. CRT assumes that each particle scatters independently of the others, that is, as if it were the only scatterer. In fact. when scatterer density exceeds about 1% by volume, the scattered fields interact in phase, resulting in reduced scattering by the individual particles. 40

41 19 I(r;s) + dl dA /0 I(r;g) Figure 4.1: Conventional radiative transfer. This section develops the form of the differential equation of transfer that Esnow uses to solve the snowpack emission problem. Section 4.3.1.1 will account for the inadequacies of the independent scattering assumption by introducing an empirical correction to the independent scattering extinction coefficient that reduces its magnitude to the level of experimental observations and more sophisticated models. Figure 4.1 diagrams the process of extinction by scattering and absorption experienced by a spectral intensity, I, as it traverses an incremental length, dr, along the ray r with a directional cosine Us = cos 0s. A plane-parallel geometry is assumed with symmetry about the z-axis. The change in intensity, dI, is given by: dI = ( -K,Idr - KIdr + KJadr + KJdr (4.1) dlh where I, and Ih are the vertically and horizontally polarized components of I, ta and Kf are the absorption and scattering coefficients (per unit length) in the medium, and Ja and Js are the absorption and scattering source functions per unit length. The differential optical depth, dT, is defined as: dr = tedr (4.2) where Ke =,a + Es is the extinction coefficient of the medium. Substituting (4.2)

42 into (4.1): d =-(1-a)I -aI + (1 -a)J, + aJs =-I+(1 -a)J, + aJ, (4.3) where a = Is/Ke is the single scattering albedo and, for independent scattering, is equal to the ratio of scattering cross section to extinction cross section for a single particle. The absorption source function can be deduced from Kirchhoff's Law: Under conditions of local thermodynamic equilibrium emission must be equal to absorption. The radiation from a differential element of the transmission medium balances that of its neighbors isotropically with a spectral intensity per polarization given by Plank's black-body radiation law: Ibb,p = c2 (ehf/ - 1(44) where p is h- or v-polarization, h is Plank's constant, f is frequency, k is Boltzmann's constant, T is temperature, c is the speed of light in vacuum, and n is the index of refraction of the medium. At microwave frequencies, the exponent is small and the Rayleigh-Jeans approximation can be used to linearize the formula: bbp 2 (1 l+hf/kT+.)-I) 2 (4.5) At 85.5 GHz-the highest frequency of interest here-the error in using the RayleighJeans approximation is +0.7% at T = 300 K. The absorption source function is then a function of absolute temperature only: Jap kf2n2T(r) (4.6) By convention, the correspondence between temperature and the source function leads to the definition of a convenient quantity, brightness temperature. The brightness

43 temperature of a source in terms of its radiant intensity measured in vacuum is: 2 Tp = I kf(4.7) Along a ray traveling through media of varying refractive index, a more appropriate quantity is the index-normalized brightness temperature which remains invariant along a ray in lossless, homogeneous media: _I c2 T -,p - k (4.8) z2 kf2' Since the problem considered here involves radiative transfer through media with continuously varying refractive index, index-normalized brightness temperature is used as the propagating term. The absorption source function for index-normalized brightness temperature is: Ja,p c2 Tn,a,p - = T(r). (4.9) n2 kf2 The scattering source function accounts for radiation from all incident directions, (/i,< i), being scattered into (/i,,,): J,(Ps) = ~ Jf Ps(/s, s; Hi fi)I(pi)dni. (4.10) where 6 is the azimuth angle around the z-axis and dQi = sin OidOid-di is the differential solid angle. Ps is the so-called phase matrix and its elements determine the scattering transformation from polarizations pi and direction (p*,co) into ps and (Ps qO). Following from (4.9), define the equivalent scattering index-normalized brightness temperature source by: Tn 2 kf2 JJ7 Ps( ~; )T()d. (4.11) Inserting (4.6) and (4.10) into the differential transfer equation, (4.3): dI + 2n2 (T(r) ) fp(; ~)I()d. -I+(1-a)k ( 2 P+4~ ff Pts1,) (4.12)

44 z=O 1 )T4 z = -d Figure 4.2: Upwelling and downwelling brightness in a snowpack layer. Replacing intensity with the definition of index-normalized brightness temperature using (4.8) and eliminating common terms, we have the differential equation of radiative transfer in terms of T,: dTn dT = -Tn + (1 -a)T,a + aTn,, (4.13) dr where Tn,a and Tn,s are given by (4.9) and (4.11), respectively. 4.3 Solution for the simulated snowpack The plane parallel Snowflow snowpack is naturally stratified by the spacing of snowfall events, with some artificial limits on layer thickness as discussed in Chapter 3. Each layer has separate values for temperature, density, grain size, and wetness. Consider the upwelling brightness, T+, and downwelling brightness, Tn-, through a layer as in Figure 4.2. From (4.13), the transfer equations are: dTd+ (i+) dTr = -TZ(y) + (1 - a)Tn,, + aTn,(+, ) (4.14) dTP

45 and dTd 8 ) = -Tn () + (1 - a)Tn,a + aTn,s(P, s ) (4.15) where drT = ~i-dz. (4.16) Tn,s alone has an azimuthal (5s) dependence through Ps in (4.11). Integrating (4.14) and (4.15) over O5: dT+ (/) = -+( ( a)T dr+* (4.17) + 9 |10 P5(+;i)T+(ui)d/i + 2 j Ps~(1+;-/i)T~ (3t)dui 2 sO (I's; Y' n T 2 so(PJO P/ n 7 and dT- (/~7) - = -T7() + (1 - a)Tn (4.18) d7 -+ 2 f Po(-si; U)T+(/i)dpi + Po(-~;- -j)T-(p)dp 2 JO 2 JO where Ps(s;pi)= 2 os dPs(s s, (4.19) We will integrate (4.19) when the functional form of PS is given in Section 4.3.1. To solve the differential equations (4.18) and (4.19), we divide each equation by attenuation from the appropriate boundary: _ _d[T+Q$)eTT+-d+z) dT.r + T(+ T ) (udz) en (-(-d,) +T )e ()] Fe+(,) (4. 2 (4.20) +1 eT(oz) d[T-(p)e7iOz)] = F) e(O dr- d7-21) 1(4.21)

46 where F+ is the sum of source terms F~(z,/u+) = (1-a)Tn,a (4.22) + a j PS(~/4; i)T+ (i)dy + - Po(~-i; — i)Tn(t )dpi 2 Jo z Jo 2 so S I Z and T (Z1, Z2)= j dTr. (4.23) z1 Integrating (4.20) over (-d, z) and (4.21) over (0, z) and substituting (4.16) for dTr, we have: T+ (/)e'+ (-,) _ -T (-d, +)e ( —d) = F+ - (-d,z') dz' P (4.24) T-(z -)eT'(OZ) -T(0 -O )eT() JF= - jF-(Oez') (4.25). J (4.25) Solving for T (z,/ ~ ): T+(z, +) = T+(-d,e ( ) + F+e e )$ dz' (4.26) TT (z,; ) = T (0, p-)e-(~oz) + j Fj-e-(z',z) dz'. (4.27) Equations (4.26) and (4.27) were derived for a Snowflow snowpack layer of thickness d. We can divide the Snowflow snowpack into layers of arbitrarily small thicknesses, 3, using interpolation from the original snowpack layer mid-points to calculate temperature, grain size, and densities. Then by choosing 3 small enough to neglect the depth dependence of ne and F~, we can make the approximations: ez'T + +z T(, ) (4.28(-) )e-'+(-,z) + +-'ldz P (4.28)

47 T (z, -) T (0, - )e — (O,) + F -e _ez/ e Ps Jz (4.29) We then solve for the upwelling brightness at the top of the layer, Tj+(O,P+), and downwelling brightness at the bottom, T- (-6, y-)): T+(0, +) T+(-, +)e —+(-6o) + F + (1-e-e6/S +) (4.30) $ n4 + F T5 (-5~;') T T (O, )ee-"-(~'-6) F- Fl e-s6/iLs( - e~es//s). bs Ke (4.31) Combining terms: T+(0,T+ ) T+ (-+, +)e-+-5,) + F+(1 - e-_e/t) (4.32) T-(-,;-) T T (0,/ )e-7 (~0-) + F (1 - e-e6/^7). (4.33) Recall that the source terms (4.23) include scattering integrals over the incident directions, fi. The method of Gauss-Legendre integration (or Gauss quadrature) allows us to numerically evaluate these integrals using a small number of Ptj and a set of N weighting coefficients, wj [48]. Gauss quadrature is accurate when the integrand is well-approximated by a polynomial and is frequently used with radiative transfer problems [49]. Applied to (4.32) and (4.33) we have: T+(0, +) m Tn-(-+( —6,o+)e-(- (1 -e-6/s )[(1 -a)T, (4.34) a N + w (Pso(+;,uij)T+ (^j ) + Po (;- )T ( ))] j=1 T (-,;-) T (0,/p-)e- (-) + (1 -e-e/ )[(1 - a)Tn, (4.35) + E (Pso(-/ -; / u)T+ (pi) + Pso (-/-; -'ii)T- 7(/ ))] 3=1 where the scattering directions, ps,: must now include the discretized incident directions, utij, as well as the direction corresponding to the observation angle of interest,

48 lobs. In this work, N = 5 was sufficient in that greater N did not significantly change the radiative transfer result. With (4.35) and (4.36) for thin layers in hand, we can use the method of invariant imbedding to find a solution for the whole snowpack as a stack of thin layers (Section 4.3.2). But first the next section deals with the form of the scattering phase matrix. 4.3.1 The scattering phase matrix and Mie scattering Grains in the Snowflow snowpack are modeled by a single grain size-the diameter, gs in each Snowflow snowpack layer. Although the grains in true snow are both indistinct and not spherical in shape, they are usually randomly oriented. Assuming this is the case, the Mie scattering phase function for spherical particles is an appropriate approximation. If the particles were small compared to the wavelengths of interest, then the simpler Rayleigh scattering function could be used. According to [50], the validity criterion for Rayleigh scattering is InrXl < 0.5 where x= 2rr (4.36) Ao6b is the size parameter, r is the particle radius, AO is the free space wavelength, eg is the background dielectric constant (the real part of the complex permittivity), and nr is the relative index of refraction of the scatterer (1.77 for ice). Table 4.1 compares Mie and Rayleigh scattering at two grain diameters-0.17 mm and 2.2 mm, the largest value from Snowflow. Although for the smaller grain diameter the Rayleigh criterion is satisfied at all three frequencies it is not valid even for the lowest frequency at the largest grain size. Here, the error in the Rayleigh scattering cross section, Qsca is 4.6% at 19 GHz and 117% at 85 GHz. Derivations of Mie scattering can be found in many texts (for example, [51]) and

49 Table 4.1: Mie-Rayleigh scattering parameter comparison. Qsa is in units of m2 will not be repeated here. Figure 4.3 diagrams the primed particle coordinate system showing the nominal scattering plane and the perpendicular and parallel polarization directions of the scattered radiation. If the incident fields are polarized vi and h as shown. the relationship between incident and scattered field components is given by: (Ej )p -s eik(r-z')SI Sh(E (437) E_ -ikbgr S~ S lh Eh where kbg is the propagation constant in the background medium and the Spq are the scattering amplitudes relating incident to scattered field at incident polarization v or h and scattered polarization 11 or 1. The scattering amplitude matrix in the particle reference frame is given by: Sp z rtc Sv Sllh S = ( 52(i5) cos 1/ -S2(s) sin' (4.38 \p Srvt Slclh - V Si (t \) sin ' Si (\u ) cos (4 ) where S2 and Si are the co-polarized and cross-polarized Mie amplitude scattering matrix elements in the plane of incident polarization and are functions of Mu only. Assuming incoherent fields, we can then write the scattering phase function for a particle in the particle's reference frame: 47, 47r- [ Sii]z2 [Sllh]2 Ps particle(:< s; pi ) = (Msca) = 4 ( S 2 [Sfh |) Qsca Qsc akg i l12v ISh2 / (4.39) where (MIsa) is the Stokes matrix for the scatterer. Figure 4.4 diagrams the transformation of the scattering amplitudes from the particle scattering plane polarizations, 11 and L, to h and v polarizations relative

50 z! Scattered radiation f/i ells Particle I Nj iss x' o Incident radiation Figure 4.3: Scattering geometry for an isolated particle.

51 z A v Scattered As radiation /, X / / "/ X, Incident radiation / z' I I ells,' \ ' hI In y — t \/ \ / / \ / \ / \/ Figure 4.4: Particle scattering geometry in the snowpack coordinate system. to the snowpack horizon. In the snowpack coordinate system, the z-axis is vertical and the x- and y-axes are arbitrarily oriented. The angle A is the azimuthal angle between the incident and scattering directions. The snowpack azimuth origin (93 = 0) is arbitrary so for convenience assume A = Oi - c5. The x'-axis from Figure 4.3 and the vi polarization are always in the same plane as the snowpack z-axis. The scattering amplitude matrix relating incident to scattered fields polarized with respect to the snowpack horizon is: Shy S Cos sin - ( snowpack = Shv ( c I' /, Saricle ( -hv Shh — s osA - - /

52 where, in terms of the snowpack scattering coordinates, sin " = sin 0i (4.41) sin Os cos d" = - cos A cos ' + sin A sin 4' cos (4.42) sin 80 sin A sin' = sin s (4.43) sin 0' cos 0' sin Oi - sin 0s cos A (444) cos' = (4.44) sin 0' cos 0i cos 0' = cos Oi cos 0s + sin 0i sin 05 cos A (4.45) A = s -X i. (4.46) The scattering phase function for a particle in the snowpack reference frame is: 47 = Ps(s, s; [Hi, Oi) = QsakSsno-wpack. (4.47) Qsca y9 Recalling (4.19) for the scattering phase matrix for the plane-parallel geometry, PSo(,s;/ti) = - j d ds( i d,iPs( 4 [,s; Hiz, (4.48) we can now eliminate one integral by z-axis symmetry and replace the other with A: - s~1 2r - Pso(s;^i)= - 2 I dAPs(/; p; A). (4.49) The remaining integral can be solved numerically using Gauss-Legendre integration as described above. The solution converges with 16 azimuthal angles. 4.3.1.1 Attenuation coefficients: Empirical scattering reduction Laboratory experiments have demonstrated that the assumption of independent scattering does not hold for particle densities in excess of about 1% of the total volume [52]. Calculations applying QCA-CP-PY-the quasi-crystalline approximation with coherent potential and Percus-Yevick pair-distribution function-reflect this result and can be used with the laboratory observations as the basis for a heuristic correction to independent scattering [53]. Figure 4.5 shows scattering coefficients from

53 independent scattering and an empirically adjusted QCA-CP-PY fit as a function of the volume fraction of ice (particles) in air. The normalized independent scattering coefficient is given by: V scat Vscat Qsct= (T~TQj = f (4.50) where Ks,ind is the scattering coefficient for independent scattering, f is the volume fraction of scatterers, Vscat is the volume occupied by a single particle, lVcat is the number of particles per unit volume, and Qscat is the scattering cross section of the particle (from Mie theory). The curve in Figure 4.5(b) fits independent scattering to empirically adjusted QCA-CP-PY results: =,QCA 7(1 - f)(|0.5 - f13 + 0.015). (4.51) Ks,ind Although experimental observations are made of total extinction not the scattering coefficient, variation in extinction in dry snow should be controlled by the scattering term. QCA is only slightly dependent on the background moisture content of the snowpack and at the frequencies of interest here attenuation by the background medium quickly exceeds attenuation by the particles [54]. Consequently, the same scattering reduction formula is used for wet and dry snow. The QCA-CP-PY scattering coefficient fit in Figure 4.5 approximates results reported in [53] and [54] for volume fractions up to 0.4. The latter used a mean particle radius of 1 mm and showed that at f = 0.4 extinction coefficients were about 6.5% of those from independent scattering over frequencies from 35 to 140 GHz. Esnow infers the shape of the QCA-CP-PY fit for f > 0.4 by assuming it is symmetric about f = 0.5. As f approaches 1, scattering becomes a function of air pockets in a mostly ice matrix and is roughly analogous to scattering by ice grains in air.

54 Q).0..0 0 C L0) 4-~ 03 C.) -a () 0) N.Z E 0 z 1.0 -0.8 -0.6 -0.4 -0.2 - (a) - - - - -...........r.* ^' 7* ---- QCA-CP-PY ^* ^* 0.0 - IIi,, II - 0.0 I I I 0.2 0.4 0.6 0.8 1.0 Volume fraction of scatterers 1.0 0.8 0.6 0.4 0.2 0.0 (b) Scattering reduction factor QCA-CP-PY/lnd. scat. 0.0 0.2 0.4 0.6 0.8 Volume fraction of scatterers 1.0 Figure 4.5: (a) Normalized scattering coefficients and (b) the independent scattering reduction factor.

55 In summary, the Esnow attenuation coefficients are: Ka = NscatQabs + bg (4.52) Ks = NscatQscat7(1 -)(]0.5-f3 + 0.015) (4.53) Ke = K = (4.54) where Qabs is the absorption cross section of a particle (from Mie theory), Kbg is the background absorption coefficient given by 2kbg(1 - f) Kbg = g(- Im{b}1, (4.55) n2bg and nbg is the complex index of refraction of the air-water background medium given in Section 4.3.4. 4.3.2 The method of invariant imbedding Equations (4.35) and (4.36) gave the upwelling and downwelling brightness from a thin snowpack layer in a set of discretized directions. The goal of this section is to solve the snowpack emission and reflection problem by iteratively adding thin layers from the bottom up and calculating upwelling emission and half-space reflectivity concurrently. This method of solving the radiative transfer problem for a layered medium has been called invariant imbedding [55]. Figure 4.6 diagrams the boundary conditions of the kth thin layer of the snowpack for which T+ and Tn satisfy (4.35) and (4.36). Consider the upwelling brightness at the top of the layer to be the following sum: T (Zk,Ps) = ak(ts)T' (Zk_-1, ) + T+o(zk, s) + T+ (Zk, s), (4.56) where N 1 T47) T+'(zk-1,ps) = T(zk-l, s) + Rk-l(us;u)T o(Z,) (4.57) j=1

56 Air m z(k+1) Snowpack layer k z(k) / / 7 0 / 0 0 Tn(k) / 0 \ T^(k-1) 0 - 0 /0......................................................................................................................................... - ------- ---- - 0 Soil --—, = Attenuated Figure 4.6: Brightness balance for the kth thin snowpack layer. T+(zk-1,,us) is the upwelling brightness leaving the k - 1 layer and ak = e-/= is the total attenuation of the kth layer independent of polarization. T+o and T-o are upwelling and downwelling absorption source terms from (4.35) and (4.36). Since + = -/, these terms are equivalent: T+ = To = (1 - a)Tn,a(l - e-'/"s1) = T (Zk, s). (4.58) RkA-l(,s; [ij) is the reflectivity matrix at the k - 1 boundary for incident radiation at the discretized direction 'ij and reflected/scattered radiation at Us. In general, Rk- (Us; pij) includes cross-polarization terms due to scattering in the volume below the boundary. The last term of (4.56), T+(zk, u4S), is the scattering source term for the layer and is given schematically by Figure 4.7. Since the scattering phase matrix for the layer is symmetric with respect to the vertical direction in the snowpack, it is convenient

57 Tnk) 41 /Pt- //P Snowpack layer k / // Tl(k-1)..... = Attenuated Figure 4.7: Components of the layer scattering source term, T+,. to define: Po = PSo(i; ij)- PSo(-~S; -i) (4.59) Po = Pso(-,; fij) = Pso(s; -fj). (4.60) In the schematic, P.+ represents scattering away from the hemisphere of the incident radiation and P.o represents back-scattering into the hemisphere of the incident radiation. The scattering source for the layer is a function of upwelling radiation only. (Downwelling emission from the layer, T-o, is reflected and included in T+'.) Upwelling radiation at all angles contributes to T+ (zk ), so to simplify the notation we can redefine the brightness vector to include all of the discretized propagation directions, the principle observation direction, and both polarizations: T, (i) i T+, (v.1) T () nh (PN ) Tnh (obs) The corresponding scattering matrix can then include the Gauss-Legendre weighting

58 functions and other terms from (4.35): 7)soQ (4.62) p4- (pl; pl)13,wl vv P,(YN; Pl)NW1 P(Yobsl Yl )obsWl Ph (pi; P0131wi Ph(PN; POONWI Ph(Pobs; 1-tl)OobsWi P..QPI (P; /PN)/31bWN..P~([LN; UN)/3NWN.. Pj~4tobs; PLN)!3obsWN 0 0 0 0 0 0 P~jbLN; /P1)!NW1 PhhjIuNb; tl)fObNW1..P~h([LN;ItN) 3NWN..P~h(,obs; 1-N)!3obsWVN..Phh(IUN;AtN>3NWN..Phhjj(Mobs; 1-N)i3 obsWN 0 0 0 0 0 0 where!3j = a(l - a,)/2. Note that the zeros correspond to incident radiation in the observation direction which do not contribute to the Gauss-Legendre approximated

59 scattering integral. Similarly, define a multi-directional volume scattering matrix: Z= -(4.63) Rvv(li; 1') vv Rv(/1;/ N) 0 Rvv(,tN; ul) - (pN'N) 0 Rvv(/obs,; t1.)..*** Rvv(obs; N) R&v(/tobs, Lobs) Rhv( 1; '1) '"' Rhv(.1; lN) 0 Rhv(lN; 1).. Rh, (lN; AN) 0 Rhv(/l(obs;I1)." Rhv(Pobs;1 N ) 0 Rvh (P1; l1) **' Rvh(Pl1; N) 0 Rvh(N; U).' Rh(N; PN) 0 Rvh(/lobs; i'1).'' Rvh (Pobs; fN) 0 Rhh(1;'1 ).'* Rhh(Pl; MN) 0 Rhh (N; 1 ) ' Rhh (UN; pN) 0 Rhh(pobs; 'l)... Rhh(,Pobs'; uN) Rhh(Mobs: obs) Additionally, we will need a diagonal matrix, ', whose elements are a(us), the layer's attenuation in direction p,. With these definitions, if we apply (4.35) up to third order then the scattering source term in (4.56) is: Tn+ (zk) = (Po k + Pkkk-lPPso + kPsotkkP-1PR l + kP P + PsTZok-1PsoRk-1Pso + * )+'(Zk-1). (4.64) Similarly, the volume reflectivity matrix for the kth boundary is given by: Rk = Pso + ckRk-1'Pso + kRk-lk + * -. (4.65) Increasing the number of terms in (4.64) or (4.65) did not change the result of the radiative transfer solution for the snowpack significantly.

60 4.3.3 Air-snow and snow-soil boundaries The only boundaries of significant dielectric contrast in the Snowflow snowpack are the soil-snow and snow-air interfaces. Section 4.3.4 discusses the calculation of the dielectric constant for these media. At the soil-snow interface, Esnow uses a specular boundary. As discussed in Chapter 3, the Snowflow soil contains a significant amount of unfrozen water even when the temperature is below the freezing point. Consequently, the soil is lossy with shallow effective emission depths and it can be modeled as a dielectric half-space. The soil-snow reflectivity matrix is: Ro = FO | h2 (4.66) where |p 12 soil snow - snow Ysoil. ir0,v\ 1 sot (4-67) 2 soil l snow + hsnow[t soil Ohh 2 soil /1 soil - snow (4 snow.68) iisoilftsoil + flsnowfisnow The direction cosine in soil is fixed by Snell's law and the discretized propagation angles in the snowpack: lsoil = cos [sin-l (nsuins (4.69) nsoil Emission from the soil-the initial condition for the snowpack emission solution (4.56)-is given by: T+(zo,/ s) = Tsoil(I - Ro) (4.70) where I is the identity matrix and Tsoil is the soil surface temperature. The snow-air reflectivity matrix is similar, with v- and h-polarized components: - o 2 Jia 12 - |snowp/air -airP'snow (4.71) air(vv -4.72) Insnowslair +,air/rsnow 2 Fair hh 2 - snow/snow - airair (4.72) nsnow/Isnow + airiair

61 When the air-snow boundary is reached using invariant embedding, Esnow calculates the snowpack brightness iteratively in a manner similar to (4.64). First, downwelling sky brightness is transmitted through the air-snow boundary, reflected, and depolarized by the volume, contributing to the topmost upwelling vector: 7(Zn+l) = iRm7 (4.73) Then multiple reflections off the boundary and scattering by the volume are accounted for, the upwelling brightness is transmitted through the snow-air interface, and the sky brightness reflected off the surface is added: T7+(Zn+l:, Is) = (I + RmRair + RmRZairRmlazri +.. )T+ (Zn+l) + (I- air)Trnky. (4.74) 4.3.4 Dielectric properties of water, ice, snow, and soil The complex permittivity of free water, w,, is well known as a function of temperature and frequency in the microwave spectrum [50]. Ice is generally accepted to have dielectric constant independent of temperature and microwave frequency: E = 3.15. (4.75) The imaginary part of the complex permittivity of ice is small (especially relative to that of water) and. consequently, difficult to measure. Esnow uses the following empirical formula from [56]: ' + = 7Off (4.76)

62 where a = (50.4 + 62.00)1.0E-4e-22.10 0 = 300/T + 1 0.585E-4 = 1.0E-4(0.445 + 2.11E-3T) + (1 /291 (1 - T/29.1)2 T is temperature in ~C and f is frequency in GHz. The complex permittivity is i = + i<'. The complex permittivity of dry snow can be calculated from the empirical formulas [50]: e = (1+ 0.508E-37y)3 (4.77) Eds 3_ 3ds(-,, / 2e)( + c1) (4.78) where 7y is the density of ice in situ (the snowpack density) and pi is the intrinsic density of ice. The background permittivity of dry snow for Mie scattering calculations is 1. The complex permittivity of wet snow can be found by solving the Polder-Van Santen mixing formula for the water-air-ice matrix [50]: ws = ds + (w - s) 1+ A ---- (4 79) 3 =,b, A(4.79) where xy is the volume fraction of water in the snowpack. The coefficients are: Aa = Ab = 0.475 and Ac = 0.05 [57]. Equation (4.79) is cubic in e~s and can be solved using a numerical complex root finding routine. The wet snow background dielectric is a water-air matrix whose permittivity can be found by replacing the host permittivity, fds, with 1 in (4.79): 9xw c + A X ) ] ( bg = + (.- (4.80) 3Lbg

63 Esnow calculates the soil complex permittivity based on a dielectric mixing of free water, bound water, dry soil, and ice. The rotation of water molecules that are adsorbed (bound) to soil grains is inhibited by weak chemical attractions [39]. Esnow assumes that only water in excess of 7% by dry weight is free while the remainder is bound. Then the soil complex permittivity is given by: Esoil = Z7% + XfreeP7c(flEw + (1 - f)i)i)/Pw (4.81) where 67% = 3.3 + i.4 f Xu - Xbound if = Xfree f xw if Tsoil > Tfpd CWU P(To -Tsoil)iwupb/Pw otherwise Tfpd = To- xuPi zwu Pb Pb P7% 0.93 XwPw Pb + XwPw p7%(mw - 0.07) pr(1 - m.) Xf ree Pw( - m ) Xbound = Zw- Xfree. 67% is the complex permittivity of soil with 7% moisture by weight (from [50]); fi is the fraction of the free water (unbound) that is in liquid form; xu is the unfrozen soil moisture volume fraction; Tfpd is called the freezing point depression and is the temperature at which x = xu; P7% is the bulk density of soil with 7% moisture by weight; m, is the weight fraction of all water with respect to the moist soil; and Xfree and Xbound are the volume fractions of free water (including ice) and bound water, respectively. From Snowflow, pb is the bulk (in situ) dry soil density, p, is the intrinsic

64 density of water, x, is the total soil moisture volume fraction, and au, and /3,u are parameters for calculating unfrozen water content. 4.4 Discussion The combined Snowflow-Esnow model operates as a responsive system with atmospheric data as the driving force. Consequently, complete testing and validation of Esnow can only be done in conjunction with Snowflow and in comparison to a natural snowpack interacting with the same atmosphere. Chapter 7 discusses these comparisons after a description of the REBEX-1 radiobrightness experiment in Chapter 5. Esnow's self-consistency can be tested independently by taking a typical Snowflow snowpack profile and setting the soil and snow temperatures to Ttest, and setting the sky brightness from all directions at the top of the snowpack to: TnkY = Testtsky (4.82) where tsky represents the unit vector of 7Tky. Then the Esnow snowpack brightness must also be equal to Ttest at all angles and frequencies according to Kirchhoff's Law. In fact, Esnow's brightness was within 0.02 K at the observation angle of 53.1~ and frequencies of 19, 37, and 85 GHz.

CHAPTER 5 GROUND BASED RADIOBRIGHTNESS OBSERVATIONS IN THE NORTHERN GREAT PLAINS: THE FIRST RADIOBRIGHTNESS ENERGY BALANCE EXPERIMENT 5.1 Introduction This chapter describes the experimental apparatus, methodology, and measurements from the first Radiobrightness Energy Balance Experiment (REBEX-1). REBEX-1 was the first in a series of field studies designed to track the microwave radiometric response of terrain to antecedent weather. The purpose of REBEX-1 was to examine the link between radiobrightness and land-atmosphere energy fluxes in the northern Great Plains through the course of wintertime freezing and spring thaw. During REBEX-1, three microwave radiometers measured the apparent radiobrightness of a grassy site at 19, 37, and 85 GHz. Augmenting these data were measurements of sky radiobrightnesses, terrain and sky infrared radiometric temperatures, net and global radiation, soil temperatures, soil heat flux, rainfall, air temperature, relative humidity, and wind speed. In addition, a video camera and digitization hardware acquired 100 images of the radiometer observation area during the experiment for later use in evaluating snowcover conditions. REBEX-1 ran from October, 1992 through April, 1993, making 17200 observation 65

66 cycles encompassing vegetation senescence, snowpack formation, soil freezing, and thaw. The study site was on the grounds of the EROS Data Center (EDC), U. S. Geological Survey, near Sioux Falls, South Dakota at 43043/ N latitude and 96030' W longitude. This chapter describes the experimental apparatus, installation, and the data collected and discusses post-experiment error handling. 5.2 Apparatus REBEX-1 had two major instrument subsystems: the Tower Mounted Radiometer System (TMRS), with microwave radiometers designed and built as a part of this thesis, and the micrometeorological subsystem (MMS), a collection of commercially available instruments for monitoring local weather conditions. The instrument subsystems were integrated around a computer-automated data acquisition and control system. More detailed descriptions of these systems can be found in [58]. 5.2.1 Micrometeorological Subsystem and system integration Table 5.1 lists the specifications of each MMS instrument along with the parameters measured. The Infrared Temperature Transducer-hereafter referred to as an IR radiometer-produced a temperature output computed from a 15~ field-of-view thermal infrared radiometric measurement and an assumed target emissivity of 0.95. The IR radiometer and the other MMS instruments were factory-calibrated. I was able to confirm only a few of these calibrations independently. I checked the soil thermistor calibration and signal conditioning circuitry by immersing the probes in an ice water bath. All thermistor channels reported the ice bath temperature to be 273.15 K to within 0.1 K. I confirmed rain gage and anemometer operation by manually actuating switch closures in each. During the experiment, nighttime global radiation values

67 Figure 5.1: Interior of the trailer on site sheltering data acquisition and device control electronics and the Macintosh computer running the FluxMon HyperCard stack which controlled the experiment. were between 0 and 2 W/m2 and the maximum relative humidity value was 101.5%, providing indirect confirmation of the calibration accuracy of the pyranonmeter and humidity probe, respectively. Figure 5.1 shows the interior of the small heated trailer (1.5 m x 2.4 m floor dimensions) on site that sheltered the data acquisition and experiment control electronics. An Apple Macintosh II computer controlled all aspects of the experiment and provided a modem linking the experiment to offices at the University of Michigan. A custom designed program called FluxMon-operating in the HyperCard software development environment-managed data acquisition from all devices except the video camera, automatically controlled power to the instruments and heaters, and provided a graphical interface for control parameter adjustment and manual radiometer calibration. FluxMon communicated with the IR radiometer via formatted character

Instrument Model Parameters Accuracy Range REBS Net Radiometer Q-6 Net radiation1 N/A N/A Eppley Black and 8-84 Global radiation ~1.0% 0-1400 W/m2 White Pyranometer (shortwave)' Thornthwaite 610 Soil heat flow N/A N/A Soil Heat Flow Disk at 2 cm depth Met-One Anemometer 014A Wind speed at 10 m height'1 1.5% or 0.11 m/s 0.447-45 m/s Texas Electronics 525 Rainfall ~1.0% 0-5.1 cm/hr Tipping Bucket Rain Gage Vaisala Relative HMP35AC Relative humidity ~ 2% RI over 0-90% RH 0-100% RH Humidity Probe + ~3% RH over 90-100% RH| Campbell Scientific 107 Air temperature at 1.8 m ~0.2 K 240-321 K Thermistor Campbell Scientific 107B Soil temperatures at 2, 4, 8, ~0.2 K 240-321 K Thermistor 16, 32, and 64 cm depths Everest Interscience Infrared 4000ALCS Thermal IR surface and ~0.5 K 243-373 K Temperature Transducer sky temperatures 'Both instantaneous and experiment cycle average values of these parameters were recorded. o0 00 Table 5.1: Microrneteorological instruments and parameters. N/A indicates data not available.

69 strings transmitted through one of two Macintosh serial communications ports. A National Instruments NB-MIO-16 board with an AMUX-64T multiplexer provided 32 differential analog to digital conversion (ADC) channels and a TTL (transistortransistor logic) counter/timer channel. NB-DIO-24 and NB-TIO-10 boards provided TTL input/output and additional counter/timer channels. All of the boards fit into internal NuBus expansion slots in the Macintosh. The NB-MIO-16 digitized the signals from all of the instruments with voltage outputs: the microwave radiometers, the internal radiometer thermistors, and most of the MMS instruments. I used thermistors in bridge circuits for all temperature measurements because of their accuracy and ease of use. Instrumentation amplifiers conditioned the rest of the MMS voltage signals for ADC. The NB-DIO-24 TTL output channels controlled power relays for the radiometers, the radiometer heaters, the radiometer fans, the humidity probe, and the motor that opened and closed the radiometer housing door. NB-DIO-24 TTL input channels read signals from three microswitches indicated the door's position: fully opened, fully closed, or opened to the sky reflection position. The NB-TIO-10 counter/timer channels generated TTL square wave signals for setting radiometer heater power levels. FluxMon determined heater duty cycle settings approximately once per minute and reset the timer channel outputs based on a ten second total period. Counter channels counted switch closures from the anemometer and rain gage. FluxMon calculated wind speed as a function of number of switch closures over an elapsed time. For the rain gage, each switch closure was equivalent to 0.245 mm (0.1 in) of rain. Timbuktu/Remote software enabled remote control of the experiment from Michigan. Remote control procedures included trouble-shooting observations and control software changes, data file and video image down-loading to Michigan, and manual

70 Frequency (GHz) 19.35 37.0 85.5 Wavelength in air (mm) 15.5 8.1 3.5 IF bandwidth (MHz) 10-250 100-1000 100-1500 Radiometric resolution (K) 0.61 0.82 N/A Mixer operation Double-sideband Polarization Horizontal Integration time 6 s Antenna 3 dB beamwidth 10~ Incidence angle 530 (terrain brightnesses) Nominal zenith angle 450 (sky brightnesses) Table 5.2: Microwave radiometer specifications. N/A indicates data not available. See Appendix B for radiometric resolution calculations. control of the video image acquisition software. 5.2.2 Design of the microwave radiometers and the Tower Mounted Radiometer System Table 5.2 lists the specifications of the TMRS microwave radiometers. The radiometers simulated the observation angle, bandwidths, and three of the four frequencies of the spaceborne Special Sensor Microwave/Imager (SSM/I), a Defense Meteorological Satellite Program instrument. TMRS measured both terrain apparent radiobrightnesses and sky radiobrightnesses at 19, 37, and 85 GHz. I manually calibrated the microwave radiometers at the beginning of the experiment using ambient and liquid nitrogen temperature microwave absorbers (Eccosorb). In addition, the system automatically executed gain recalibration of the radiometers during the experiment using internal noise reference sources (matched microwave loads). Figure 5.2 shows the TMRS housing at the top of the REBEX tower. The towerbased electronics were divided into five modules: one each for the 19, 37, and 85 GHz radiometers, one for the IR radiometer and video camera, and one in the back of the housing's center compartment for the electrical bus. In addition, a motor

71 19 GHz Radiometer 85 GHz Radiometer Infrared Radiometer Door & Sky Reflector Figure 5.2: The TMRS-1 radiometer housing.

72 Figure 5.3: The 85 GHz radiometer. The 19 and 37 GHz radiometers have layouts which are comparable component by component. and screw drive mechanism in the center compartment positioned the door for sky brightness measurements during each experiment cycle. Figure 5.3 shows the 85 GHz radiometer module. The 19 and 37 GHz radiometers have similar component layouts. In each radiometer module, a mixer down-converted the RF signal to IF which then passed two amplifier stages and a bandpass filter. Three 12.2 m coaxial cables carried the IF signals from the tower to the trailer. The three radiometer modules differed only in the frequencies of their RF (radio frequency) and IF (intermediate frequency) components and in the voltage level of the regulators for their local oscillators. For each radiometer, a square law detector converted the IF signal from the tower to AF (audio frequency, in this case 0-20 kHz). AF amplifiers in a temperature controlled compartment then conditioned these signals for ADC by the NB-MIO-16. The NB-MIO-16 sampled the AF radiometer signals separately on three ADC channels at 40 ksamples/s for 6 seconds. FluxMon then calculated radiometer output values in instrument counts depending on which of two possible radiometer modes was activated-total power mode or Dicke mode. In total power radiometer (TPR)

73 mode. FluxMon calculated radiometer outputs by simply averaging each data stream. In Dicke radiometer mode, a 1250 Hz TTL signal generated by the NB-MIO-16 counter/timer modulated the RF input between the antenna and the internal reference load. This TTL signal also triggered ADC sampling, synchronizing it to the RF input modulation. FluxMon calculated the Dicke-mode radiometer output values by first numerically demodulating the data streams and then averaging. During automatic radiometer operation, FluxMon used TPR-mode for gain recalibration measurements off the internal reference loads and Dicke-mode for brightness measurements. Complete manual calibration of the microwave radiometers required TPR-mode and Dicke-mode measurements (in ADC instrument counts) of (i) a microwave absorber soaked in liquid nitrogen (LN2) and (ii) an absorber at ambient temperature. TPR measurements of the internal reference noise sources also made at calibration time established a baseline value for radiometer gain drift corrections. The radiometers were always under computer-controlled temperature stabilization during calibration-that is, FluxMon automatically measured temperatures and set radiometer heater duty cycles about once every minute. The calibration data acquisition procedures built into FluxMon also triggered measurements of the radiometer internal reference load and antenna temperatures with each data sample. A thermistor embedded in the ambient temperature absorber registered its temperature. The temperature of the LN2 soaked absorber is fixed at about 80 K by the liquid-gaseous phase change. 5.2.2.1 Microwave radiometer calibration Figure 5.4 is a block diagram showing the components of the microwave radiometers, all of which followed the same basic design. During the experiment, the calibra

74 Figure 5.4: Microwave radiometer block diagram.

75 tion parameterization included estimates of transmission line losses from the antenna to the receiver and from the reference load to the receiver. Arbitrarily defining the receiver as beginning at the output port (port 1) of the RF switch (a latching ferrite circulator) we have the following forward radiometer equation for the radiometer output signal: VD = CsD (L TAP + -r)Tpa- Doff) (5.1) L21 L21 where VD is the measured Dicke radiometer output in counts, TAP is the apparentor radiometric-temperature being measured, L21 is transmission line loss from the antenna to the receiver, rji is the antenna radiation efficiency, Tpa is the physical temperature of the antenna, and CsD and Doff are the Dicke-mode gain and offset parameters, respectively. Inverting (5.1) for apparent radiobrightness: TAP = - (1 - L)Tpa + Doff) (5.2) 711 CsD L21 CsD and Doff were found through a two point Dicke-mode calibration: Doff = - Dof f =1 VDO(TAP1l ~ (1 - r1)Tpal) - Vbl(TAPOnlz + (1 77i)TpaQ) L21 VDO - VD1 (5.3) - VDO _ _ (5.4) D TAPO/IL21 + (I1 - 771)Tpao/L2i - Doff () where subscripts 1 and 0 indicate data from the ambient and LN2 temperature sources. I estimated values for L21, L31, and qr/, listed in Table 5.3, from manufacturers specifications. Precise radiometric measurements require a calibration curve established at a time as close as possible to the time of measurement. This is because radiometer outputs are sensitive to gain variations. The TMRS radiometers were most vulnerable to gain

76 Parameter 119 GHz | 37 GHz _ 85 GHz ri1 0.9 0.9 0.9 L21 1.1001 1.034854 1.13186 L31 1.18 1.034982 1.13494 Table 5.3: Estimated loss parameters. drift through the inevitable change in IF coaxial cable temperature and the commensurate change in loss through the cable. These cables were directly exposed to weather over most of their 12.2 m lengths. IF and AF amplifiers were also subject to gain drift over the seven month length of the experiment. By using the parameterization in (5.1), gain variation may be isolated to the parameter CsD. Reference load temperature and RF transmission line losses and temperatures are then the primary determinants of the radiometer offset parameter, Doff, and I assume these terms to be constant. To track gain drift, FluxMon automatically made TPR-mode measurements of the internal reference loads and their temperatures during each experiment cycle. Assuming that receiver noise temperature and DC offsets remained constant, the reference load gain parameter, CsREF, is: F = VREF (5E5) CsREF ref /L31 + Trec where VREF is the TPR-mode output in counts when switched to the reference load, Tref is the temperature of the reference load, L31 is transmission line loss from the reference load to the receiver, and Trec is the TPR-mode receiver noise temperature parameter which includes DC offsets in the AF amplifier. Each manual calibration determined Trec using LN2 and ambient temperature absorber data: 1 VANTO(TAP1rlI + (1 - rl)Tpal) - VANT1(TAPOr77 + (1 - rji)Tpa0) T Lre - VANT1 - VANTO (5.6) where the parameters are defined as in (5.3) except that VANT is the TPR-mode

77 output in counts when the RF switch is set to the antenna input. Dicke-mode and TPR-mode gains differ due to the differing mismatches and losses in the antenna and reference load transmission lines. The RF components preceding the receiver are passive and their losses and mismatches are constant. Provided that internal radiometer temperatures and, consequently, Trec are stable, the ratio of CsD to CsREF will remain constant in time. That is, relative gain variation in the TPR reference load radiometer can be used to track variation in Dicke-mode gain using the relationship: CsD(t) - CsREF(t) CsD(O) CsREF(O) where t indicates the time of the experiment cycle and t = 0 is the calibration time. Figure B.8 in Appendix B is a graph of the TPR gain factors, CsREF, measured during each experiment cycle using (5.5). In each experiment cycle, FluxMon calculated CsD from CsREF using: CsD(t) = CsREF(t CsD () (5.8) CsREF(O) I manually calibrated the microwave radiometers and changed the calibration parameters accordingly on the dates listed in Table 5.4. Appendix B contains data used to evaluate the accuracy and precision of the microwave radiometers. In summary, the radiometric resolutions (repeatabilities) of the 19 and 37 GHz radiometers were 0.61 and 0.82 K, respectively, based on 20 measurements of known sources. Insufficient data are available to give the radiometric resolution of the 85 GHz radiometer. The average calibration accuracies were 0.24, -0.61, and -0.53 K for the 19, 37, and 85 GHz radiometers, respectively.

78 Radiometer Day Time [ Doff Trec CSD/CSREF 19 GHz 279 1800 272.552 146.796 0.91918 309 1600 271.099 84.4456 0.91065 403 0100 270.710 68.48 0.91162 37 GHz 279 1800 310.754 2.76098 0.93106 288 1900 310.055 -10.007 0.93520 309 1600 309.465 -53.472 0.93452 403 0100 303.357 -38.874 0.94690 85 GHz 309 1600 292.432 225.183 0.91427 403 0100 276.149 122.348 0.92667 Table 5.4: Calibration parameters used during REBEX-1 from day indicated to day of next calibration. Some 37 and 85 GHz parameters were later modified. See Sections 5.5.3 and 5.5.5. 5.3 Installation Figure 5.5 shows the EDC site as seen from the east. The TMRS radiometer housing is positioned atop the REBEX tower. The housing was attached to the 9.14 m (30 ft) aluminum tower via a winched shuttle. I installed and calibrated the radiometers with the housing lowered and left it at full height during all experiment cycles. The housing was made from aluminum sheet welded to a tube frame and the bottom hinged door was stainless steel. Figure 5.6 shows the housing with its back cover removed, revealing the housing power bus, the protruding door motor mechanism, and the housing's center module. The bracket mounting the housing to the shuttle permits rotation of the housing into a vertical position for servicing and module removal. Figure 5.7 shows the distribution of the other instruments to the southeast of the tower. The radiometer observation area was kept undisturbed through the installation process. I chose an observation area to the southeast of the trailer and the MMS instruments so that wintertime prevailing winds from the northwest would not cause snow drifts on the site. A graduated 5.1 cm (2 in) diameter PVC pipe with alternating

79 Radiometer Housing Guy Wire Anemometer f Net Radiometer - eter Soil Temperature Probes (Buried) Pyranom Trailer \ Air Temperature & Soil Heat Flow Relative Humidity Probes Disk (Buried) Figure 5.5: View of the REBEX-1 site from the east.

80 Figure 5.6: View of the radiometer housing with its back cover removed. Inserted into the housing are (from the far side) the 85 GHz radiometer, the IR radiometer and video camera module, the center connector box, the 37 GHz radiometer, and the 19 GHz radiometer.

81 N 5 — m Guy wire Trailer Gu W,iji Air temperature &... relative humidity ly wire probes..~ 'Tower 57~ Guy wire Rain gage Net Pyranometer / radiometer * Soil 12 m temperature probes S (buried) Soil heat flow disk (buried) Radiometer observation area 4m Snow gage * /......-201-ft"'-'-, Figure 5.7: Plan view of the REBEX-1 site.

82 Figure 5.8: Graduated PVC pipe with alternating 1.27 cm (0.5 in) black and white stripes used for gaging snow depth from video stills of the REBEX-1 site. 1.27 cm (0.5 in) black and white stripes was installed within the video camera fieldof-view, as in Figure 5.8. Video images of this gauge were used to make the snow depth estimates shown in Figure B.9. I installed the soil heat flow disk at 2 cm depth below the soil surface and the soil temperature thermistors at depths of 2, 4, 8, 16, 32, and 64 cm, as shown in Figure 5.9. I chose an undisturbed area under which to install the temperature probes. Because the soil was obscured by grass roots and litter, identifying the surface was possible to within only about 0.5 cm. The temperature probes themselves were 1 cm in diameter and were inserted into the soil through six 46 cm long horizontal holes made in the side of a trench, which was then back-filled. I installed the soil heat flux disk under a separate undisturbed area by cutting into the sod and soil with a knife and inserting the disk horizontally.

83 Figure 5.9: Insertion of the soil temperature probes. At the time of this picture. I had already inserted the probes in the side of the trench and refilled it, burying the 64, 32, and 16 cm probe cables. Cables leading to the 2, 4, and 8 cm deep probes protrude from the side of the trench.

84 5.4 Experiment log The REBEX-1 field data report gives a detailed account of the experiment log [58]. This section summarizes key qualitative information from the experiment. All experiment dates in this thesis are given as day numbers from January 1, 1992 (day 1). For example, January 1, 1993 is day 367 since 1992 was a leap year. A REBEX-1 day number vs. calendar date chart is given in Appendix B. All times are Universal Time (UT) which is six hours ahead of Central Standard Time (CST) at the site. While the experiment was operational, FluxMon initiated measurements at preset times-initially at every 10 minute mark of the hour and then later every 15 minutes. Data sets were time stamped at the end of each experiment run which lasted approximately 5 minutes. I acquired video images via the telephone link infrequently at first but then almost every day when there was snow on the ground. Approximately 100 video images were recorded over the course of the experiment. The frontispiece of this thesis shows one such image from day 420 (February 23, 1993). The REBEX-1 field data report [58] contains copies of all the images. Setup of the experiment began on day 269 and installation was completed by day 271. Several equipment failures forced a delay in experiment commencement until day 279, when data taking began at 1805 UT. Chief among these were the failures of the IR radiometer and one of the three microwave radiometer IF detector units and general electrical bus noise. Rewiring the infrared temperature transducer onto an independent DC power supply circuit resolved that failure and rewiring one of the 85 GHz circuits resolved the DC problems. There were no spare IF detector units so I left the two good detectors in the more reliable 19 and 37 GHz radiometers and placed the faulty diode with the 85 GHz radiometer. Although the 85 GHz output was unusable, I left the 85 GHz radiometer installed to monitor the performance

85 of the complete electrical system. This was the configuration when the experiment began. Over the course of the experiment, cold weather periodically affected instrument performance. The rain gage was not heated and so did not record snowfall accurately, if at all, and only worked reliably in warm weather. Dew, frost, and snow interfered with operation of the net radiometer and pyranometer, covering the instrument domes and blocking radiation., A heavy frost in early February, 1993 apparently caused the seizure of the anemometer lasting from day 404 to 409. On day 289, the computer clock stopped at 0934 UT and did not resume until 1318 UT when a worker manually disturbed the computer keyboard or mouse. Similar clock stoppages occurred several times during the course of the experiment and were resolved by manual means each time. The problem seemed to be attributable to a suspension of normal computer time-keeping interrupt generation when large data streams were collected from the microwave radiometers by the NB-MIO-16. The problem was resolved in later implementations of the software by initiating a query of the system clock after each large data acquisition run. Table 5.5 summarizes hardware problems that lasted for significant portions of the experiment and divides the experiment into seven time periods by general site condition. More detailed descriptions of most entries can be found in [58]. Section 5.5.5 discusses the problems with the 85 GHz radiometer and Section 5.5.2 discusses TMRS housing door problems. The experiment itself lasted from day 279 to day 471. 5.5 Correcting radiobrightness errors The following sections describe post-experiment processing of the radiometer data. The analysis covers deletion of out of range values, identification of sky brightnesses

86 Days Site coverage conditions 279-306 No snow, green grass cover, unfrozen ground 307-315 Snow cover 316-337 No snow, dormant grass cover, unfrozen ground 338-401 Snow cover 403-407 Snow cover cleared manually exposing grass over frozen ground 408-452 Snow cover 453-471 Mostly snow free, dormant grass cover Days Hardware problems 279-295 IR radiometer: Some data missing due to serial communications error 295-394 IR radiometer: Data drop-outs of temperatures less than 249.9 K (-9.9~F) 279-309 85 GHz radiometer: Not installed due to missing IF detector 309-347 85 GHz radiometer: Installed but malfunctioning 347-403 85 GHz radiometer: Functioning but not manually calibrated 403-471 85 GHz radiometer: Functioning and calibrated 350-471 TMRS housing door: Frequently fails to open to sky reflection position 415-471 TMRS housing door: Frequently fails to close past sky reflection position 404-409 Anemometer: Not spinning due to frost Table 5.5: Summary of REBEX-1 site coverage conditions and hardware problems. corrupted by reflector positioning errors, 37 GHz radiometer calibration errors on day 403, estimation of actual sky brightnesses, 85 GHz recalibration for days 347 -403, an alternative calibration parameterization, and radiobrightness sensitivity to 77i assumptions. 5.5.1 Removing out of range radiobrightness values Both terrain and sky radiobrightnesses exhibited occasional out of range values, typically near 296, 320, and 332 K for the 19, 37, and 85 GHz radiometers, respectively. The spurious points occurred singly in either the sky or terrain measurement and were not accompanied by spurious TPR gain factor readings. The out of range values occurred about once every other day in each instrument. The spurious ra

87 diometric measurements were always greater than the corresponding infrared surface temperatures, so they were distinguishable as non-physical. I manually deleted these values from the data set. 5.5.2 Sky reflector positioning errors Sky radiobrightness measurements suffered from two major sources of error during REBEX-1-incorrect sky reflector position and inadequate reflector size. Section 5.5.4 discusses the later of these errors in detail. Motion of the TMRS housing door-which also served as a sky reflector-became erratic on day 351 due to wear in the screw drive mechanism. At that time, only door closures were affected, with incomplete closures occurring periodically and usually in groups. Commensurately, sky brightness measurements were at times unusually high but I did not tie them to the door problems until day 414. From day 351 to 414, if an experiment cycle began when the door was not completely closed, FluxMon would make sky measurements without moving the door from its current, partially closed position. Many flawed sky brightnesses from this period were greater than the concurrent terrain brightnesses and could be easily identified. However, the sporadic occurrences of sky-brightening clouds or precipitation make sky brightness inherently highly variable, corrmplicating the job of weeding out faulty data. Since terrain brightnesses are a combination of emitted radiation and reflected radiation from the sky, it is possible to use terrain brightness variation as an indicator of sky brightness validity. For example, if sky brightness jumps by 100 K from one measurement cycle to the next, one would expect terrain brightness to increase by 10 K if the reflectivity of the terrain were 10c%. Manual editing of sky brightnesses in the day 351 to 414 period used the following criteria: sky brightnesses were deleted when they exceeded terrain brightnesses or

88 when a sudden, large change in sky brightness was not accompanied by a commensurate change in terrain radiobrightness. When the data were ambiguous, I chose to remove the points rather than risk using bad data. In most cases, a group of removed points would include at least one very high value, or a group of low values would be bounded by large jumps not mirrored in the terrain brightnesses. It is likely that faulty data remain, however, since zenith angles between 12~ and the the nominal sky measurement angle of 45~ were possible and would have yielded cold sky brightnesses that were indistinguishable by our criteria from valid 45~ measurements. The door mechanism failed in a new way on day 414 when it would not open completely. Software changes at this time eliminated the need for subjective editing of sky measurements. During the sky measurement, the reflector was either in the correct sky measurement position or fully open. Erroneous sky measurements were close to terrain measurements and easily identified. From day 414 to the end of the experiment, I deleted a sky measurement set if both its 19 and 37 GHz sky brightnesses were greater than 0.95 times the respective terrain brightnesses. I checked this criterion against the data from before day 350 and found no matching points, indirectly confirming that no valid points from days 414-471 were removed. 5.5.3 Revision of the day 403 37 GHz radiometer calibration Following the experiment, the t 37 GHz radiometer physical antenna temperature (Tpa) values recorded during the day 403 calibration were found to be erroneous. Since no modifications had been made to the system following this calibration, the erroneous value is likely to have persisted through the end of the experiment. Physical antenna temperature affects apparent radiobrightness calculations only when it changes from its value at calibration time. Over the course of the experiment, the range of 37 GHz

89 physical antenna temperature was approximately 9 K. This variation is modified by the factor of (1 - ril)/L2 in (5.2) so the error could have a maximum effect of about 0.9 K on TAP. FluxMon did not save the physical antenna or reference load temperatures used in the calculation of TAP during the experiment. In order to correct the 37 GHz calibration for the period from day 403 to the end of the experiment, I estimated these temperatures using outside air temperature. The following steps gave corrected 37 GHz data: * Using data available from calibrations on days 279, 288, and 309, perform linear regressions between the 37 GHz Tref and Tpa and air temperature. * Calculate Tref and Tpa for each experiment cycle after the day 403 calibration using the factors from the linear regression. * Using the erroneous Tpa value (282.73 K), the original calibration factors from day 403, and CSREF values saved from each experiment cycle, estimate the radiometer output, VD, for each experiment cycle after the day 403 calibration by applying (5.1). * Use the raw calibration data from the day 403 calibration with the estimate for the correct TAP to calculate a new set of calibration factors (Doff, Trec, and CSD/CSREF). * Calculate 1VREF for each experiment cycle by using the original Trec value and solving (5.5). * Calculate new CSREF and TAP values for each experiment cycle from (5.5) and (5.2) using the new calibration factors and previously estimated Tref and Tpa.

90 Radiometer Day Time Doff Trec CSD/CsREF 37 GHz 403 0100 305.936 -41.7428 0.93759 Table 5.6: Corrected calibration parameters for the period from days 403 to 471. Table 5.6 lists the corrected calibration parameters. Of the approximately 5700 experiment cycles between day 403 and the end of the experiment on day 471 the largest correction by this procedure was -0.914 K and the smallest was -0.124 K. The average correction was -0.538 K. 5.5.4 Estimation of actual sky radiobrightnesses High REBEX-1 sky radiobrightnesses suggest that these measurements were corrupted when radiation from the terrain reached the radiometers. As seen in Figure 5.2, the 19 and 85 GHz radiometers were close to the edges of TMRS housing door. Although the size of the door accommodated 10~ field-of-view main antenna beams, it did not account for side lobes and near-field diffraction effects. Sky brightnesses are used to calculate terrain emission from apparent radiobrightness and modeled terrain reflectivity. Since reflected sky brightness is usually small in comparison to emission from the surface, in most cases an approximate sky brightness will introduce only second order errors in terrain emission calculations. It is therefore valuable to attempt to estimate sky radiobrightnesses from the available flawed sky and terrain measurements. To make a sky brightness estimate, consider the hypothesis that for each radiometer the sky reflector had a fixed efficiency, Trrfi, such that: TRFL = 7TrflTSKY + (1 - T7rfl)TTER (5.9) where TRFL is sky brightness measured with the reflector, TSKY is actual sky brightness, and TTER is apparent terrain radiobrightness. The hypothesis suffers from some

91 obvious flaws since it assumes that radiation corrupting the sky measurement originates from the same area and at the same incident angle as TTER. If we assume that the radiation source corrupting the sky measurements is of terrestrial origin, then TTER is at least a valid surrogate. If, however, the corrupting source is either emission by the TMRS housing components-including the reflector itself which at times was wet or ice covered-or radiation from angles near the horizon, then TTER is at best the same order of magnitude as the source and the estimate of TsiKy will be invalid. The 37 GHz radiometer is most likely to fall into this second category since it was positioned near the center of the TMRS housing and did not have direct line of site to the terrain except at high incidence angles. I estimated r7rfi for each radiometer by calculating TSKY from Huron, SD rawinsonde data using radiative transfer. (See Section B.3 for a description of the data.) Rawinsonde profiles from Huron were available twice daily at 1100 and 2300 UT. I selected 20 profiles corresponding to clear sky periods as reported in both the Sioux Falls and Huron LCD's. TSKY is the sum of TDN, the downwelling radiobrightness from the atmosphere, and cosmic radiation attenuated by the atmosphere. Calculation of TDN is based on solution of the equation of radiative transfer [50]: dB +B =J (5.10) dT where B is brightness (W/m2sr) and J is the effective total source function at a point r in direction r, dr is incremental optical depth defined as dr = e dr, (5.11) and Ke is the total extinction coefficient. At 19, 37, and 85 GHz, the atmosphere is scatter-free in the absence of clouds and precipitation, i, = 9g, the coefficient for absorption by atmospheric gases, and J = Ja, the isotropic absorption source

92 function. Using the Rayleigh-Jeans approximation to Plank's blackbody radiation law, we can express Ja in terms of air temperature, T(r): JJ(r) = A2kT(r) (5.12) where k is Boltzmann's constant, A is wavelength, and Af defines the bandwidth. (In general, A is a function of r but its variance is negligible in the atmosphere.) Similarly, we can use Rayleigh-Jeans to define apparent radiometric temperature, TAP, in terms of B: 2k B(r)rf= A/|TAp(r)r. (5.13) Solving for downwelling apparent radiometric temperature, TDN, at the surface from zenith angle 0 in a plane-stratified atmosphere, we have: TDN(O) = sec 0/ j 9,(z')T(z')e-T(Oz')sec dz (5.14) where r(0, z')= j g(z)dz (5.15) and z is vertical height in the atmosphere. Oxygen and water vapor are the predominant absorbers in the microwave region so _g(z) is a function of temperature, pressure, and atmospheric water vapor content. The semi-empirical formulations in [50] give Kg. The solution to (5.14) can then be found by integrating numerically over the range of altitudes available in the Huron atmospheric profiles. Before integrating, I increased the vertical sampling density of the rawinsonde data by 20 times using cubic spline interpolation. Without the denser data, the non-linearity of (5.14) produced erroneous results when tested with a modified isothermal Huron atmospheric profile.

93 19 GHz 37 GHz 85 GHz rlrfi 0.701 0.723 0.636 Number of samples 20 20 7 cr 0.02 0.06 0.12 Percent TSKY < 0 8.1 17.7 0.1 Table 5.7: Reflector efficiencies from (5.18), the number of clear sky profiles used, and the standard deviation of rlrfl. Also shown are the percentage of TSKY values from (5.19) that were less than zero. I calculated TSKY for each of the 20 profiles using: TSKY = TDN + TCOS/L(0) (5.16) where Tcos is the cosmic background radiation (2.7 K) and L is the atmospheric loss at zenith angle 0: L(O) = e(Oo0o)sec. (5.17) Inverting (5.9) for 77rf:, we have: TRFL - TTER 7rfri = (5.18) TSKY - TTER I then calculated r1rfl using all 20 clear sky points for the 19 and 37 GHz radiometers and only the last seven for the 85 GHz radiometer, which was not calibrated until day 403. Table 5.7 gives the resultant reflector efficiencies. Having calculated the reflector efficiencies, I applied the following equation for estimated sky brightness, TSKY- to the entire REBEX-1 data set: TSKY = (TRFL - (1 - r7rfI)TTER)/17rfi (5.19) For each radiometer there were a percentage-given in Table 5.7-of TSK1-Y values less than zero, a non-physical result. This percentage was greatest at 37 C0Hz. suggesting that the reflector efficiency hypothesis expressed by (5.9) was least applicable

94 to the the 37 GHz radiometer. TSKY and TRFL are graphed in Figures B.3 through B.5 of Appendix B. Figure B.4 shows that the negative 37 GHz TSKY values are mostly concentrated in the fall and spring whereas the later chapters of this thesis focus on the winter months. The remainder of this thesis uses TSKY as a surrogate for true sky brightness measurements with the knowledge that TSKY is a first order approximation. 5.5.5 Calibration of the 85 GHz radiometer for days 347-403 As discussed in the experiment log (Section 5.4), the 85 GHz radiometer was unused from the beginning of the experiment until day 309 when I installed a new IF detector unit and calibrated the radiometer. The instrument quickly drifted out of calibration as seen in the TPR gain factor plot (Appendix B, Figure B.8). On day 347 the TPR gain factor data reveal an apparent transition to a functioning state. Since there is no explanation for either the initial drift or the day 347 transition, the day 309 calibration is not valid for the operational period after day 347. This section describes a means of determining an 85 GHz radiometer calibration set for this period. Equation (5.1) describes the forward radiometer parameterization for Dicke-mode output: (r (1 - /)Tpa _ ) VD = CsD (TAP -- + L -Doff) (0.20) L 21 L21 and (5.5) gives the TPR-mode reference load parameterization: Tref /REF = CsREF( + Tre). (5.21) L31 Before recalibrating the radiometer. I first evaluated (5.20) and (5.21) to recover the radiometer outputs originally measured during each experiment cycle. The 85 GHz

95 radiometer took data from days 347 to 403 using the day 309 calibration parameters given in Tables 5.3 and 5.4. The remaining terms in (5.20) and (5.21)-the temperature variables Tpa and Tref-were derived from air temperatures using the method described in Section 5.5.3. Having recovered VD and VREF, a new calibration set for the day 347-403 period was needed to recalculate TAP. I assumed that the ratio CsD(O)/CsREF(O) in (5.8) remained constant between this period and the day 403 calibration and that CSREF could be found for each experiment cycle by inverting (5.21) once Trec was found. Of the remaining parameters-Doff and Trec-the data in Table 5.4 shows that Trec varied between calibrations more than Doff for the 19 and 37 GHz radiometers. Consequently, I chose to use Doff from the day 403 calibration in recalibrating the day 347-403 data, leaving Trec to be determined. From day 347 through 403, the only possible calibration sources were the eleven sky brightnesses, TSKY, calculated from clear-sky Huron rawinsondes as described in Section 5.5.4. The inefficient sky reflector-also described in Section 5.5.4 -complicated the calculation of Trec using these TSKY. We can solve for Trec by inverting the TPR-mode reference load parameterization, (5.21): Tr e VREF Tref (.22) CsREF L31 where (from (5.8)): CsREF(0) CsREF = CsD Cs(O) (5.23) To calculate Trec from (5.22), we must find CsD using TSKY values from the Huron rawinsondes. We begin with the hypothesized formulation for reflector-measured sky brightness, (5.9): TRFL = T/rflTsIKY + (1 - r]rfl)TTER. (5.24)

96 For each rawinsonde-estimated TSKY, we have Dicke-mode measurements of TRFL and TTER from the coincident experiment cycles that we can use in (5.20): v r (-r)l T p -Drll VRFL = CsD (21TRFL + 21 Tpa- Doff (5.25) 21 L21 5.2 VTE TER = CD (LTER + 1 Tpa - Doff (5.26) LV21 L21 Solving (5.24), (5.25), and (5.26) for CsD, we have: VRFL - (1 - Tlrfl)VTER CsD -- (5.27) 7rll7rfIlTsKY/L21 + (1 - I77)r rflTpa/L21 - r rflDofff and we can then calculate Te,, from (5.22) and (5.23). The average Trec for the eleven rawinsondes was 116.359 K with values ranging from 92.195 to 139.247 K. Trec was 122.348 K in the day 403 calibration (Table 5.4), for a difference of 6 K. For comparison, the 19 and 37 GHz Trec's changed by -16 and 14.6 K, respectively, between their day 309 and 403 calibrations. Table 5.8 summarizes the revised 85 GHz calibration parameters. I used the new calibration parameters and the 85 GHz radiometer outputs-from (5.20) and (5.21)-to calculate TTER, TRFL, and TSKy data for days 347-403. Plots of these parameters appearing in Appendix B use the recalibrated data only. Figure B.8 of Appendix B shows both the original and recalculated TPR gain factor, CsREF, between days 347 and 403. The day 347-403 recalculated CsREF values are on average somewhat higher than those after day 403 but there is some crossover. Since there is no check of the recalibrated 85 GHz radiobrightnesses except by the same model used to derive the calibration, the recalibrated 85 GHz data should be considered approximate and of lesser reliability than the 19 and 37 GHz data over the same interval.

97 Radiometer Days Doff Trec CSD/CSREF 85 GHz 347-403 276.149 116.359 0.92667 Table 5.8: Revised 85 GHz calibration parameters for days 347-403. 5.5.6 An alternative calibration parameterization The REBEX-1 calibration parameterization used approximate values for the transmission line loss factors, L21 and L31, and ignored the effects of mismatch differences between the antenna and reference load radiometer transmission lines. Because radiometer output is linear in apparent radiobrightness, only two parameters need be determined through calibration in order to make an immediate radiobrightness measurement of some unknown source. In this case a very simple parameterization for the radiometer equation, (5.1), could be used, for example: V = ATAp+ B (5.28) where A differs from the gain parameter, CsD, in (5.1) by constant terms only while B includes variable gain, physical antenna temperature, and receiver noise terms. I used the parameterization in (5.1) instead of (5.28) because in practice each of the terms in B can be treated separately. This effectively isolated the most variable factor-gain-into one term so that it could be corrected through experimnent-time TPR reference load measurements. It was not necessary, however, to include the loss parameters L21 and L31 in the parameterization. L21 may be eliminated from the Dicke-mode calculations in (5.2) through (5.4) by redefining C/D = CSD/L21 and D'1 = DoffL2. Since L21 is constant, it has no time-varying effect on (CsD, that is, the ratio CsD(t)/CsD(O) in (5.7) will remain the same despite the transformation. L21 and L31 may be similarly removed from the TPR calibration equations (5.5) and (5.6). Because the transmission lines for the antenna and reference load were at nearly

98 Radiometer Day Time Do Te CSD/CSREF 19 GHz 279 1800 299.982 161.896 0.9635 309 1600 298.227 96.678 0.96266 403 0100 297.097 72.500 0.96763 37 GHz 279 1800 321.333 2.62806 0.93155 288 1900 320.862 -10.3558 0.93532 309 1600 320.196 -55.563 0.93704 403 0100 316.599 -41.7428 0.93759 85 GHz 403 0100 313.057 139.472 0.92653 Table 5.9: Alternative calibration parameters. the same temperature as the reference load itself, self emission by the transmission lines compensated for loss of signal from the load. Consequently, the alternative parameterization: C'/EF VREF (5.29) CRF = Tref + rec29) is a better model of the TPR reference load mode of the radiometer, although the (false) assumption is made that Tre, is the same for the antenna and reference load radiometer modes. We may justify this assumption-made in this and the original calibration-by taking (5.29) to be an index of gain variation but not a true radiometer equation for the TPR reference load mode. The effect of the alternative calibration on the other calibration equations is simply as if the values of L21 and L31 parameters were set to unity, so these equations will not be restated here. The alternative calibration parameters were calculated from each set of manual calibration data and the results are listed in Table 5.9. Applying the alternative calibration changed the calculated value of TAP by at most 0.016 K. Its advantage is simplicity and the elimination of arbitrary constants.

99 5.5.7 Sensitivity of radiobrightness to antenna efficiency assumptions To examine the sensitivity of TAP to rl1, I recalibrated and recalculated 37 GHz TAP using alternative rli values of 0.95 and 0.85 and data from day 403 to day 471. With r/1 = 0.95, the average deviation from the 0.9 case was -0.024, the maximum deviation was 0.165, and the minimum deviation was -0.253. With r71 = 0.85, the average deviation from the 0.9 case was 0.028, the maximum was 0.283 and the minimum was -0.184. These deviations are less than the radiometric resolution of about 1 K which justifies the initial assumption of r77 = 0.9. 5.6 Discussion Appendix B contains both overview and month-by-month plots of all the REBEX-1 radiobrightness and micrometeorological measurements. Discussion of the data is left to Chapter 6, which compares REBEX-1 radiobrightnesses to coincident SSM/I measurements, and Chapter 7, which uses the Snowflow and Esnow models from Chapters 3 and 4 to analyze REBEX-1 observations of the snowpack.

CHAPTER 6 MICROWAVE RADIOMETRY FROM SPACE: THE SPECIAL SENSOR MICROWAVE/IMAGER 6.1 Introduction The last three chapters have presented the theoretical foundations for a dynamic snowpack radiobrightness simulation and experimental observations linking radiobrightness to atmospheric conditions. This chapter adds a third elementobservations of radiobrightness from a space-borne instrument, the Special Sensor Microwave/Imager (SSM/I). It is important to present SSM/I observations before comparing theory with experiment for two reasons. Firstly, the validity of any remote sensing model must be confirmed for the intended instrument-that is, a space-borne one. The SSM/I field of view encompasses a range of surface cover conditions that are not necessarily well-characterized by conditions at the REBEX-1 site. Consequently, there will be times when differences between REBEX-1 measured brightnesses and SSM/I brightnesses are large, and understanding these differences will help determine the generality of the REBEX-1 observations and the simulations based on them. Also, the space-borne instrument measures radiation from the earth and an intervening atmosphere whose effects must be minimized in studying emission from the terrain. Secondly, the SSM/I has both vertical and horizontal polarizations, complementing 100

101 Channel Pol. Pass-band 3 dB beamwidth (deg) EFOV on earth frequency (V/H) (MHz) E-plane H-plane H-plane cross-. along(GHz) IFOV IFOV EFOV scan scan 19.35 V 10-250 1.86 1.87 1.93 69 43 19.35 H 10-250 1.88 1.87 1.93 69 43 22.235 V 10-250 1.6 1.65 1.83 60 40 37.0 V 100-1000 1.00 1.10 1.27 37 28 37.0 H 100-1000 1.00 1.10 1.31 37 29 85.5 V 100-1500 0.41 0.43 0.60 15 13 85.5 H 100-1500 0.42 0.45 0.60 15 13 Table 6.1: SSM/I sensor specifications [59]. EFOV on earth is in kmn. the REBEX-1 radiometers which had only horizontal polarization. The addition of the v-polarized channels will contribute to evaluation of the snowpack simulation in Chapter 7. 6.2 Data from the SSM/I The SSM/I instrument is a part of the Defense Meteorological Satellite Program (DMSP) and four DMSP platforms have carried versions of it. SSM/I data are archived for civilian use by the National Environmental Satellite, Data, and Information Service, an agency of the Department of Commerce's Nation Oceanic and Atmospheric Administration (NOAA). Beginning in 1990, a joint NOAA-NASA (National Aeronautics and Space Administration) program called Pathfinder has been archiving the data for global change research at NASA's Marshall Space Flight Center (MSFC) Distributed Active Archive Center (DAAC). The MSFC-DAAC is an element in NASA's Earth Observing System (EOS) Data and Information System (EOSDIS). The MSFC DAAC provided all the SSM/I data discussed in this thesis. Table 6.1 gives the characteristics of the SSM/I sensor channels. The antenna beamwidth data are from antenna range measurements with the first SSM/I on the DMSP F08 platform but apply to later instruments as well [60]. The SSM/I is a

102 conically scanning instrument with a fixed incidence angle of 53.1~ at the earth's surface. Samples are taken in a 102~ arc centered on the ground track of the satellite. The 85 GHz radiometer takes 128 samples in each scan and the other radiometers take 64 samples in every other scan. The radiometer integration time, T, is 7.95 ms at 19, 22, and 37 GHz and 3.89 ms at 85 GHz. The F08 SSM/I channels had pre-launch radiometric resolutions of 0.45, 0.73, 0.38, and 0.73 K at 19, 22, 37, and 85 GHz, respectively. In Table 6.1, the instantaneous field of view (IFOV) is the beamwidth of the stationary instrument and the effective field of view (EFOV) on the earth's surface takes into account integration time and movement of the beam along the scan. The antenna H-plane corresponds to the along-scan direction. This analysis uses data from the DMSP F11 platform. During the REBEX-1 period, the F08 and F10 platforms also carried SSM/I instruments but the F08 85 GHz channel was inoperative and the F10 was in a high-eccentricity orbit, distorting the SSM/I field of view. The F11 crosses the equator at 17:04 (ascending) and 05:04 (descending) local time, and measurements at Sioux Falls, SD fall in the ranges 22:31 -23:51 and 11:42-13:03 UTC (Universal Time). SSM/I processing included (a) subsetting the data to a local region from the global data provided, (b) calculation of antenna temperatures from sensor counts, and (c) resampling of the data to a common earth-registered grid at a common spatial resolution. Appendix C gives the details of these processes. Before satellite and ground-based measurements can be compared, the effect of atmospheric interference must be accounted for. This is the subject of the next section.

103 6.3 Compensating for atmospheric attenuation and emission The SSM/I measures radiobrightnesses from above an absorbing and emitting atmosphere whose characteristics are continuously changing. Section 5.5.4 discussed calculation of the downwelling brightness temperature from the atmosphere, TDN(0), and the expression for upwelling brightness, Tup(O), is similar: Tup(0) = sec O t K(zT')T(z')e-T(z') sec dz' (6.1) where z is vertical height in the atmosphere, H is the effective top of the atmosphere, and T(Z', H)= j, (z) dz. (6.2) As discussed in section 5.5.4, the predominant absorbers in a cloudless atmosphere at SSM/I frequencies are oxygen and water vapor. Ulaby, et al. [50] give the total gas absorption coefficient, Kg as a semi-empirical function of air temperature, pressure, and water vapor. (Clouds and precipitation also contribute to atmospheric attenuation but will not be included in this analysis.) The apparent (radiometric) temperature, TAP, in the sensor field of view is given by:1 TTER TA + TUP = TA (6.3) La(O) where TTER is the idealized apparent temperature of the terrain without atmospheric interference and La is the loss factor of the atmosphere: La,() = e(0oZ)secO. (6.4) 1The apparent temperature, TAP, is strictly defined along a ray incident upon the sensor antenna and antenna temperature TA is the integral of TAP over 47r solid angle. In this analysis, assume TAP is constant over the main beam of the antenna and that the antenna sidelobes are negligible such that TA = TAP.

104 TTER is the quantity we seek for remote sensing of the land surface. Given estimates of TUP and La, TTER can be estimated from SSM/I antenna temperatures by: TTER = (TA - Tup)La(O). (6.5) Atmospheric temperature and water vapor pressure are parameters that can vary over location and time throughout the year, especially near the surface where high pressure and water vapor density increase the attenuation. There are several methods by which the atmosphere can be characterized for a particular time and place, including weather balloons (rawinsondes), atmospheric models, and remote sensing techniques. The analysis of TDN in section 5.5.4 used rawinsondes from Huron, South Dakota to characterize 21 clear-sky atmospheres near the REBEX-1 site. Now we use these same rawinsondes to estimate TTER from SSM/I TA measurements using (6.5). 6.3.1 Using the rawinsonde atmospheric profiles Figure 6.1 compares three types of brightness temperature measurements-terrain radiobrightnesses measured by the ground-based REBEX-1 radiometers (REBEX-1 TTER), SSM/I antenna temperatures with no atmospheric compensation (SSM/I TA), and SSM/I antenna temperatures with the Huron rawinsonde atmosphere removed via (6.5) (SSM/I TA - rawin. atm. = TTER). All the measurements were made under clear sky conditions within two hours of the launch time of the Huron rawinsonde at either 1100 or 2300 UTC. Because SSM/I coverage is incomplete, there are just 12 times at which SSM/I and rawinsonde measurements correspond. As discussed in Appendix C, resampling of the SSM/I 85 GHz channel produces both high resolution-that is, the original 85 GHz EFOV given in Table 6.1-and low resolution samples approximating the 19 GHz EFOV. 19 and 37 GHz resample to low resolution only. This section examines both low and high 85 GHz resolution levels.

105 300 - -250 - I-I i 200 -Co CO150 150 - 19 GHz. o'/ + +y4,,'+ 300 - v E 250 -(C.I 200 - CO 150 -co 0I 300 19 GHz,, +,,,~l,,o+, -. I I 150 200 REBEX-1 250 TTER (K) I I I 150 200 250 REBEX-1 TTER (K) 300 300 (a) (b) 37 GHz 300 - 250 -!2 200 -co CO 150 -,~-7 "-\ l - o/ unL 1-.. +. 7+o 300 - - 250-.d " 200-:Z HCo 150 -co.~'+ +t,.... I I I 150 200 250 REBEX-1 TTER (K) 300 I I 150 200 REBEX-1 250 TTER (K) 300 300 (c) (d) Figure 6.1: Comparison of SSM/I antenna temperatures and REBEX-1 terrain radiobrightnesses (on the abscissa) at times coincident with rawinsonde measurements. All brightnesses are h-pol. (a) 19 GHz SSM/I antenna temperatures, (b) 19 GHz SSM/I antenna temperatures with rawinsonde atmosphere removed (TTER) (see text), (c) 37 GHz SSM/I antenna temperatures, (d) 37 GHz SSM/I antenna temperatures with rawinsonde atmosphere removed (TTER)- (Continued on following page.)

106 300 - 250 - I2 200 -CD CO 150 - 85 GHz, low resolution..4.7 + + -^ +... 300 - v - 250-.C I" 200 -F:-, Cn 150 -cn 85 GHz, low resolution / /+.4+,.. 4-/,., /,, I I I 150 200 250 REBEX-1 TTER (K) 300 I I I 150 200 250 REBEX-1 TTER (K) 3 300 (e) (f) 300 -,250 - -- i 200 -0) CD 150 - 85 GHz, high resolution 300 - v I 250-.C " 200 -C1 -0 150 -0CD t'+ + -W + +,, 85 GHz, high resolution.." +.~,,+,, +-..-."+,,,'.,| I I i 150 200 250 REBEX-1 TTER (K) 300 I I I 150 200 250 REBEX-1 TTER (K) 300 (g) (h) Figure 6.1: (Continued from previous page.) (e) 85 GHz low resolution SSM/I antenna temperatures low resolution, (f) 85 GHz low resolution SSM/I antenna temperatures with rawinsonde atmosphere removed (TTER), (g) 85 GHz high resolution SSM/I antenna temperatures, (h) 85 GHz high resolution SSM/I antenna temperatures with rawinsonde atmosphere removed (TTER).

107 Channel e = TA - TTER e = TTER - TTER (GHz) e ce |e ce 19 -2.0 11.0 -3.5 10.9 37 -0.9 11.1 -4.4 10.8 85 (low res.) 7.4 13.6 -11.0 12.1 85 (high res.) 2.2 13.5 -16.8 11.2 Table 6.2: Average, e, and standard deviation, ae, of difference between SSM/I (TA and TTER) and REBEX-1 brightnesses. All brightnesses are h-pol. Table 6.2 gives the average difference between REBEX-1 (TTER) and SSMI/I (TTER or TA) brightnesses. Note that TTER is always less than TA. This occurs because the atmosphere emits in proportion to its temperature and to its rate of absorption, while radiation from the ground-initially scaled down from the surface temperature by an emissivity factor-is absorbed in the atmosphere at an equal rate. An analogous process explains why TTER(SSM/I) is often less than TTER(REBEX-1). The REBEX-1 site differed from most of the surrounding farmland in that it was covered with a thick mat of grass while the wintertime farms were mostly bare. Although perhaps physically colder than the underlying soil, the grass layer is an absorbing cloud with a negligible dielectric contrast to air. Consequently, radiometrically cold emission from the REBEX-1 ground-reduced by the high ground-air dielectric contrast-is attenuated in the grass layer while the grass emits proportionately to its temperature and the same attenuation/emission function. Table 6.2 indicates that, when the atmosphere is taken into account, the terrain is on average radiometrically colder in the SSM/I field of view than at the REBEX-1 site. But the average is by no means a complete description of the difference. Figure 6.2 plots the difference TTER - TTER over the duration of REBEX-1. As described above, in October and mid-December-before there is significant snowcover —grass cover makes the REBEX-1 site radiometrically warmer than the SSM/I field of view.

108 11/1/92 1 Date 1/1/93 3/1/93 I 5/1/93 I 20 0 -20 - -40 -L - 20 w 0 m w LU, -20 a: -40 CO U) 20 0 -20 -40 19 GHz ___ +, + -____, + + -, ' i I +I I -- 37 GHz + + + + + - + +.::;, + I 85 GHz | + I I + ' +* + + 85,,=,=. j., t + +-',..=,, +......... I I I I I I I I I ' I I I I 260 280 300 320 340 360 380 400 420 440 460 480 Day from Jan. 1, 1992 Figure 6.2: Variation with time of the difference between TTER from SSM/I and TTER from REBEX-1.

109 But snow accumulating in the second half of December seems to create the opposite effect-that is, the REBEX-1 site is radiometrically colder than SSM/I-but only at the lower frequencies (19 and 37 GHz). Section 6.4 uses additional data to address these questions. 6.3.2 Compensation without a priori information Rawinsondes are the most reliable way to characterize the atmosphere for estimating TTER but they are rarely available. Without a priori information about the atmosphere, this section attempts to derive a first-order correction method that uses only the data available from the SSM/I. For this purpose, the 19 and 22 GHz v-polarized channels are good candidates because 22.235 GHz is a water vapor absorption line but is close enough to 19.35 GHz that we can assume the brightness temperature of the terrain, TB, is the same for both: TB(19, V) = TB(22, V). (6.6) Applying (6.3), we have: (TA(19, 1V) - TUP(19, V))La(19, V) - RTsKY(19, V) = (TA(22, V) - Tp(22, V))L,(22, V) - RTsK(22, V) (6.7) where RTsKY is reflected sky brightness and TB = TTER - RTDN. (6.8) Since reflected sky brightness is a second order effect, we further assume that R = 0, such that for the true values: A19,22 = (TA(19, V) - Tup(19, V))La(19, V) - (TA(22, V) - Tup(22, V))La(22, V) _ 0. (6.9)

110 To find the unknowns in (6.9)-19 and 22 GHz TUP and LA-we seek a model atmosphere that minimizes A19,22. Using approximate expressions for the U.S. Standard Atmosphere from [50], we have the profiles: T(z) T()- z < 11km, T(z) T(11) 11km < z < 20km, (6.10) P(z) = Poe-/H3, (6.11) P = poe-z/H4, (6.12) where 7 is the temperature lapse rate (6.5 K/km), Po is sea-level atmospheric pressure (1013 mbar), po is surface water vapor density, and the scale heights are taken to be H3 = 7.7 km and H4 = 2.3 km. Of the two remaining unknowns, T(0) and po, atmospheric attenuation is most sensitive to po at these frequencies. Consequently, for our simple correction set T(0) = 270K. This is a reasonable first order approximation for wintertime temperatures in South Dakota. Now the model atmosphere is characterized only by the po that minimizes A19,22 in (6.9). But can a minimum always be found? Figures 6.3 shows the sensitivity of TUP to po for the model atmosphere with T(0) = 270K. As mentioned above, the difference Tup(22) - Tup(19) increases with po such that-if TUP were independently measured-the difference could be used to infer po.. Figure 6.4 shows a simulation of what is actually measured-antenna temperatures-given terrain brightnesses, TTER, of 200, 235, and 270 K. The simple correction method proposes to use the information in the difference TA(22) -TA(19) to infer po, but the simulation shows that when TTER is around 235 K the method fails because TA(22) - TA(19) is insensitive to po-that is, there is no information in the difference that can be used to determine po. Without po information, a reasonable first order assumption is that po = 0. This will yield the minimum atmospheric correction-one based on atmospheric oxygen

Ill. - 19 GHz 200-......... 22 GHz -— "'...... 37 GHz —. — - -- 85 GHz.." 150 — 85GHZ,~ ' ~.,,...ooO.' -., 10,.....* } 100 i........ 50............ 0 5 10 15 3 20 25 30 Pv(O) (gm/m ) Figure 6.3: Brightness temperatures upwelling from the atmosphere as a function of surface water vapor density. 270......'.......... _ 260 TTER = 270 K................. -: — = 19 GHz,% 15 ~ T - 5- * 22 G.z I- 230 -220 -............. + 240 - TTER=...................... t-. 230 - oO....o................... 2210- T^ = 200,................... < ~ T 200:. —210 -.T....... 200 -........... 0 5 10 15 20 25 30 Pv(O) (gm/m3) Figure 6.4: Simulated antenna temperatures as a function of surface water vapor density given terrain brightnesses of 200. 235, and 270 K.

112 Channel ATTER e = TTER- TA e = TTER- TA (GHz) rawin.-dry atm. (rawinsonde) (dry atmosphere) e 'e e O'e e oe 19 1.1 0.7 -1.5 0.33 -0.35 0.33 37 1.5 1.0 -3.5 1.8 -2.0 1.8 85 8.1 5.7 -16.8 8.5 -8.7 6.2 Table 6.3: Comparison of TTER from the rawinsonde and dry standard atmosphere methods showing average, e, and standard deviation, oe,, of the difference e = TTER(rawin.) - TTER(dry atm.). Also given are these statistics for the magnitude of the correction, that is e = TTER - TA. All brightnesses are h-pol. content only-and is unlikely to over-compensate for the atmosphere. Figure 6.5 compares SSM/I TTER calculated using the rawinsonde method in section 6.3.1 and the dry standard atmosphere method (T(0) = 270K). As in section 6.3.1, each graph has about 12 points where rawinsondes and SSM/I samples correspond. Table 6.3 summarizes the statistics of each correction. The relative magnitudes of the rawinsonde corrections are about 3%, 5%, and 16% of the range of TTER at 19, 37, and 85 GHz respectively. The dry atmosphere corrections are about a factor of two smaller, that is, a dry atmosphere can model about 50% of the brightness change imparted by an intervening atmosphere under clear sky conditions during winter. In warmer months and when there are clouds the dry atmosphere correction will under-estimate the correction to a greater degree, especially at 85 GHz. The rest of this thesis uses SSM/I TTER calculated by the dry atmosphere method with the caveat that the true terrain brightness may be-and usually will be-lower.

113 300 - -o E 250 - -i0 -' 200 -< 150 -Co 150 - 19 GHz I;1 I I I 150 SSM/I I I 200 250 300 TA - rawin. atm. (K) (a) 300 -^-o E 250 -<200 --- co 150 - 37 GHz 300 - 85 GHz /",..,,+.+ v E 250 -'200 -Co 0C 150 - +,,' + f-" - - - I I I I 150 200 250 300 SSM/I TA - rawin. atm. (K) I I I 1 150 200 250 300 SSM/I TA - rawin. atm. (K) (b) (c) Figure 6.5: Comparisons between SSM/I terrain brightness estimates, TTER, made with a dry standard atmosphere (TA - dry atm.) and a rawinsonde atmosphere (TA - rawin. atm.). All brightnesses are h-pol.

114 6.4 Comparison of SSM/I and REBEX-1 radiobrightnesses Figures 6.6 through 6.8 compare SSM/I and REBEX-1 h-polarized terrain brightnesses at 19, 37, and 85 GHz, respectively. The clear trend at 19 GHz is for better correspondence between REBEX-1 and SSM/I during the snow season (late December through March) than in spring and fall. This trend is consistent with snowfall correlation lengths at least on the order of the SSM/I EFOV (Table 6.1) such that the snow cover at the REBEX-1 site is characteristic of the surrounding region. In spring and fall, when there is no snow cover, the grass-covered REBEX-1 site is brighter than predominantly bare farmland of the SSM/I region. Section 6.3.1 described this brightening, which is due to emission from the absorbing grass cover. Grass brightening is seen to a lesser degree at 37 GHz (Figure 6.7) and not at all at 85 GHz (Figure 6.8), although no fall REBEX-1 data are available at that frequency. Figure 6.9 graphs the differences between SSM/I estimated terrain brightnesses and REBEX-1 terrain brightnesses, TTER(SSM/I) - TTER(REBEX). Section 6.3.1 presented similar graphs for the few times at which rawinsonde data were available and noted that the snowpack-which is well established by January-makes the REBEX-1 site radiometrically colder than the SSM/I EFOV at 19 and 37 GHz. For example, Table 6.4 summarizes the difference statistics for January indicating that the REBEX-1 site is on average 6.6 K colder than SSM/I at 19 GHz. The same effect occurs briefly even in early November when snow covered the REBEX-1 site for about 5 days. Figure 6.6 showed that increasing snow is associated with the coldest SSM/I brightnesses, and we know from REBEX-1 that the site was uniformly snow-covered through all of January. Consequently, these observations suggest that the SSM/I viewed a mixed scene with less snow cover on average than the REBEX-1 site.

115 CO a3) CO 0) ---------------------------------------- ----------------------------------------------------------------------------,in............. ------------------------------.............................................................. - --------------................. --------------------------------------------- ----------- -------------------- -jr; 0 CO 0) 0Y) CY) 0) 1-. 11 -C'j 0) 1-1 C'j 0 Col o C~ 0* 0 _ Co co co co CM C\j 0) 11 -03) 0 0 0 0 0 0 0 0 0 0 0 co (0 "I ~ 0 co CO "t CO C'J C~J C~J C'J Cj T rI T (>i) sseauLjt~uqoipej uei~ia I

116. -. C). - 0. = = —.......... --..... C)00...... Ale '- -.............................. f-....................................3 - ^g J," CL^ O X /,,,-.,., - 0).... - -- --...................... 9-<:'-. '"7 L - - ---------------------- = ---------------------------------------- -----------------. --- —----------------------.' O.,, NC 0 0 0 0 0 0 0 0 0 O CO (0 C\j O O0 0CD CY) C\J C CM C\ C\M - - 9 -(>I) ssaut6I!Jqo!peJ uiejjai

117 C) (O............................... - -- -- ---- ------ ------- - -- - --- --—......- ---- - - - -—.................................................................. 0)I..p. -O C *.::::-' ~ d - TI -.........T, (S.>. - 0.. ~~~~~~~~~~~~~~~~~~- — *- *-;w- t --- -- ------ -- -- - - ----------— * --- —-----— * --- —-------— * ----*,I............... -t; M.......:..... 0..................... o. CO C a)3 0......,,.,.,.,.,.*...-............................. ------------- -; i-WRR9 --- —— @- 0 -O CO::D:-C ~ I - -d Chi 2 ^ *,;;o CJ uv Ej 0:........... -- - -------------------------------------------------------------------------------- ------------------------------------ L c -..................c... - C o- o 0 0 C ) C\ C (> — ) SS-U- qO-p- - —!-jj.;.......... 0...::,,. 0 0 j CM 0~ 0) 00 0

118 Date 1/1/93 11/1/92 3/1/93 5/1/93 40 20 0 -20 I y LU IJ m cr a: LUJ I -40 40 20 0 -20 -40 40 20 -20 -40 44I I I480 440 460 480 260 280 300 320 340 360 380 400 420 Day from Jan. 1, 1992 Figure 6.9: TTER(SSM/I) polarized. Difference between REBEX-1 and SSM/I terrain brightnesses, - TTER(REBEX), at 19, 37, and 85 GHz. All brightnesses are h

119 Channel TTER(SSM/I) - TTER(REBEX) (GHz) e oe crae/DR 19 6.6 6.4 0.08 37 1.7 7.7 0.09 85 1.0 13.9 0.12 Table 6.4: Statistics of difference between REBEX-1 and SSM/I terrain brightnesses in January. All brightnesses are h-pol. DR is the January dynamic range of the REBEX-1 brightnesses. Other features of Figure 6.9 include: * Negative SSM/I-REBEX-1 differences in February at all frequencies. Recall from chapter 5 that the REBEX-1 site was cleared of all snow in early February, so even after new snow fell the SSM/I viewed on average a deeper-and radiometrically colder-snowpack than that at the site. * High variance in the differences in March. March brightnesses are wildly varying (at least on the time scale shown here) due to melting and refreezing in the snowpack-usually on a diurnal cycle. (Chapter 2 discussed how snowpack brightness increases with wetness.) The variation in the difference is probably a function of melt/freeze timing which is unlikely to be synchronized over the entire SSM/I EFOV. * Generally high variability in the 85 GHz difference. 85 GHz January brightness temperatures are usually the lowest of the three frequencies, and, consequently, 85 GHz snowpack reflectivity is commensurately high. Add to this the natural variability in 85 GHz upwelling and downwelling atmospheric brightness (see Figure 6.3) and the result is amplification of the variability seen at the other frequencies. Figures 6.10 through 6.12 at the end of the chapter compare SSM/I brightnesses

120 at v- and h-polarization. Since the REBEX-1 radiometers were single polarization (hpol.), SSM/I is the only source of v-pol. brightnesses available to us. The minimum polarization difference occurs in late September as thick vegetation cover obscures the soil surface. Polarization difference at 85 GHz is always small, suggesting that scattering and emission by the cover medium dominate the signal. Figures 6.10 through 6.12 show that polarization difference is always largest at 19 GHz. Wet soil dielectric constants are higher at 19 GHz and scattering by the cover medium is minimized by the longer wavelength. For many terrain types, the SSM/I incidence angle is close to the point of v-polarization total transmission-that is, the Brewster angle, OB, in the emitting medium. In air, the propagation angle corresponding to the Brewster angle is: tan Oair = c1/2 (6.13) where e is the dielectric constant of the emitting terrain and is equal to 1.77 for total transmission at 53.1~. At the SSM/I frequencies, the dielectric constant of soils may range from 3 for frozen soils to 15 for warm, wet soils [50], for total transmission angles of 63.4~ and 75.5~, respectively. But even at high dielectric constants, v-pol. emissivity at 53.1~ remains significantly higher than h-pol. The very large 19 GHz polarization difference during the snowpack season (days 350 through about 450) is more difficult to explain by this reasoning. The next chapter will address this issue further when modeling results are presented.

121 CY) a) 10. CO~ 0) c-) 0) 11.1 CO) 0) c'J 0) 1-1 COJ 0) 1-. CN 01) C15 Cf4 0) 1 -0rl 0 0 (0 0 0 CMl 0 0) 1~ E 0 0 0 Ct CO) 0 CO) 0 (0 cM 0 (0 CMj ct 4CD 0 -0~ 0 0 0 0 0 0 0 0 0 0 co CO I'l CMl 0 CCO (0 11 ClF) CMl CMl CMl CMj CMl?- I(>I) ssaujq46!Jqoipe4 uieaUg 0 0 0 Iq CMl (>i) H - A

122 C1) Lf)~ CL) 00 a) (0 > 0 I............~..... _..~.......................................... - --------- --- - CO O \] AC... C1) *..:CM CC' CC Ca *. aV) -..C-.. a C'..y.. 0... -D CC.) 0 COJ Co C\J N0 CM C~l \J 0 CO _ _ _ _ _ _ _ _ _ _ _ _ ( 0 0 0,-0 0 00 0 0 0 0 0 0 0 ( 0)0() ca 1 ~~~~ - - --- -- --- -- -- --- --- -- ------ -- --- --- - - - -- C)~ Cr)j 0 0 0 0 0 0 0 0 0 0 co(O t CM0 c (O t q I Co CVI C\J CVj C1J C' - I(>i) sseu)t6ijqo!pei uiejjejL (>) H -A

123 Co 0) r~ LO co C') Cl) a) CO (1) CM) 0 CO 0 (C 0T 01 0 CM\ 0 C\j To co n E 0 cOC C CO CM C' NZ CM 0) IrrCM C) C TT C 0 CO Cj (0 (NJ o o 0 o o o o o o o o o (C C\j C\J C\M C\M T- 2 (>4) ssuj46iuqoipej ulejjei (>) H -A

CHAPTER 7 COMPARING MODELS TO OBSERVATIONS 7.1 Introduction Chapters 3 and 4 presented the theoretical groundwork for a SVAT-linked snowpack radiobrightness model and Chapters 5 and 6 described an observational data set comprising continuous terrain brightness and micrometeorological measurements spanning 193 days. This chapter connects model with observation and, through the comparison, examines the physical processes that tie microwave emission to snowpackatmosphere fluxes. 7.2 Model initialization and inputs This section briefly describes how REBEX-1 data drive the Snowflow snowpack SVAT and the Esnow brightness simulation. Since Snowflow's inputs are taken directly from the REBEX-1 micrometeorological record and Esnow uses measured sky brightness to calculate terrain brightness, it is important to clarify which parameters are modeled and which measured. Some information here is repeated from Chapters 3 and 4 and the reader should turn to those chapters for more detail about the models. Also. Appendix A lists fixed Snowflow constitutive parameters. their symbols, and their values used in these simulations. Figure 7.1 diagrams inputs to the snowpack simulation. Snowflow's simulation 124

125 Initial state Soil temperature profile Fixed soil moisture No Snow Atmospheric inputs from REBEX-1 U3, Tair, RH, P, Qsw, Qlw Snowpack profile Ts, ps, Ws, gs, ds Sky brightness from REBEX-1...... - Esnow (TSKY) 19, 37, 85 GHz Terrain brightnesses (TTER) 19, 37, 85 GHz v- and h-polarizations Figure 7.1: Schematic of the Snowflow and Esnow inputs and products.

126 always begins with no snowcover and soil temperatures are the only input constraint on the initial state. For the results in this chapter, Snowflow initialized soil temperatures using REBEX-1 measurements at six soil depths, interpolating to the fixed depths of the model soil layers. For layers below the deepest REBEX-1 temperature measurement, Snowflow extrapolated temperatures to a maximum depth of 1.5 m. Initializing temperatures of deeper layers had minimal effect on surface temperatures in short term trials (2-3 model days.) As discussed in section 3.2, deep temperatures and heat flux will have an effect in the longer term. Of the fixed Snowflow parameters listed in Appendix A, soil moisture requires special attention. Snowflow treats total soil moisture as an invariant constitutive parameter of the soil when in fact it varies dynamically under the influence of moisture gradients, gravity, temperature, and freezing zones. This chapter discusses Snowflow model results with two soil moisture treatments that bracket the observed range. "Wet" soil has volumetric moisture content before freezing of x, = 0.43, which is the highest value measured during REBEX-1 and just under saturation at 0.45. "Dry" soil has xw = 0.20, which is far drier than the lowest REBEX-1 moisture (0.30). This discussion takes wet soil results to be the more realistic simulation and uses dry soil results to examine soil moisture sensitivity. In each model time step, Snowflow's only inputs are a set of atmospheric forcing variables-wind speed, air temperature, relative humidity, precipitation, shortwave radiation, and longwave radiation. All of these variables are taken directly from the REBEX-1 micrometeorology except longwave radiation which is estimated from air temperature and RH (see section 3.4.4.1). When gaps in the record longer than one hour occur, Snowflow uses micrometeorology from 24 hours before. All Snowflow model results have a nominal 15 minute interval.

127 Output from the Snowflow snowpack simulation includes soil temperatures at the soil surface and six depths and snow density. temperature, liquid water content. avxerage grain size. and thickness in every modeled snowpack layer. A maximum of 43 snowpack model layers were produced. Esnow uses these data to calculate snowpack emission and reflectivityr and then uses input REBEX-1 sky brightnesses to produce terrain brightnesses. In this chapter, Esnow brightness results have been calculated about every four hours. The model-experiment comparisons here are based on a Snowflow test period running from December 2, 1992 through February 4, 1993 (REBEX-1 days 337-401). First snowfall occurred on Dec. 3 but the snowpack was not permanent until Dec. 14. Terrain brightness comparisons here begin on Dec. 13. As described in Chapter 5. snow was cleared from the REBEX-1 site on Feb. 6 so comparisons to the model are invalid after this point. Focusing on the test period, the following sections evaluate first the snowpack simulation then the brightness simulation-which, of course. cannot be analyzed independently of the model snowpack. The last section discusses the link between radiobrightness and snowpack structure, grain size, and wetness. 7.3 Evaluation of the snowpack SVAT This section compares the Snowflow wet-soil snowpack simulation to observations of four variables from REBEX-1, all of '- ich are measured independently of the Snowflow inputs. These are snow depth (from Sioux Falls, SD LCD), surface temperature (from the REBEX-1 IR radiometer), and soil temperature and heat flux at 2 cm depth. The discussion also covers the sensitivities of soil temperature and unfrozen water content to total soil moisture. Table 7.1 summarizes the modeled-observed difference statistics for all four com

128 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 _,,,,,,,,,,,,,1,,,,,,,J,,,,,,!,,,,1,,,, ll,,,,,,,,!,,,,:- Snow depth... 0.2 Sioux Falls LCD?... E........... Snowflow model (wet soil) |. 0.1 A. -—. |. 0.1 Difference -0.1 - 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.2: Modeled and observed snow depths with wet model soil. parison variables. Figure 7.2 compares modeled and observed snow depths. The average aow depth difference is 0.02 m with a standard deviation of 0.027 m. For comparison, the measurement resolution of the LCD snow depth data is 0.025 m. The prominent feature of this comparison is the sharp difference on day 364. Before the snowfall here, the model already has 0.05 m more snow than observed. Then the model under-estimates the density of the new snow and diverges further from ooser _ VariableA or 0A/DR Snow depth (m) 0.02 0.027 0.12 Surf. temp. (K) -5.5 5.0 0.17 2 cm soil temp. (K) -1.5 0.96 0.64 ~2 cm heat flow (W/m2) -3.6 15.9 0.59 Table 7.1: Comparisons of modeled snowpack variables to REBEX-1 measurements. \ = model - observation. a is standard deviation, and DR is the dynamic range of the measured data over the test period.

129 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 290 I I I I IIII I I II I I I I I I I I I I,, I,,, I, I... Snowpack surface temperature REBEX-1 0.. Snowflow model (wet soil) 270 - 7280 E 260 111, V, irl. i '!.^ I. I I, I, I I I I I i I, ~ I I i I, I I p, i I: I I,.I i I' I I I I, H20 -Difference a 10 - '' ~ 1,, l..... 0 - -20 r4 J i l, I I I I Il,, I I [l, I,,,i,i l,i,,i, 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.3: Modeled and observed snowpack surface temperature with wet model soil. vations. New snow density estimates are based solely on air temperatures, reaching their minimum value at Tair = 258 K. Figure 7.3 shows that the snowpack surface temperature on day 364 was around this value. The snow depths are in closer agreement after the snowpack deepens on day 377, with the subsequent ablation occurring at comparable rates. Figure 7.3 compares surface temperatures from the Snowflow wet soil model to REBEX-1 observations. The average difference is -5.5 K with 5.0 K standard deviation. The negative bias could be caused by under-estimated values of Qlw. the longwave radiation from the atmosphere. Since this input parameter is estimated from air temperature and relative humidity only, the warming contribution of clouds

130 is systematically neglected. The other possible source of surface temperature error is the latent and sensible heat transfer parameterization, but there is nothing to suggest that such an error would have a negative bias. The remaining directly comparable Snowflow variables are soil temperature (Figure 7.4) and soil heat flux (Figure 7.5), both measured at 2 cm depth in the soil. For soil temperature, the average difference is -1.5 K with 0.96 K standard deviation. Although the magnitudes of these values are small they are significant nevertheless because the soil temperature is near the freezing point. As discussed in section 3.3.1, soil water freezes over a range of temperatures and this elevates the heat capacity of the soil and its thermal inertia over that range. Consequently, the soil temperature solution is sensitive to errors in the temperature dependence of the model soil's heat capacity. The average model soil temperature error is a direct consequence of negative bias in the model snowpack surface temperature, discussed above, but higher model variability may be caused by differences in either heat capacity or the balance of fluxes at the soil surface. The rapid divergence of model soil temperatures from measurements in the first model day-before REBEX-1 soil reached freezing-suggests that the balance of fluxes at the soil surface is the larger source of error. Figure 7.5 compares modeled and observed heat flux at 2 cm depth and shows the components of model flux at the soil surface: the net flux into the soil (Fg = Fsg + Fswg), the quasi-conductive flux between the snow and ground from section 3.4.2 (Fsg), and the shortwave radiative flux absorbed by the soil surface from section 3.4.3 (Fswg). The average difference between the modeled and observed 2 cm heat flux is -3.6 W/m2 with 15.9 W/m2 standard deviation. The strong diurnal cycle of Fs,,,g is the largest component of F9 and is the driving force behind diurnal temperature variation in the

131 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 278 Soil temperature at 2 cm REBEX-1 276-...... Snowflow model (wet soil) < 274 - i- ~- - - - -- - - - ----------------- - - - - - - - - - -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 272 - 270 _. ^.. 268 -266 -6 -4 - Difference Ia -2 (3D -2 ) - -6 --8 -8 - -', lJ[ ] i [ [ [ i,, i, I I [, i 1 i, I j il l,,, 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.4: Modeled and observed soil temperature at 2 cm depth with wet model soil.

132 Date 12/1 1/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 cs E x Ia I 60 40 20 0O -20 -40 c^s E x I i-s CZg NE I -1 - 100 - 80 -60 -40 -20 -0 - Shortwave flux at soil surface -—........... Snowflow model (wet soil) i!: L.; I_. -6.i %...P.Z I- -;- J t.- - - %.; .j Lz %.. ... -: 1- %..;;.- .;..; -.4 a.; - A.- ww................................................................. I rod ll I 1 -l 1 1 1 1 1 1 1 1 1 1 1 1 1 1 j 1 1 1 1 1 1 I 1 1 I 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 I 1 1 1 1 1 1 1 1 1 40 20 0 -20 -40 Snow-soil heat flux........... Snowflow model (wet soil) 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.5: Modeled and observed soil heat flux at 2 cm depth with wet model soil. Positive values indicate heat flow into the soil.

133 model soil (Figure 7.4). As discussed in section 3.4.2.1, Snowflow accounted for the grass mat by increasing the thermal resistance at the snow-soil interface but it neglected shortwave attenuation by grass inside and at the bottom of the snowpack. The lack of an observable diurnal cycle in soil temperatures suggests that little shortwave radiation penetrated to the REBEX-1 soil surface. The sensitivity of soil temperature to moisture content is demonstrated in Figure 7.6 which plots 2 cm soil temperatures from the dry soil model. The average difference with respect to the observed temperature is -2.8 K with a standard deviation of 2.4 K. Figure 7.7 plots the wet and dry soil unfrozen water content calculated as a function of temperature as in section 3.3.1. Although dry soil temperatures are up to 4.7 K colder than wet soil temperatures, all the moisture in the drier soil remains unfrozen through most of the test period. This is because the freezing point of the soil water is reduced to 268.0 K in the dry soil while it is 272.5 K in wet soil. Consequently, the heat capacity of the drier soil is comparatively small and its temperature reacts more quickly to changes in heat flux. 7.4 Radiobrightness comparisons This section compares the radiobrightness simulations of the linked snowpack SVAT and emission models to observed terrain brightnesses from the REBEX-1 ground-based radiometers and SSM/I. The discussion also covers brightnesses from the wet and dry soil models, the relationship between snow wetness and brightness, the relative effects of soil moisture and the snowpack at 19 GHz, and the sensitivity of brightness to snowpack layering. Figures 7.8 through 7.10 compare model terrain brightnesses at h-pol. with groundbased observations from REBEX-1. Table 7.2 summarizes the difference statistics

134 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 a) a) Q. E a, 278 276 274 272 270 268 266 6 4 2 0 -2 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.6: Modeled and observed soil temperature at 2 cm depth with dry model soil.

135 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 m E CO E E 0 > 0.4 0.3 0.2 0.1 0.0 Co E cso, CO E E > 0.4 - 0.3 - Modeled unfrozen soil water (dry soil) 0.2 - V I v V\N^AArr^,.-.. 0.1 - 0.0 I.I I. I.... I.....I. r..... 340 350 I.... II.... I....... I.... I.' 360 370 380 Day from Jan. 1, 1992 390 400 Figure 7.7: Modeled unfrozen soil water content for model soils with initial moisture volume fractions of 0.43 (top) and 0.20 (bottom) before freezing.

136 along with those for the SSM/I data, discussed below. The usefulness of a brightness comparison in a remote sensing context is the information content of the modelobservation difference. Model brightnesses at 19 GHz h-pol. (Figure 7.8) are at first colder than observed, then briefly warmer, and at the end of the test period the difference is highly variable. Since interaction with the dry snowpack is weak at 19 GHz, coldness in the 19 GHz model is likely caused by the dielectric contrast at the soil surface-that is, either the contrast is too high in the model or the un-modeled grass cover increases the observed brightness. A sharp dip in the observed brightness around day 370 may be the dielectric signal of wetness at the top of the soil or the effect of dielectric contrasts in the snowpack layers. Lastly, the high-variability period, which is distinct at 37 GHz as well, is a result of partial melt and refreezing cycles in the snowpack combined with wetness at the top of the model soil. These phenomena are all discussed in more detail below. Figure 7.9 compares 37 GHz brightnesses from the models and REBEX-1. Except for the middle period (days 365-385), there is little in the comparison to suggest that any particular process is unmodeled, although the simulation fails at times to match the details of the processes. In the early period (days 347-364) during which 19 GHz modeled brightness were cold, 37 GHz modeled brightnesses are unbiased. Since at 273 K the dielectric constant of water at 19 GHz is about twice that at Channel REBEX-1, h-pol. SSM/I, h-pol. SSM/I, v-pol. (GHz) /A UA oA/DR A CA A/DR A CA a/DR 19 -11.5 22.0 0.29 -11.0 19.4 0.38 -15.1 11.6 0.72 37 8.1 14.5 0.16 10.8 14.8 0.23 3.4 10.7 0.20 85 10.6 14.7 0.13 8.3 18.8 0.16 3.7 18.8 0.17 Table 7.2: Comparison of modeled brightnesses (K) to h-pol. observations from REBEX-1 and h-pol. and v-pol. observations from SSM/I. A = model - observation, a is standard deviation, and DR is the dynamic range of the measured data over the test period.

137 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 t 180- 19 GHz, h-pol. brightness. A........... Models (wet soil) p 'V 140 — 60 - 40 340 350 360 370 380 390 40 200 Figure 7.8: Modeled and observed 19 GHz h-pol. terrain brightness with wet model 37 GHz, the modeled dielectric contrast at the soil surface and sensitivity to soil 340moisture are reduced at 37 388 to 40 the 37 GHz signature in both the simulation and the measured data. Where the two diverge strongly, information about the timing of partial snowmelt may be inferred, as discussed below. Figure 7.10 compares 85 GHz modeled and observed brightnesses. The best match-and the least residual information-occurs between day 370 and 385 where 19 and 37 GHz brightnesses were least precisely modeled. The outstanding features of the 85 GHz data are that brightnesses are very low and variability very high. These

138 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 v 0) ((.) c~.cCQ c.r <d 260 240 220 200 180 160 140 60 40 20 0 -20 -40 -60 340 350 360 370 380 Day from Jan. 1, 1992 Figure 7.9: Modeled and observed 37 GHz h-pol. terrain brightness soil and ground-based observations. 390 400 with wet model features are related because, by Kirchhoff's law, low emissivity means high reflectivity, and the reflected source-the sky-is highly variable at 85 GHz, as discussed in Chapter 6. Emissivity is low at 85 GHz because of scattering in the dry snowpack that elevates extinction and decreases the effective depth of emission. Low emission depths mean that snow thicknesses beyond a threshold value have little effect on emission and neither do soil conditions. Figure 7.11 uses SSM/I terrain brightnesses (TTER in Chapter 6) to evaluate the v-pol. model performance. The model-SSM/I difference statistics (Table 7.2) indicate that the model deviates least from observations when compared to 37 and 85 GHz v

139 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 260 - 7 240 160- REBEX-1.Models (wet soil) 60 - Difference C' 20 - 1 -20 < 40 -6 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.10: Modeled and observed 85 GHz h-pol. terrain brightness with wet model soil and ground-based observations. pol. SSM/I brightnesses. In contrast, modeled 19 GHz v-pol. brightness is on average 15 K colder than SSM/I which is the worst model-observation mismatch. Since v-pol. radiation is least sensitive to dielectric contrasts (the Brewster angle effect discussed in section 6.4), the weak contrasts in dry snow due to layering are unlikely to affect it strongly. The closest 19 GHz v-pol. model-SSM/I brightness match for dry snow occurs around day 375, and this is when modeled unfrozen soil water content was at its lowest value (about 25% by volume in Figure 7.7). The sensitivity of all three frequencies to soil moisture is examined in Figure 7.12 which plots h-pol. brightnesses from the wet and dry soil models. The plots show the

140 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 7260 220 I -~~ 19 Gr-z, v-pol. brightness 200 -SSM/l co....... Models (wet soil) ~260 240 -C)........Models (wet soil) 260 -7240 -220 -200 -180- 85 GHz, v-pol. brightness.....Models (wet soil) 1 4 0 - I I I I 1 1 I I 1 1 1 1 1 I UI I I I I I I I I I III III I r 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.11: Modeled and observed v-pol. terrain brightness at 19, 37, and 85 GHz wi'th wet model soil and space-based observations from SSM/I.

141 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93,1,i1,,, I I III L I I I, ll I I I I I,,I l, 1, I Ill1 I. I,,,,l. v CD (/n (/) a, c 0).m a) c C0) m y 1-1.-0 C Cm 260 240 -220 -200 -180 -160 - 140 - 260 240 220 -200 -180 -160 140-r 260 240 -220 200 -180 -160 - 140 - 19 GHz, h-pol. brightness Models (dry soil)........... Models (wet soil) II 151111111111 ~111111111111111111111111 IlIllilli 111111111111111 37 GHz, h-pol. brightness Models (dry soil) —...... Models (wet soil) ri I 1 I I I I111 1 11 I i I I, I I I' I I 1 II I i, I 1 I II I I i I 85 GHz, h-pol. brightness - Models (dry soil) -.. —.. Models (wet soil) -T,ll1 II ll l,,,l I X, I * I X! ' I 1111 111,,' lii 340 350 360 370 380 390 Day from Jan. 1, 1992 400 Figure 7.12: Modeled h-pol. terrain brightness at 19, 37, and 85 GHz with both wet and dry model soil.

142 brightening effect of lower soil moisture and also the overall decreased influence of soil conditions with increasing frequency. For example, the average difference between wet and dry soil brightnesses is 15.4 K at 19 GHz but only 9.7 K at 37 GHz. The two cases are almost indistinguishable at 85 GHz except after day 393 when the snowpack has thinned and the wet soil is mostly thawed. Figures 7.13 through 7.16 present additional data on model-observation mismatches at 19 GHz. The discussion above suggested that 19 GHz emission is primarily a function of the amount and state of soil moisture. The top two graphs of Figure 7.13 demonstrate that this is in fact the case for modeled brightness, which rises and falls in synchronization with partial freezing and thawing of the model soil. REBEX-1 19 GHz brightness with 2 cm soil temperature observations are also plotted and the following argument suggests that its variation is also the result of liquid water content. On day 337 at the beginning of the plot, the soil at 2 cm is just above freezing and 19 GHz brightness is at 270 K. Compared with the IR radiometric temperature (see Appendix B), the 19 GHz emissivity was 0.97 or higher from day 324-the last time rain was recorded at the REBEX-1 site until day 344 when it began to decline. Video images of the REBEX-1 site indicate light snow cover as early as day 338 with surface temperatures remaining below freezing from day 338 to 344. From day 344 to 349, the surface temperature was often at or above freezing. These data suggest that light sub-freezing snow fell on a dry soil on day 338 and had little effect on emissivity until some melting occurred between day 344 to 349. On day 349 a more substantial snowpack formed covering and insulating the still warm soil from the cold air. And as Figure 7.13 shows, the 2 cm soil temperature rose slightly under the snowpack from day 349 through 355. Some melting at the snowpack base likely occurred, adding to the soil moisture.

143 v U) U) 4-' - n0) 260 - 240 - 220 - 200 -180 -160 - Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 1 1!f lI 11 l l l ll lJl l lI I f l ll l ll l I l l II l I l l I ll Il Il l;.~.19 GHz, h-pot. bi. rghtness. -pModels..(wet soil). i,.........Mdls(e ol i'~ "'i' 140............................................................... I CT co E Cv) E 0.4 - 0.3 - rT 1111rl 1111 11717111 1 11111717111111 11 II 1 ~I ~~ C.. ~I ~~ ~~ ri.. .. ~'; ~. I r I~ ~ ~r ~. r r~ * I ~~ 1. 3 Z I ~~ r c ~" 1R Modeled unfrozen soil water (wet soil) ~., 3 z x .., ~, ~ '.~ ~~ i ~~::......5.. '' r ~' ~' ~ I~ ~~,I.. i i:I " "''' '',, ~~ II ii\j.. ~ ~~ 4 ~~ I'*'' ~.. ~~.... r ~s r ~ ~~.. ~ ~r ~r~..:,.. 0.2 -............................. y () Q) - cn a, m v - -_ E aD Q FE 4 -C~ 260 -240 - 220 - f''l^flff''l''I''l ff'l'l'fl I flffl f 'lij 11f I''' 111111' 19 GHz, h-pol. brightness REBEX-1 200 - 180 I, I I i l I 274 272 0.2 - 0.1 -0.0 - Snow depth,.., Sioux Falls LCD' —.. ---—.... Snowflow model (wet soil), -..: 1............................................................. T I I I I I I I I I I I I II I I I I I I I fIII I I I I I I I I I I UI I I I I I I I I 1 1 340 350 360 370 380 390 Day from Jan. 1, 1992 II I 1 400 Figure 7.13: Plots of modeled 19 GHz h-pol. terrain brightness, modeled soil surface unfrozen water content, observed 19 GHz h-pol. terrain brightness, observed 2 cm soil temperature, and modeled and observed snowpack depths.

144 As the 2 cm temperature dropped below freezing around day 361, some portion of the soil moisture began to freeze. At the same time, 19 GHz brightness reached a temporary minimum at 223 K, which stood until day 368. A two-day decline in brightness began on day 368 leading to the lowest 19 GHz brightnesses around 190 K on day 370. Soil temperature at 2 cm rose to a maximum of 272.9 K on day 270 commensurate with the declining brightnesses and consistent with melting at the soil surface. If soil moisture beneath the snowpack is as dynamic a quantity as this sequence suggests then a coupled moisture and energy transfer model for freezing soil is required to achieve an accurate brightness simulation. Simple modification of the fixed soil moisture in the Snowflow model would fail to match the dynamic variations described in this scenario. A possible alternative to the soil moisture variation scenario is that snowpack stratification is the source of observed 19 GHz dynamics. The surface of a snowpack is subjected to compaction from wind and atmospheric conditions, and this effect is not modeled by Snowflow. When surfaces are buried by new-fallen snow, interfaces may develop with moderate dielectric contrasts. A series of moderately sharp contrasts can have the same darkening effect on emission as a single strong boundary. Creation of such boundaries and their subsequent disappearance through compression of the snowpack could explain the variability in REBEX-1 19 GHz brightness. Figures 7.14 and 7.15 test the enhanced layering hypothesis by comparing the effects of enhanced dielectric contrasts to those of the dry soil model. For these graphs, the Esnow snowpack emission model was modified such that the dielectric contrast at each snowpack layer boundary was increased by a factor of two by modifying the greater of the two dielectric constants. The upper limit on the modified dielectric constants was that of pure ice (3.15). Figure 7.14 shows the effect of enhanced layering

145 for h-pol. brightnesses which are reduced by about 15 K at 19 and 37 GHz. The change is minimal at 85 GHz and when the snowpack is wet because in these cases emissions originate near the top of the snowpack. The effect of layering is also small at vpol. for all frequencies because of the Brewster angle effect, with v-pol. brightnesses decreasing by 1.4, 3.1, and 3.6 K on average at 19, 37, and 85 GHz, respectively. This test shows that layering may be the cause of 19 GHz h-pol. darkening in general. But with respect to the temporal signature of 19 GHz brightness at the REBEX-1 site (and from coincident SSM/I measurements), the model with enhanced layering still fails to simulate the initial gradual decline in brightness and the strong dip around day 370. Figure 7.15 shows 19 GHz SSM/I brightnesses at v- and h-pol. with plots of brightnesses from the three model treatments: wet soil (the baseline case), dry soil, and wet soil with enhanced layering. What this data shows is that in all cases there is dynamic information missing from the modeled brightnesses. The wet soil baseline model matches v- and h-pol. SSM/I brightnesses simultaneously for a short period around day 374 and the dry soil model matches some points between day 350 and 355. These results are consistent with the hypothesis that soil moisture variation is driving brightness dynamics. And since the net effect of enhanced layering is to cool h-pol. brightnesses at all times, there is no indication that layering in the natural snowpack is responsible for increased brightness variation. Another test of layering is to examine SSM/I data for another site where temperatures are colder and the snowpack is deeper. Figure 7.16 shows the 19 GHz vand h-pol. SSM/I terrain brightness near the Langdon Experiment Station in North Dakota (48.75N latitude, 98.33W longitude). Also given are the polarization difference and snow depths from Local Climatological Data. In January, the snow depth at Langdon is up to 0.6 m deep compared to only 0.2 m at the REBEX-1 site. Yet

146 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93.,......... I.....................I..........I.... 1/30/93.1.. I... I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I c 0).-7 co m Cl) co - e. 0).C/ c C-.mo 260 - 240 - 220 - 200 - 180 - 160 - 140 - 260 - 240 - 220 - 200 - 180 - 160 - 140 - 260 - 240 - 220 - 200 - 180 - 160 - 140 - 19 GHz, h-pol. brightness........... Models (wet soil) -- Models (wet soil, layering) 'I.............................................................. 37 GHz, h-pol. brightness........... Models (wet soil) - Models (wet soil, layering) h.k *11 f a I I I I a I I I I I f I I I I I I I II - - 85 GHz, h-pol. brightness........... Models (wet soil) Models (wet soil, layering) ' I....... a I I I I III I I I I I I I I I I I I I I I I I I I........... 340 I i 350 I i 350 360 370 380 Day from Jan. 1, 1992 390 400 Figure 7.14: Modeled h-pol. terrain brightnesses with dielectric layering artificially enhanced.

147 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 11/30/93 280 260 220 19 GHz brightness 200 -v-pal. SSM/l........h-pol. SSM/l 280 -v-pot, models (wet soil). --- —-h-pol. models (wet soil) 260........ ~240 C 220-.2 200 -180 -160- hfI 111 111T11111111II 7 111jl 280 -v-pol. models (dry soil). --- —-h-pol. models (dry soil) ~260 Cn 240 *. 220 -200 -280 -v-pal, models (wet soil, layering)......h-pal. models (wet soil, layering) 260 -240 220 -180 -16 0 -IIIIIIIII IIII i III III* I 340 350 360 370 380 390 400 Day from Jan. 1, 1992 Figure 7.15: SSM/I 19 GHz v- and h-pol. brightnesses (top graph and overlays) compared to model results with wet soil, dry soil, and wet soil plus enhanced dielectric layering.

148 co............................................................ -- - - ---- --- ---- --- ---- --- ---- --- ---- --- 0 ~...~ (0 -) - - - -- 1~bz lba. t CO %. CO CC C) cm) -- - -------------------—, -- ----------- - -------— ~-C ---........ - -------- -- ~-~i~ --- —--— ~ —. — - ------ 4 0Co CZ-0 C CZ ~ Co 0)0 C~Co c' '.. C) 0 00. ~3.... I92 C~ Co I I 0~ 0' 0 0) ZZ Z), C', Co co1 CV 0 0) I, I I I. C'Je 0 0 0 0 0 0 0 0 0 0 (CDLO 'tCY) C'J 0 0 co CD "t C~l 0 co ~ 'c C\j 666 66 Co (' ' ' (4 (4 ' CY) sse CbjqC\p3i u ij (>)sCjtbjope n a N -A ()L4dpm u

149 at 20 K, the average polarization difference for the month is less than that for Sioux Falls, which was 28 K (wet snow periods excluded). In addition, the polarization difference at Langdon decreases slightly with increasing snow depth in December as the snow-soil boundary is obscured. Average air temperatures at Langdon in December and January were well below freezing (257 and 255 K, respectively). If we assume that such cold temperatures imply very low liquid water content in the soil, then these data suggest that layering contributes some to decreased h-pol. brightness but that soil conditions are the more dynamic driver. 7.5 Comments on snowpack structure, grain size, and wetness The strongest temporal snowpack radiobrightness signature is caused by partial snowmelt. Figure 7.17 plots the correspondence between model snowpack liquid water content and modeled and observed 37 GHz h-pol. brightness. Note that modeled brightness spikes on days 390, 393, and 394 each produce similar maximum brightness temperatures-that is, once a threshold amount of moisture is present the brightness temperature saturates at about 265 K. The four arrows point to modeled snowmelt events that correspond to observed brightness jumps. Where model and observation are in disagreement, the contrast may be used to adjust parameters of the snowpack SVAT model that are otherwise fixed. For example, the snowmelt event on day 390 reaches a maximum wetness of 1.84 mm. A 25% reduction in the absorbed shortwave radiation-representing, for example, an increase in the snowpack shortwave albedoreduces maximum wetness by 50% to 0.9 mm. Of course, more detailed study is needed to determine the appropriate feedback parameters and gains but this example demonstrates the process.

150 Date 12/11/92 12/21/92 12/31/92 1/10/93 1/20/93 1/30/93 Y(I) CD c.mt__ 260 - 240 - 220 - 200 - 180 -160 - 37 GHz, h-pol. brightness --- REBEX-1 140 m 1 1 11 1[ I I I I 1 I I I I I I I I I 1I11 I I I 1 I I 1I11 1 I 11 1 I I I I I I I1I I1I I I I I 1 I I 1 1I II1'1' 1''1''1'1'1IIaIIAIIIIIIIIIII I....A II'' I''1' 1''1, I. III I IIII I II ~ (I),) Co -C 260 - 240 - 220 -200 - 180 -160 - 37 GHz, h-pol. brightness..... — Models (wet soil) E E -' Co Q.-. LU 140 3.0 2.0 1.0 0.0..-.... - - - - - - - - - - - - - -. - - - - - - - -...... I 11111 ""1111111 it.jiljiIliil 1i111111111111 11 jiIl 111Elj Snowflow modeled snowpack wetness -- 340 350 360 370 380 390 Day from Jan. 1, 1992 400 Figure 7.17: Modeled and observed h-pol. 37 GHz terrain brightnesses and modeled total snowpack liquid water content. The arrows indicate modeled partial snowmelt events that correspond to observed brightness jumps.

151 E 0 - cn.C.0) c ET.I - E 3: 0 c cca) I 0.20 - Day 392.817 0.15 -- \ 0.10 0.05 0.00- I ' 0.0 1.0 2.0 0.20 Day 393.003 0.15 -0.10 0.05 0.00- - I - I 0.0 1.0 2.0 0.20 X Day 393.189 0.15 0.10 0.05 0.00 I - 0.0 1.0 2.0 Day 392.879 0.0 1.0 2.0 Day 393.065 — i- | —i — 0.0 1.0 2.0 Day 393.251 - -/ ' j I I I 0.0 1.0 2.0 Day 392.941 - 0.0 1.0 2.0 Day 393.127 0.0 1.0 2.0 Day 393.313 I I I I I 0.0 1.0 2.0 Snow liquid water content (m /m x 100) Figure 7.18: Snow liquid water content profiles during a cycle of partial melt and refreeze. The cycle of snowpack wetting is shown in detail in Figure 7.18 which chronicles partial melt and refreeze on day 392-393.1 First melt occurs at 12:15 CST and begins near but not at the surface of the snowpack, which is 0.18 m deep at this time. The maximum wetness occurs at 17:00 CST, and the snowpack is completely refrozen by 02:00 CST on day 393. There are two implications of this data. The first is that satellite measurements-which for a sun-synchronous orbit have a nominal 12 hour spacing-need to be synchronized such that one of the daily measurements occurs during the coldest part of the day-ideally, 04:00-06:00 local time. The second 1Fractional day is in Universal Time (UTC) which is 6 hours ahead of local time at Sioux Falls (CST). 393.0 fractional day corresponds to 18:00 local time on day 392.

152 0.20 - E 0.20 Day 355 Day 365 Day 375 0.15 0 |" 010 s 0.05 '0.00 I o.! o- - i- I I I! I I I I 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 Snow grain diameter (mm) 0.20 Day 385 " Day 395 | 0.15 0 0.10 r 0.05 ' oo —\ -- - ' — \-v I 0.00 I I 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 Snow grain diameter (mm) Figure 7.19: Snow grain diameter profiles from the Snowflow model. implication concerns the shape of the melt profile. If the melt region extends to and wets the soil surface, then residual moisture may be left at the soil surface, which is separated from the cold air and sky by the insulating snowpack. Consequently, when the snow completely refreezes but the soil stays wet, the brightness after the meltrefreeze cycle will be colder than before at frequencies which penetrate the snowpack. Another snowpack property affected by partial melt-refreeze cycles is the grain size. As discussed in section 3.4.7, grain growth is accelerated by the presence of liquid water. Figure 7.19 shows example grain size profiles in 10-day intervals through the test period. On day 375 the grain size at all levels is about 1.0 mm or greater. At the same time, both model and REBEX-1 85 GHz brightnesses are at their lowest point to date-around 180 K-indicating strong scattering (see Figure 7.10). Even lower brightnesses (155 K) were observed around day 395 after several melt-freeze cycles. Model and observation diverge at this point, with the 85 GHz model about 30 K

153 warmer than observations. A probable cause is that model grain size did not develop as quickly as grains in the natural snowpack. In contrast, at 37 GHz (Figure 7.9 the model matches observations to within 5 K during the re-freeze portion of the cycle on days 395-396 while the 19 GHz model is about 40 K too cold. Since the model soil is almost completely thawed here (see Figure 7.7), 19 GHz model darkening is directly attributable to the high dielectric contrast. Then the 37 GHz darkening is probably due to separate effects: larger grains in the natural snowpack and wet soil in the model. Figures 7.20 through 7.22 show pairs of SSM/I brightness images that demonstrate how the melt-freeze process correlates over spatial scales. The data are from days 401 and 405 near the end of the period of melt-freeze cycling which lasted fromnt day 397 to 407.2 The day 401 image was acquired at about 17:24 local time and the day 405 image at about 06:44 local time and probably correspond with times of high snowpack wetness and refreeze, respectively, at the REBEX-1 site. The strongest melt-freeze signal can be seen at 85 GHz in Figure 7.22. In the bottom image, a darkened region extends from about 43N to 45N latitude and 95W to 97W including the REBEX-1 site at 43.72N, 96.5W. In the top image the same region is of comparable brightness to the surrounding terrain. The consistently brighter region to the west may have a wet snow cover that extends over into the REBEX-1 region during daylight on day 401. 19 and 37 GHz images in Figures 7.20 and 7.21 show the same pattern but with the contrast between wet and frozen snow decreasing with frequency. This decrease is consistent with the REBEX-1 observations of melt-freeze cycles. By day 405 the REBEX-1 site had undergone several melt-freeze cycles with commensurate opportunities for accelerated grain growth. Consequently, it is likely that the REBEX-1 region is 2See the figures in Appendix B for REBEX-1 data extending past day 401.

154 radiometrically darker at 85 GHz than other regions that did not undergo melt-freeze cycles.

155 -1,13 19 GHz Horizontal Pol. (Atm. Removed) 1992:401:23:24 UTC -1.05 -97 1 pixel - 03 deg, 160 180 200 220 240 260 280 300 Rodiobrightness (K) 19 GHz -105 Horizontal Pol. (Atm. Removed) 1992:405:12:44 UTC -97 -89.38 / 1 pixel = 0.3 deg. 160 180 200 220 240 260 280 300.................... Rodiobrightness (K) Figure 7.20: SSM/I 19 GHz h-pol. images from February 4. 1992 at 23:24 1UT( (top) and February S. 1992 at 12:44 UTC. See text for description.

156 37 GHz Horizontal Pol. (Atm. Removed) 1992:401:23:24 UTC -1.05 -97 -t13 -89 t- t 1 pixel = 0.3 deg. 160 180 200 220 240 260 280 300 Rodiobrightness (K) -1.13 37 GHz Horizontal Pol. (Atm. Removed) 1992:405:12:44 UTC -1.05 -97 -89 /......... 1 pixel. 0.3 deg. 160 180 200 220 240 260 280 300 Rodiobrightness (K) Figure 7.21: SSM/I 37 GHz h-pol. images from February 4. 1992 at 23:24 UTCr (top) and February S. 1992 at 12:44 ITTC. See text for description.

157 85 GHz Horizontal Pol. (Atm. Removed), High Res. 199 gg2:401:23:24 UTC -1.05 -97 -113 -89 1 pixel = 0.2 deg. 160 180 200 220 240 260 280 300 Rodiobrightness (K) -1.13 85 GHz Horizontal Pol. (Atm. Removed), High Res. 1992:405:12:44 UTC -1.05 -97: -89 -- 1 pixel = 0.2 deg. 38...... 160 180 200 220 240 260 280 300............. Rodiobringhtness (K) Figure 7.22: SSM/I 85 (Hz h-pol. images from February 4., 1992 at 23:24 and February 8, 1992 at 12:44 UTC. See text for description. IT'TC (tojp)

CHAPTER 8 CONCLUSIONS, CONTRIBUTIONS, AND RECOMMENDATIONS This chapter summarizes the major results and implications of the radiobrightness data analysis in Chapters 6 and 7. The chapter also lists the contributions of the thesis. discusses some lessons learned from experimental work. and makes reconimmendations for further study. 8.1 Conclusions The major results of this thesis come from the comparisons of simulated and observed snowpack radiobrightnesses in Chapter 7. The principle implications of these comparisons include: (a) 19 and 37 GHz radiobrightnesses are sensitive to the d( ynamically-varving unfrozen soil moisture beneath the snowpack, (b) snowpack stratification may reduce these brightnesses but is not as dynamic a, driver as soil wetness, (c) dry snow 85 GHz brightness is sensitive to grain size development during snowpack partial melt and sky- brightness is the main driver of its variability, and (d) the sharp contrast between the brightness of dry and partially melted snowpacks may be explained by a lossy-cloud model of wet snowpack emission. The careful treatment of grass may be critical to producing a simulation that accurately models energy exchange between the snowpack and ground-especially regarding shortwave 158

159 radiation penetration. Soil wetness modeling requires the handling of melt-.water from the snowpack and a coupled moisture-heat transfer model for freezing soil. The measurement of microwave brightness by ground-based instruments has been shown to be a viable way of simulating long-term space-borne monitoring, although the large field-of-view of satellite instruments must be taken into account. Comparisons between ground-based and satellite radiobrightnesses indicate that brightnesses correlate better when a snowpack is present than at other times because of the relative homogeneity of regional snowfall compared to vegetation cover. There is less agreement between ground- and space-based measurements during episodes of snowpack partial melt and refreeze because of the magnitude of the wet-dry radiobrightness difference and the fact that the process is not usually synchronized across the SSM/I field-of-view. Also, it is likely that variability in atmospheric water vapor content and temperature at times adds to the difference between ground- and space-based measurements-especially at 85 GHz. Yet with only simple atmospheric corrections of SSM/I brightnesses, the differences between ground- and space-based measurements were within about 10% with respect to the January snowpack's radiometric dynamic range at all frequencies. The use of REBEX-1 data with models suggests several ways in which future experiments of this type may be improved.1 REBEX-1 data were used to drive the snowpack SVAT in the same way that an atmospheric model would drive its landatmosphere boundary condition parameterization. The key parameters missing in this context were downwelling longwave flux and precipitation from a heated gauge. As for model diagnostics, the most critical missing parameter was unfrozen soil water content, which may be measured automatically through time domain reflectometry 'Some of these lessons have already been incorporated in implementation of REBEX-3 in Alaska, which was begun in August, 1994.

160 techniques. Other measurements that would have been useful-and are possible with existing methods-include soil surface temperature under the snowpack, temperatures inside the snowpack, albedo (shortwave reflectivity), and Bowen ratio (the ratio of sensible to latent heat transfer rates) for closure of the land-atmosphere energy balance. It would also aid modeling if one of the key sources of measurement errorfrost and snow obstructing radiation instrument domes-could be eliminated. Also, Chapter 5 spends a great deal of time correcting REBEX-1 sky br-ightnesses for the effects of the reflector used in their measurement. The accuracy of these measurements could be improved if the brightnesses could be measured directly without a reflector. The results presented in this thesis suggest that a SVAT model simulating radiobrightness may be a useful diagnostic tool in atmospheric modeling. Consider an atmospheric model with a spatial resolution comparable to that of the EASE-Grid radiobrightnesses and run for a historical test period like that of Chapter 7-that is, December, 1992 and January, 1993. On a grid-point by grid-point basis, modeled brightnesses from the atmosphere-SVAT-emission model could be evaluated against SSM/I measurements of the actual terrain, and differences could be used to verify aspects of the model's performance. For example, the distinct temporal radiobrightness signal of melt-freeze cycles are a sensitive indicator of the amount of energy incident on and absorbed by the snowpack. Differences in the timing of melt-freeze episodes between model and measurements could be used to judge the model's net land-atmosphere energy exchange rate. Another possible feedback parameter is snow grain size. Divergence between modeled and observed 85 GHz brightness would track grain size development which is closely tied to snowpack temperature gradients and internal heat flux in dry snow. Lower modeled brightnesses would indicate larger

161 grain sizes and greater flux rates than in the historical record. Of course, the model relating grain size development to temperature gradients would need improvement over the one presented here to make this difference meaningful. Lastly, since the wetness of frozen soil is closely tied to the deposition of energy in the snow and soil, 19 GHz darkening due to increased soil wetness might also be a sensitive indicator of energy exchange rate differences. The link between radiative balance and snowpack state may be exploited by a well-tuned snowpack SVAT-emission model to evaluate the shortwave and longwave fluxes from an atmospheric model. A snowpack SVAT's parameterizations of sensible and latent heat flux and albedo can be tuned through extended measurements of the complete snowpack-atmosphere energy balance. With a well-tuned snowpack SVAT, radiobrightness would be most sensitive to radiative forcing from a model atmosphere. As described above, 19 GHz brightness is sensitive to shortwave radiation penetrating to the soil and 37 GHz brightness is sensitive to net radiation through snowpack melt-refreeze cycles. Both these processes are tied to time-averaged conditions so once or twice daily satellite measurements would be adequate samples. Since radiative forcing from the atmosphere is a direct function of the percent cover and height of clouds, radiometric monitoring of snowpacks may be a good way to determine the accuracy of the cloud component of an atmospheric model. This is demonstrated by the example in section 7.5 that showed that-on a day when the snowpack melted partially an artificial 25% reduction in the total diurnal shortwave flux reduced the maximum snowpack wetness by 50%. Consequently, 37 GHz brightness would be sensitive to cloud cover and could be used to find the percent cover range corresponding to observed snowpack conditions. Linking models of atmosphere, snowpack, and radiobrightness may also present

162 a way to improve the monitoring of snowpack parameters from space. Snowpack equivalent water content and onset of melt are the primary parameters of interest in snowpack monitoring. With the aid of an atmospheric model that simulates the current state of the atmosphere (often called a "nowcast"), a snowpack SVAT could be used to infer snowpack history and soil conditions. Since a particular radiobrightness signature may be indicative of multiple snowpack states, a priori information from a measurement-adjusted model can help eliminate possible interpretations of new data. For example, increasing snowpack depth and grain size both lead to darkening of 37 GHz brightnesses. If snowpack melt and refreeze occur, then grain size will increase while snow depth decreases. If the model accurately simulates the meltrefreeze cycle, then the darkened brightnesses would not be mistaken for a signal of increased snow depth. 8.2 Contributions This thesis presented a record of concurrent microwave radiometric and meteorological measurements at a northern Great Plains site (REBEX-1) and developed a SVAT-linked microwave radiobrightness model that simulated the observed snowpack. The major contribution was the REBEX-1 data set which includes 192 days worth of continuous h-pol. radiobrightness measurements at 19, 37, and 85 GHz. The snowpack SVAT and radiobrightness models are implementations of existing knowledge but their combination and use with the unique REBEX-1 record represents a novel approach in the development of radiometric models and their testing. The thesis also demonstrates the resampling of a long-term SSM/I brightness record to a common resolution at a single earth surface coordinate. Successful comparisons of this data with terrain radiobrightnesses from REBEX-1 show that the resampling technique is

163 sound and useful and that ground-based radiometric measurements are meaningful surrogates for high-temporal space-based observations. 8.3 Recommendations The successful implementation of REBEX-1 suggests the usefulness of further long-term ground-based studies of the link between radiobrightness and antecedent weather. The analysis here covered only a portion of the REBEX-1 data set. Further analysis of the data concerning the radiobrightness of grass-covered terrain in fall (lasting about 60 days) and during spring thaw (about 46 days) is suggested. Also, a repetition of the experiment with improvements to the instrumentation recommended above would be valuable in resolving some of the ambiguities in the REBEX-1 data. A more northerly site may be more representative of the Great Plains as a whole. Also, the site's vegetation cover should be consistent with that of the surrounding terrain to facilitate comparisons to space-based instruments. The SVAT-linked emission model may be improved by the addition of a coupled heat and moisture transfer model of the soil, although further empirical work on the relationship between temperature below freezing and the unfrozen water content of organic soils needs to be done. Current research on the growth of snow grains and a parameterization of new snow density accounting for wind packing should also be incorporated. The total albedo and distribution of shortwave radiation il shallow snowpacks with grass need to be further examined although improvements to the snowpack SVAT may be made immediately by increasing the number of shortwave absorption bands. Lastly, the empirical correction to independent scattering discussed in Chapter 4 needs to be re-examined in the context of finding a comprehensive but computationally efficient solution to the problem of microwave radiative transfer and

164 scattering in snowpacks. The usefulness of dynamic SVAT-linked radiobrightness models can only be exploited if satellite observations are frequent enough to capture the events of interest. The minimum useful measurement periodicity for capturing events like partial melt and refreezing of the snowpack is twelve hours. Each additional measurement could be used to determine the timing of the melt event. The coverage of each SSM/I instrument is incomplete at the latitudes of the Great Plains and two sun-synchronous satellites are required to guarantee diurnal coverage. Based on the results in this thesis, the preferred timing of the first two satellites would be at 04:00-06:00 local equator crossing time. The maximum benefit from additional satellites could be gained with a 12:00-14:00 local equator crossing time. In addition, a radiometer channel at 10 GHz would be useful as a probe of soil surface conditions since radiation at this frequency would interact less with the dry snowpack and be more more responsive to soil water.

APPENDICES 165

APPENDIX A SNOWFLOW PARAMETERIZATIONS This appendix gives formulas for variables used in Snowflow (Chapter 3). Section A.6 lists values for the parameters used in these formulas and those in Chapter 3. A.1 Saturation vapor pressure and density The saturation vapor pressure, Pws, is found by solving the Clausius-Clapeyron differential equation: d In pw lk(T) dT RwT2 where lk is either ly in the case of vaporization or Is for sublimation. Colbeck [61] suggests the following solution for use in snow modeling: ls(T-To) pwsoe RWTTo T < T, Pws(T) = (A.2) I l (T-To) pwoe R TTo T > To. where Pwso is empirically determined at T = To. Then the saturation vapor density is: f(T) vs = RP- (A.3) RT' 166

167 Snowflow requires the first, second, and third derivatives of f(T) and these are given by: f'(T)- Pw [TS -k i] (A.4) f (T)= R 2 R +T ( RT, (A.5) RwT3 RwT RwI p s P(T) 181k 31 2 k 3( f"'(T) 4- T -6 + R ) ] (A.6) wTT4 lwT wT. wT A.2 Wet bulb temperature The wet bulb temperature, Twb, is given by recursive evaluation of the relationship: wb = Tair IkPws (Twb) R Twb = Tar - l k-(TweR(1 - RH) (A.7) Cp,airP where Ik is the latent heat of vaporization 1, or sublimation Is: pS(T) is the vapor pressure of saturated air at temperature T (from (A.2), ER is the ratio of the dry air and water vapor gas constants (Rd/Rw), cp,air is the specific heat capacity of air, p is the air pressure, and RH is the relative humidity. A.3 Water vapor diffusion coefficient The effective water vapor diffusion coefficient is given by: De(T) = DeoP (1) T P T"o/ = CDeTD (A.8) where po is atmospheric pressure at sea-level and CDe and nD are constants. A.4 Composite conductivity factor, Fl The composite conductivity factor is: F1(T) = K + lkDe(T)f'(T) = K + lkCDeT D f'(T) (A.9)

168 where K is thermal conductivity and lk is either lv or I,. A.5 Freezing of water flowing into dry snow Water at T = TO that flows into a dry layer will freeze until an equilibrium is reached. Snowflow finds the amount of freezing by first assuming the equilibrium temperature is TO. Then the water content of the equilibrium state would be (in m): q = o (dst ((CH + H2To)To - (CH1 + T)T If P \ 2 2 (A.10) where Wo is the depth of water before freezing occurs, (d,p,)t is the mass of ice in the layer, T is the temperature of the layer before the water was added, and CH1 and CH2 are coefficients of the temperature dependent heat capacity of ice, Cice(T) = CH1 + CH2T. If Weq is less than zero, then all the water that seeped in eventually froze and the equilibrium temperature of the layer is found from: Te CV I + 2CH2Qfreeze - CH1 (A.1) Te - (A.I1) CH2 where Wopw (lf + (CH1 + T)) + (dss)t(CH1 + 2rT)T Qfreeze (d p) + Wo (A (A.12) A.6 Summary of Snowflow parameters The following tables list Snowflow parameter values that are (a) set by the user, (b) empirically determined but fixed, and (c) fundamental physical constants. All parameter values are given in SI units.

169 Table A.1: Snowflow parameters set by the user. Symbol Description Value A snowpack albedo (range 0.4-0.95 [62]) 0.7 dgo top soil layer thickness 0.001 m ds,opt optimal initial snow layer thickness 0.02 m p air pressure 1.013E5 N/m2 Qg heat from below the bottom soil layer 0 J/m2 Wm weight fraction of non-quartz matter in dry soil 0.75 Wq weight fraction of quartz in dry soil 0.25 Xw initial soil water volume fraction 0.20 or 0.43 Table A.2: Snowflow empirical parameters and their sources. Symbol Description Value Source _ _ _ Cpav average dry soil heat capacity C5 snow compaction coefficient c6 snow compaction coefficient Cgsl snow grain growth coefficient Cgs2 snow grain growth coefficient CH1 ice heat capacity coefficient CH2 ice heat capacity coefficient Deog water vapor diffusion coefficient in ground Deos water vapor diffusion coefficient in snow d zero displacement height for snow eIW longwave emissivity of snow ga,p size factor of soil solids ka von Karman constant nD,g water vapor diffusion temperature exponent in ground nD,s water vapor diffusion temperature exponent in snow Sr irreducible water saturation constant xads volume fraction of adsorbed soil water Xfld volume fraction of water at field capacity xv volume fraction of soil voids Zo roughness length for snow Z2 screen height (air temperature) z3 wind speed height 900 J/kgK 0.08 1/K 0.021 m3/kg 5.0E-7 m4/kg 4.0E-12 m2/kg 92.95 J/kgK 7.369 J/kgK2 x,-1.61E-5 m2/s 0.92E-4 m2/s 0 m 0.97 0.125 0.4 2.3 6 0.04 0.07 0.1 0.45 2.0E-3 m 1.8 m 10 m [39] [38] [38] [38] [38] [37] [37] [38] [38] [46] [38] [42] [46] [38] [38] [38] [42] [42] REBEX-1 [46] REBEX-1 REBEX-1 continued on next page

170 continued from previous page Symbol Description Value Source acu unfrozen soil moisture factor 0.238 [41] 3Wu unfrozen soil moisture factor -0.360 [41] r/o snow viscosity coefficient 3.6E6 Ns/m2 [38] Aair,d thermal conductivity of dry air at To 0.0241 W/mK [63] Ai thermal conductivity of ice 2.18 W/mK [42] Am mean thermal conductivity of 2.93 W/mK [42] minerals and organics Aorganic thermal conductivity of organic matter 0.25 W/mK [42] Aq thermal conductivity of quartz 8.16 W/mK [63] XA thermal conductivity of water 0.561 W/mK [63] pb soil bulk density 972 kg/m3 REBEX-1 pq intrinsic density of quartz 2.66E3 kg/m3 [42] pm mean intrinsic density of 2.8E3 kg/m3 [42] non-quartz soil minerals vo snow shortwave attenuation factor 0.084 [37] vi ice shortwave attenuation 1.72E3 m [37] Table A.3: Other physical parameters used in Snowflow. Symbol Description Value Cp,air Cpi Cpw g if Is Iv PwsO Rd Rw To 7d Pi Pw Pw ar specific heat capacity of dry air (at 273 K) specific heat capacity of ice (at 273 K) specific heat capacity of water (at 273 K) acceleration of gravity latent heat of fusion (at 273 K) latent heat of sublimation (at 273 K) latent heat of vaporization (at 273 K) saturation water vapor pressure at 273 K dry air gas constant water vapor gas constant freezing point of water dry adiabatic lapse rate dynamic viscosity of water density of air at sea level intrinsic density of ice intrinsic density of water Stefan-Boltzmann constant 1005 J/kgK 2106 J/kgK 4218 J/kgK 9.8 m/s2 0.334E6 J/kg 2.834E6 J/kg 2.501E6 J/kg 610.5 N/m2 287.05 J/kgK 461.51 J/kgK 273.15 K 9.8E-3 K/m 1.78E-3 Ns/m2 1.225 kg/m3 917 kg/m3 1000 kg/m3 5.6696E-8 W/m2K4

171 APPENDIX B DATA FROM REBEX-1 B.1 Microwave radiometer performance evaluation To test the calibration accuracy and radiometric resolution of the microwave radiometers, I placed a known source (microwave absorber) in front of the radiometer during regular experiment cycles. The use of the complete experiment cycle for these tests was necessary to simulate experiment conditions including internal load TPRmode gain factor recalibration. Table B.1 gives the results of these tests. Here, radiometric resolution is the standard deviation of the errors. I did not include the day 403 1448 and 1504 UT calibration checks in the radiometric resolution calculations because the Eccosorb temperature was more than 3 K over the ambient air temperature at that time. The discrepancy suggests that there was a temperature gradient between the Eccosorb surface and the thermistor embedded in it. For these times, Table B.1 compares the 19 GHz apparent brightness to the Eccosorb thermistor temperature and the 37 and 85 GHz brightnesses to air temperature as approximations to actual emitting temperatures. In the future, the Eccosorb should be shielded from wind and direct sunlight and should not be brought into direct contact with the surface of the heated radiometers.

172 - - -o-III 19 GHz 37 GHz 85 GHz Day Time TEcco TAP A TECCO TAP A TECCO TAP 279 1813 296.7 296.9 0.2 296.7 297.0 0.3 296.7 295.1 -1.1 296.7 296.3 -0.4__ 1823 297.4 297.6 0.2 297.4 297.6 0.2_ _297.4 296.8 -0.6 297.4 297.3 -0.1 288 1743 279.9 281.2 1.3 279.9 278.0 -1.9 279.8 280.0 0.2 279.8 277.9 -1.9 1753 279.9 280.2 0.3 279.9 278.6 -1.3 279.8 280.2 0.4 1854 281.8 281.8 0.0 281.8 280.6 -1.2 281.2 280.3 -0.9 281.2 279.9 -1.3 1903 80 80.8 0.8 1913 281.1 281.7 0.6 281.1 280.5 -0.6_ 281.5 281.7 0.2 281.5 281.2 -0.3 289 1644 278.1 278.5 0.4 278.1 276.8 -1.3_ 278.3 278.5 0.2 278.3 276.4 -1.9 1658 277.6 277 -0.6___ 1703 277.9 279.3 1.4 277.9 277.2 -0.7 277.9 278.3 0.4 277.9 277.3 -0.6 403 0238 273.4 273.5 0.1 273.4 273.7 0.3 273.4 272.9 -0.5 273.7 274.5 0.8 273.7 273.9 0.2 273.7 273.4 -0.3 0243 273.3 273.3 0 _273.3 272.5 -0.8 273.6 274.3 0.7 273.6 273.8 0.2 1448 272.3 271.7 -0.61 269.82 266.2 -3.61 269.82 271.3 1.51 272.7 273.3 0.61 269.82 269.0 -0.8 269.82 271.4 1.61 1504 273.0 271.2 -1.81 269.82 267.7 -2.11 269.82 271.4 2.61 _ 273.2 273.0 -0.21 269.82 268.7 -1.11 269.8 268.6 -1.21 Average A 0.24 Radiometric resolution 0.61 -0.61 0.82 -0.53. 1These data were not included in radiometric resolution calculations as noted in the text. 2Air temperature used in place of measured Eccosorb temperature. Table B.1: Calibration check results. TECCO is the temperature of the Eccosorb target at the time the TAP were measured. There is insufficient data to calculate 85 GHz radiometric resolution.

173 B.2 Soil and canopy moisture samples I took soil samples from the REBEX-1 site on five days in October and November, 1992, using an 8.255 cm (3.25 in) diameter cylindrical coring tool. A 6.35 cm (2.5 in) high bit screwed on to the end of the coring tool, holding in place a set of rings which lined the cylinder. From the top of the bit, we inserted rings of 2, 2, 2, 8, 2, and 8 cm heights. To take samples, we drove the coring tool into the soil to a depth of 12.35 cm. We twisted the tool to break the soil at the bottom of the bit and extracted the core. We then sliced off the soil extending from the bottom of the bit, unscrewed the bit from the cylinder, and removed the core from the cylinder surrounded by the bit and the three lowest 2 cm rings. We sliced between the rings and between the bottom ring and the bit with a knife to make three samples of 2 cm thickness. We then slid the soil remaining in the bit into a 4 cm ring and removed the excess from the bottom. The soil samples then spanned 0-2, 2-4, 4-6, and 6-10 cm depths. We usually took grass samples before coring by cutting the grass over the sampling position as close to the soil horizon as possible. I first dried the grass and most of the soil samples in a 70~C oven until the mass had reached a constant value. The drying time was 9 days for the grass samples and 3-3.5 days for the soil samples. I then dried all of the soil samples at 105~C for an additional 24 hours, as is the usual practice. The gravimetric calculation of soil and grass moisture as a mixing ratio to dry matter is: w = Mdry (B.1) Mdry - Mtare where w stands for gravimetric moisture content, Mwet is the measured mass of the wet sample, Mdry is the measured mass of the dried sample, and Mtare is the container (tare) mass. To calculate volumetric moisture content, I first estimated the soil bulk

174 density by taking the maximum ratio of dry soil mass to core volume. By this method the soil bulk density was 0.972 gm/cm3, which is within the wide range of possible organic soil densities. Volumetric moisture content, 0, is then given as: W0 wb (B.2) Pw where Pb is the soil bulk density, and Pw is the density of water. The soil and grass data are given in Table B.2. B.3 Data from the National Weather Service To augment data obtained by TMRS and MMS during REBEX-1, I acquired two data sets from the National Weather Service: the October, 1992 through April, 1993 monthly Local Climatological Data (LCD) for Sioux Falls, SD and selected Radiosonde/Rawinsonde Observations and the LCD's from Huron, SD. The LCD's include daily summaries of snow depth, which could only be roughly estimated from REBEX-1 video stills, and precipitation, which was unavailable from the REBEX-1 rain gage during cold weather. The LCD's also give reports of sky cover, ceiling, and weather at six times per day. These data were used in conjunction with the Huron rawinsondes to determine when the sky was cloud free and to calculate the clear-sky downwelling radiobrightnesses discussed in Sections 5.5.4 and 5.5.5. The rawinsonde data set includes atmospheric profiles of pressure and temperature from ground level to about the 10 mbar level and dew point and relative humidity to the height at which dew point falls to about -50~C.

175 Day Time Depth Tare Wet Dry Dry Gravimetric Vol. Mass Mass Mass Mass Moisture Mois. 70~C 105~C 700C 105~C 105~C 287 2100 Grass 12.26 22.17 18.12 - 0.69 0-2 cm 8.99 152.12 115.75 113.01 0.341 0.376 0.365 2-4 cm 9.00 110.85 85.75 84.69 0.327 0.346 0.336 4-6 cm 8.99 100.67 78.77 77.89 0.314 0.331 0.322 6-10 cm 9.04 260.40 204.55 201.81 0.286 0.304 0.295 288 1400 Grass 12.24 20.59 17.88 0.48 0-2 cm 9.02 126.17 96.90 95.65 0.333 0.352 0.342 2-4 cm 9.00 122.15 95.07 93.89 0.315 0.333 0.324 4-6 cm 9.03 124.04 97.36 96.25 0.302 0.319 0.310 6-10 cm 9.03 237.95 186.79 184.42 0.288 0.305 0.296 2100 Grass 12.26 26.91 22.06 - 0.49 0-2 cm 9.00 141.10 108.73 107.27 0.325 0.344 0.334 2-4 cm 8.98 127.36 98.42 97.15 0.324 0.343 0.333 4-6 cm 9.02 119.24 92.39 91.34 0.322 0.339 0.330 6-10 cm 9.00 223.16 172.34 171.02 0.311 0.322 0.313 289 1400 Grass 12.23 35.72 25.64 - 0.75 0-2 cm 8.97 117.81 90.37 89.15 0.337 0.357 0.347 2-4 cm 9.03 121.19 93.93 92.59 0.321 0.342 0.332 4-6 cm 9.01 130.73 101.51 100.32 0.316 0.333 0.324 6-10 cm 9.05 234.07 182.88 180.42 0.294 0.313 0.304 2100 Grass 12.29 25.84 20.68 0.62 - - 0-2 cm 9.06 138.99 108.36 106.72 0.308 0.330 0.321 2-4 cm 9.02 126.34 98.38 97.00 0.313 0.333 0.324 4-6 cm 9.08 129.38 101.55 100.21 0.301 0.320 0.311 6-10 cm 8.99 239.31 188.64 185.78 0.282 0.303 0.295 290 1400 Grass 12.27 32.30 24.89 - 0.59 - 0-2 cm 9.08 101.61 81.06 79.95 0.285 0.306 0.297 2-4 cm 9.07 120.95 95.31 93.92 0.297 0.319 0.310 4-6 cm 9.02 123.04 97.71 96.48 0.286 0.304 0.295 6-10 cm 8.99 227.72 181.01 178.08 0.272 0.294 0.286 307 2100 0-2 cm 7.25 125.74 89.38 - 0.443 0.431 2-4 cm 7.20 127.98 - 94.87 - 0.378 0.367 4-6 cm 7.22 133.84 - 100.37 - 0.359 0.349 6-10 cm 7.20 233.01 - 176.43 - 0.334 0.325 Table B.2: Soil moisture sampling data. Masses are in grams and include tare mass. Drying was first done at 700C then at 1050C. Gravimetric moisture content is tabulated at both temperatures and volumetric moisture content is tabulated only at 105lC.

176 B.4 Overview plots of continuously measured parameters The graphs in this section summarize the measured parameters as well as several derived parameters from the entire period of REBEX-1. The parameters are 19, 37, and 85 GHz terrain apparent radiobrightnesses, reflector-measured sky radiobrightnesses, estimated sky radiobrightnesses (see the main text for a description of the estimation process), thermal infrared sky temperature, global radiation, net radiation, wind speed at 10 m, soil temperatures, heat flux in the soil, relative humidity, rainfall, thermal infrared surface temperature, air temperature, terrain radiobrightness spectral gradients, TPR gain factors, estimated 85 GHz TPR gain factor (for the uncalibrated period discussed in the main text), and snow depths. These figures contain terrain and sky radiobrightnesses edited for flawed data, with points removed that were either out of range or taken when the internal radiometer temperatures drifted from their design values. Sky radiobrightnesses and IR sky temperatures that were in error due to incorrect reflector positioning have also been removed. Chapter 5 describes in detail reflector position errors and the editing methodology applied to the sky data. Because of a communications error described in Section 5.4, the IR radiometer sky data also contains additional incorrect points but these have not been removed here. The terrain radiobrightness spectral gradients, S, shown in Figure B.8 are given by: Sg (fl, f2) TTER(f2) - TTER(fl) (B.3) f2 - fil where fi and f2 are radiometer frequencies and fi is less than f2.

177 ~~~~~...................~ ~~i:%,i; C o....................................................................................... ~,.......................-:.......... -......'"E" '.:.:.7.;.'.... '"'-.....,..,,...:;:::L;~:::~;i:;:,,~L;L;L;..:.; ---;:;?~,;,~..'.'. ~............................................................................................................................:......~.....~?i:i-~,,,~:.........................................................,______.:.:.:...._........ 0+,................................................... ~ -I _'^__________^^'.^^f-'...... -............................................................................................................................. 0..... -^ —..-. —.. —. —..-. ---................... O:- ^'^ ICOC\ -.......................................... ~.. _..........................................;.................... C) C) l...............................~'.'"-',.:-.-::......... -— ~.... 0 ~~~~ 0 CO >'.7........................................ - - - - - - - -- -R...........................................................,............._. ".':-.:.."".'.......':......................................... CZ III).....................(0......................................................................;i "...'.................................................... O0r 0 CO) NNN.................................. co 0 00 0D 0 0 0 0o 0D 0t 0... (0... C\. 0 (0 (0 CV) ~" (>1) sseulL4Buqoipej lueJiddV

178 -- - - -- - - - -- - - -..._...... - - - - - - - - - -- - - -- - - -................................... V 1111 IjIjIII. II 3111 o- 0 --- —-- o......... 0... ----------.......................................... 7.. 7. '. T.. .-..!:. '. '.- 7. 7 7 1................. NAw.................................................................................................................. NOW....................................................................................................................... - ------------ - --- - ------............................................................................................. -------------............ wig-.............................. 0 0 -N 10 CMj Cl 0+ -0 0 II o D 0) 03) 11) E 0 >1O 0 C co0 0 N1 04-. 0 -At b0 17 -0t 0d 2i......... -....... -- ------- - ------------------- i".-M i q A I I I I I I I I 1 1 1 1 1 1 1................:.!.:. - -.=... - "I!.,=,, =M-................................... (5 M00 a) L r-rl IrII IITmII f 0 co N~ 0 0 0 0 0 0 0 co CD "CT N~ 0 co CD )dwel jo ss;9j4uLpqoipui ANS 0 0 qt N~ 0 (>I) einjeje

179..... 0. (0............................ -..............................................0...................... ___ IC)......................................... 0 +..... - --— 0-............................0..........0 C\.................................. ci 0 0 C)f -0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.(0...J0...~C~ 0( (. C'O~( C~) N N N N......... (>i) ssquLjt6!~jqoip-ej

180 - -.......... ----------------------------.......... ---------------------------------- --- ft.................................................................................................................................................................................................................................................. ---- --------- --- ---------------- ------- ---- Mm............................................................................................................. ---------------------------................................................... -------------........................................................................................................................................................................................ 0 -CD 0t 0 -'t 0+ - C ol CDj - a) 0 00. 4l) -Cl), QCl)Z 0+ C/) (I) 1a, -C0 0) -0 C/) (I) 01) E (4i) uJ 00 00 0 0 D (0 't N C1 N N N N 0 1 1 1 1 1 1 1 A I I I I....... 0 00 0 000 C 0 co (0 It NO CD0 N T- '.- ' — i(>) sseu1L45ijqoip-eb 0 0 0 000 C 0 I'll N N~ IIC CD I

181............ ----------------------------------------................................................................................................................................................................................................................................................................. -------------------.......................................................................................................................................................................................................................................................... --------------------..................................................................................................... 77. - ---- -- - - -....................................................................................................................................................................................................................................................................................................................................................................................................... 0 -C'j 0+ - CO) 0 co coC~ C)j 0d cc 0E 04 -0n 0d 0c 0 -So U) a1) 4-*-.2) a) LE5a) N LO -a) E -C.) ul..................................................................... 0 0 O~ ' 0 O (0 0 C' 0 0(0 0 C J CO C\J C\ C\1 C~ C\i T- -- '- ' — ~-r 0 00 C lt CD

182 Wind speed at 1 0 m (mis) 0 0 0 0 0 LO It co) N T Heat flux (W/m2) 0 0 0 0 co (0 '1i C\IJ0 0 0 0) T 0 C\j 0 d (0 It 0 0 o U COJ 0) 0 0 It C 0) CO 0 CI) 0~ z 0i 00000 0 0 0 0 00000 0 0 0 0 (0 C~ '~t ( W /Ml5uo!ju'peu 0 L0 0 10 0 a~ co~ co cN cr (>I) ajnleaidwe~1 LO (.0 C'J 0 (0

183 Rain (mm/cycle) LO) 't co) N~ T- 0............................................................................................. AW........................... CZ ----------------............... ---------------— _ ---....................................................................................................... --------------........................... - - - - - - - - - - - - - - - - - - - - - - - - - - - - -............... ------- ------- -----------....... -. t - - - - - ------ -- 0d 0. Q 0 + C) 0D 0) C ) C E ~ 0> 0 C') 0 02D 0 CO) (0 c'J................... I1II1 I 1 II1 II I'I' I I I I 0 0 0 0 0 0 co 0 I'l C\J (%) Al~piwnH 9A!4~leUj IIIr0 0 0 0 0 0 03) (0 CY) C~J N~ C\j III0 (0 C'%j 'I Il II' 0 LO I C\ (>I) ainjuiedwe 1

184........ 2........... J - --- --- -- 44~ --- — a1) C/) 0 10 ~0 -E. C') U, 0 ' N - 00 N CT) 0 0 -'t - 't 0+ - Cl) 0 U ) a ) -CO -. CO O\ 0) -0 Co 0 -0 N\ Cd 9)c N 0)............................................................................................. --------- - ---................................ I- 0 C\j Co) CD 't C'.] 0 co co I'l C\J C (ZHO/N) IU9.IPBJf5 leilogdS iolo3el u!b UdcjI

185 0 ON (o 0 0+) CC 4 0 co0 0 co * 0 0 0 0 0 0 (wo) Ljdeq] 0

186 B.5 Monthly plots of continuously measured parameters The graphs in this section are monthly compilations of most of the physical parameters measured during REBEX-1. Each month from October, 1992 through April, 1993 is covered by two sets of graphs in consecutive figures. The parameters graphed in the first set for each month are terrain apparent radiobrightnesses, reflector-measured sky radiobrightnesses (no corrections applied), global and net radiation, and wind speed at 10 m. The parameters in the second set for each month are subsurface soil temperatures, vertical heat flux at 2 cm depth in the soil, relative humidity, rainfall, air temperature, and thermal infrared surface temperature. Some terrain and sky radiobrightness data have been edited out as discussed in Section B.4.

187 Radiation (W/m2 0 0 0 000 00 0 00 0 00 00 0 M DItN 0N 't0......S..... (..................... A............................... 2................ I................N N. I x...........................0................. 0...... I.................... t......................... N................................................................................... C. -- -..... 0 0 0 0 00....................... 0.. N...... N.....N. (>.. )....... u......... q........ 0 C0 CY) 0 -Cr, 0 CO 04 0) + 0) LC) C) N 0 co CZ 0 0c d N0 -4-c 0M u co 0~ Nt No 0-;NM 0~ 0) NC C'j 0) 0) 0 I - IC-) 0 IIIaII IIIIII 0 CO I I I I I I I I I I I 0 0 0 0 0 0 0 0 Irl N CO N~ N~ N~ N4 (>i) sseujt46iqoipaj 7Thf0 0 0 0o NM v (s/w) peeds puiM

188 Heat flux (W/m2) 0 0 0 0 0 co co o Oj 0 CJ CO CO ~' C 0 J Rain (mm/cycle)........................................................................................... ----- ---................................................ 77'.................... I......................................................................................................................................... ---------------------------------------.......................... 7 7.................................................. -----------------................................................................................ CZ 0, C-Z cr-.............................................. iii.111111~~~.ii "" '....".~ 1...... 'i/ 1..:'_.................................................. - ---—......... '.........:::.... "........:'..,.............. -..-..........::.::::.:............................:.................................?...i..............................~: ~......... Ca Ca) C) ---- ----..... i... K..................<.. (D - 0 -o 0 CO 'CMj -Co CO -o ' 0 -0) co - cJ C) CM -0) LO - oc-1 - E 0 - CM II -vCc CM cq co c. CM -g - I,,. 0 -0 CO o -4 C)g -CO - I',M JCM -CCO -CMl '0 -CMl CMj rCM 0;-4 4- 0 - 0)... a:... 0 0 II III aI I 0 L 0 LO 0 CO 0 0 0 0 0 0 0 0 0 0 0 0 0) o CO -N N ( CO 0 CO (0 N 0 0) CO N- CD C' CN CM CN C C -r- CO C CM CJ CMJ (>4) ainiaedLuew (%) HU (>I) ainijedwuaeI

189 Radiation (W/m2 0 00 0 00 000 0000 0 000 O M D It 0 N I'CD....... 1,................. -I.-......: I'.'::....... -...................N N...........I.I..........................>. a )............. w......................... I.....................(..........)....................................................... - l..... 3...................... r n - r n -i......h.....T....................... 0.......... 0 ---- 0.....0 0 0 0 0 0.............. 0................. 0.......................-..................s....u................... C)O 0 CD C CY) CY) CMj CY) CY) 0 7 C') CM CO CM~ C) CMH C'f) CMO CM) C CM. CZ~ C')0 C)~ CD~ CI Cf) ~C) CM CO) C) ct CD 0 co) LLI I0 z 2-Il'.......... 0 CMj CY) I ' I I I I I I I 0 0 0 0 0 CD CD 113 co) CMj CM CM (>s) ssault45uqc 0 0 0 0 0 0 (s/uw) peads puiMA

190 HaflxWm2) Rain (mm/cycle) 0 0 0 0 0 co CD It N 0 N LC ~t CO) N '- 0 ------------------------------....................................................................!............................/....................)...(........................................... CO CV) CV) 0 CV) CY) -4 — CC) CM CO C\o N - Cm 0) CZ (o M CMC~ N CY) O1-D CDO ot 0 1J................................................. / C\J (1C).r...fI~CD.Ii~ - i~j~ ClCZ.rM i 0 LO 0 LO 0 LO 0 00 0 0 0 0 0 0 0 0 0 C) CD CD r-N CD CD 0 CD C m 0 C) CD r- CD N N~ N~ N~ N~ N~ N\ ~-co NM N N Nq (>i);ainl~ejdwej (% HU (>4) einjiejdw9_

191 Radiation (W/m2 0 0 0 000 00 0 00 0 00 00 0 ((0 1' NO0 N " (0 i.................................. N N....................2........................................ I...... ' Il~u~lli I~..u........ 0 0 0 0 00..............................* 0 (...N..... N.. N........-....... (.i........ q...p.. 0: 0O (0 (0 CD 0D 0 (0) (0 0 LON LO CM~ CoN cj 0f) -)D c-o ('l0 (0 0d C 04 (0 (0 (0 (0 (0 0o -i CO;-4 0 0 0 0 0 0 N 0 (0 (0 It N (0 (0 N N N N (>t) sseujtLbiqoipeHj 0 0 0 (0 N~ T 0 0 0 0 0 (0 t (0 N r (s/w) peeds pu!MA 0

192 Heat flux (W/m2) Rain (mm/cycle) LO It CO cN 0 0 co o 0 0 CD It N~ L I I -, I 0 0 Nl............. E E Cl co. FE E. E C.) E CD...........[ 1II9 11*............ j......... nil... Tt,1 I........................... t J....... I a................................. ---------.................................................................................................................................................................................................................................................................... --------------------------------------.................................................................................................................... ----------------------- --------------------.................................................................................................................................................................................................................................................................................................................-............................................................................................ CD,CC.. C', ~(0 ~CO -CO -CD ~CY -CO co CY) -LOC - CM LO Lo co co CO CD CO CM -Co CY) CO -co CY).II CM 0) CC S Q 9) 0 0 S 0 Q 0 So 0 U 0 QQ i5 0 0 0 U 9)0 cS 0 EC) x CZ......... - VA:7.............................................. f-77.7 vii.......................................................................................................... IE e,I-J -w r*L a:..1.... 1'...I HU 0 0 0 0C co N- CD LO I NM N~ N~ N~ C\j (>i) 9inj~iedw9_............................. I 03) N~ CDj 0 CD) 0 CD co N- N CD (>4) C9mnjsadw9_

193 Radiation (W/m2 0 00 0 00 00 0 00 0 00 00 0 ' — MD ItN0N 'CD..................... N............. I............ n...... >.................p.-1C........ - -- - - -............2..............................'............. 0 0 0 0 0 0..................... 0. I...................... N...I.. I....-.... (> 4).............................................................................................................................................................................................................................................................................................................. --- -------........................ CO 0) CO 0Y) C\o 0O 0) 0) Coj 0D C) CY) Do. NM M CO+ COH () -D oo C) Eo 0 CO CO0 CD' CO O 0co Cf)Co rCo Co b0 z CZ.0 -, 0..... Z5... E.-. - 0 - CZ -.-:5 '... (z a: 0 0 0 0 0 0 co (0 1q- N 0 co N~ C N N N N - (>N) sseuIL46uqoipeH 0 0 (0 I'l 0 0 L) II 0 0 0 0 Co) N~ ' (S/W) Pee9ds pu!M

194 Rain (mm/cycle) Heat flux (W/m2 0 0 0 0 0 co CD N 0 CM......... C................ 0...........................................................l...... 0 0 0 0 0 0.................. 0............................................................................................................................................................................ I............I... -.1 I......I............I............................I............................... r T 2...........................). CZ 4Arr............................ Mz 0 0 --- —- 0 --- 0..... 0....................N.. N... CO CD Cr) CO CO Cr) C) CO) a) CY) 00 CY) CDl *CD CO) COO COE 0 CO N C 0CO CD> coo CO COl NCO CIO 4Q 04 -04 -4C) 0t (>1) ena~inaedwa i

195 C Radiation (W/m ) 0 00 00 00 0 00 000 0O0 0 Co ( IrN0 N '- o (5) C CD mr I vi MI C\) Col N C) 0 7 N~ C) co C C)j C')O CZ 1-Q C) ~CZ Co2 C~C) 0 co0 Ct Co - D CY~) CY) LL I I I1 1 I1 I I I I I 1I 0 0 0 0 0 0 0 co o (D Nt C 0 co Co0 N~ N~ N N N (>1) sseujtL6uqoipeHJ 0 0 00 0 00 0 0 0 0 0 0 II 0 Co N Co lq Co ~I M~ N N N?-r(>4) ssauwt46pqoipeH (s/W) Peads puiMA

196 Heat flux (W/m2) Rain (mm/cycle) 10 'I CO N TJ 0 o 0 0 0 co CD 11 C\i 0 0 C~J I:................................. --- ---............................................................................................................................................................................... --- -. a1) CZ) CZ).............. a) CDl C\J Col -OD - ') - CZ E -0 C~'J -o 0) CY) 7+a) Q S -a) 'a E a0) uoa) I~ i....... LL................... I 0 L0 0 10 0 10 00 0 00 00 C) co co N- _ CD CD 0 o CD "cl C J C'J C\j C\J C'J C'J CJ C\I T(>) inj~eadwaj (%) HUJ I I I II I I I 0 0 0 0 co N- CD 10 C-I C-I C-I C-I N) ainvieadwe i Cl

197 Radiation (W/m2) 0 0 0 0 0 0 00 0000 0 0 00 -0 CO (D N 0 N 'Ct (D! i!..::... -.-.-..- -:11......................................-... N N N....................................................... 0......?- ------ ~Nz _..'....'%.'..-:.Z ":!'';::'.:....11............. ~ r...........::::'...... w......;::..............::.......................... ':..,.......:. ' ' -.........................................................,,........::..::..........:. '..;:...............................................-.............. - '... "..'............................................................:. L.................................;. -.....................::,..-':: -........ N...... -....................... -........................-.................. - O........................................ -..... -........................... f.... --... r -r........................................ -,,.*:::........................ S....: Lo' 00 CO 4-z LO (0 U) 0a 0 CI0 LO h -. co E CO) 1q1. co,0 CO CMl,o c. - CM - U: 0 t ( bo...... I- -. . -.... -,J. I I 1 I I I I I I I 0 0 0 o 0 0 CO D (0 C l 0 0 Cj CM Cs CM CM (X) sseulq6uqoipea I I I 0 co 0 r 0 0 0 0 0 0 't 0 COD CM O0 t (>) ss9u4q6!jqo!Peu 0 0 0 0 0 (lO C N I(s/w) peads PU!M 0

198 Heat flux (W/m2 0 0 0 0 co CD 111 04 0 Rain (mm/cycle) 0 Nl CD "I (1 N Tr 0..................................... 0.)....................................... o...........0....0 0.......0. 0.................................... H....................................... C.).......................I.................... a.............................. E......................)- - -... H...................................................I.......................... 4CD 04 LO 0 CD -C co 0~ 0 + N ct 0) c CO 0 CZE COo c ct COO co co co 0 Q CD CM 9. 10 0 &) 04 -S) 0 0 LO 0 CD) 0 LO 0 0) co co N- N- CD cD (>i) einijejdwe 0 0 0 0 0) co Nl- CD (>) ainmiedwejL 0 LO N

199 Radiation (W/M2 00 0 00 0c00 00 0 00 00 0 M0 'I' 0N "t0....... N N N....... cCZ )DNl-LO................ O ) C....................... C) N LO C coo.............. CCC 0 CC) 9)) (:O) C) CC 0.o N c - CD~ EC*D CD,z 0C 0 co CJ 0O4 -0C 0 — 0:~ C0 c 0 — C) -j Cc.............. 0 CO, IIIT0 0 CO) I II I I I 0 0 0 0C 0 11 N1 C~1 C~l 0 CI.J C'j 0F 0 0r 0 0 0 0 C0 0 0 C'J C\J 1- -r (>4) sseulq6!iqoipLH 0 0 0 0 0 0 LO) "I CY) C'J r (s/w) paeds puiM (>4) ssau~5jt4qoipeH

200 Heat flux (W/m2) o o 0 0 co o ( - Nt C III I I I I I 0 0 I I I Rain (mm/cycle) LO It CO N 'r- 0 li II I II II II II IIIIII.............% E 0 CD E CD 0 UC o................................................................................................ --- ---- ------- -- a).................................... ----------------------............................................................................................................................................ -... -— m ------------......... I..... CZ)...............2 7 - - - - - -7.-............ CZ Co -Co -co C- 0 ~0 -Co CZ -0 co CM ) O C.) CO)....... z~..u )I i I................................ r-Ir I I I I I I I I I I lII IIIIIIIIIIII I 11 11 I I II I IIII III 0 CO 0 CO 0 CO 0 0) co co r- tr- Co Co N~ N~ N~ N N~ N~ N~ (>I) einjeiedwe 1 o 0 0 00 0 0 Co Co I'l N (%) HU 0 0 0 0 0 0 0) co I!,- Co co N~ N ~ N~ N~

201 B.6 REBEX-1 day to calendar date conversion Use Table B.3 to convert REBEX-1 day numbers to calendar date. The REBEX-1 day number is equivalent to Julian day for 1992 dates and is 366 plus the Julian day for 1993 dates (1992 was a leap year). Calculate fractional day by adding time of day to the day number. For example, 1200 UT on January, 1 1992 becomes fractional day 1.5.

202.. 1992 1993 - -, I -, I - - -., I OCT INOV DEC JAN I FEB I MAR I APR 1 275 306 336 367 398 426 457 2 276 307 337 368 399 427 458 3 277 308 338 369 400 428 459 4 278 309 339 370 401 429 460 5 2791 310 340 371 402 430 461 6 280 311 341 372 403 431 462 7 281 312 342 373 404 432 463 8 282 313 343 374 405 433 464 9 283 314 344 375 406 434 465 10 284 315 345 376 407 435 466 11 285 316 346 377 408 436 467 12 286 317 347 378 409 437 468 13 287 318 348 379 410 438 469 14 288 319 349 380 411 439 470 15 289 320 350 381 412 440 4712 16 290 321 351 382 413 441 472 17 291 322 352 383 414 442 473 18 292 323 353 384 415 443 474 19 293 324 354 385 416 444 475 20 294 325 355 386 417 445 476 21 295 326 356 387 418 446 477 22 296 327 357 388 419 447 478 23 297 328 358 389 420 448 479 24 298 329 359 390 421 449 480 25 299 330 360 391 422 450 481 26 300 331 361 392 423 451 482 27 301 332 362 393 424 452 483 28 302 333 363 394 425 453 484 29 303 334 364 395 454 485 30 304 335 365 396 455 486 31 305 366 397 456 1First day of experiment. 2Last day of experiment. Table B.3: REBEX-1 day to calendar date conversion chart.

203 APPENDIX C RESAMPLING SSM/I RADIOBRIGHTNESSES TO A COMMON GRID This appendix describes the how the SSM/I brightnesses in Chapter 6 were resampled to the geographically fixed coordinates of the Equal-Area SSM/I Earth Grid (EASE-Grid). The resampling process takes two steps: (a) Identification of a set of 16 SSM/I boresight loci that surround a grid point, and (b) interpolation of the SSM/I brightnesses to the closest of a set of predetermined interpolation loci in the sensor reference frame. The set of interpolation loci are dense enough that the error in re-mapping by nearest neighbor selection is reduced below the geolocation error of the data. Interpolating to points that are fixed with respect to the satellite reference frame enables the use of time-saving interpolation coefficients calculated once for a particular sensing geometry. The advantages of this process over conventional mapping techniques are two-fold. First, brightness variation with time can be somewhat removed from spatial variation since the earth grid is fixed for all times and sensors. Secondly, the interpolation process can include the simulation of arbitrary sensor patterns. The antenna pattern of the SSM/I varies with the frequency of its channels, with the 19 GHz channel viewing the largest area of earth (see Table 6.1). To make inter-channel comparisons measurements must represent the brightness of the same geographical area. This can

204 be achieved by choosing a single optimal interpolated pattern and using it for all the channels. Here, the chosen pattern is that of the 19 GHz V-polarized channel. The first section of this appendix details the optimizing interpolation scheme for geophysical data first developed by Stogryn [64] and first adapted for SSM/I by Poe [65]. Section C.1.1 describes implementation of the method for the available SSM/I data stream. Section C.1.2 discusses the use of a common interpolated pattern and section C.2 describes remapping to the EASE-Grid system. C.1 Optimal interpolation An interpolation filter for geophysical data may be described as optimal if interpolation to a particular location gives the same value as would have been measured had the sensor boresight been pointed there. Poe [65] presents a method for closely approximating this optimal interpolation for SSM/I data based on the application of Backus-Gilbert inversion methodology [66] [64]. The SSM/I views the earth with a conical scanning geometry as shown in Figure C.1. In terms of an instantaneous antenna gain function, (o(t'),(t')), the antenna temperature measured in the direction SA is given by: TA(SA) = J d h(t') JJ dQ g (o(t'), A(t'))TB (P, ~(t'), t') (C.1) -00 E where So(t') is the instantaneous boresight direction of the antenna, s(t') is the unit vector from the antenna in the direction of the solid angle dQ, p is the position vector of a spherical earth surface coordinate system, and h is the impulse response of the radiometer receiver low-pass filter. TB(P,S.(t'),t') stands for the brightness temperature at a point p, time t', and in the direction s(t'). The region of integration, E, encompasses the entire surface of the earth intercepted by s(t'). Note that both h and g are normalized such that fJ+~ dt h = 1 and fE dQ = 1. The antenna boresight

205 - -- -- - Figure C.1: SSM/I (F08 platform) scan geometry showing earth and satellite coordinate system reference vectors. p and PA are normal to the earth's surface. The F11 platform's geometry differs only in the direction of the ground track. (Modified from [59].)

206 pointing direction, SA, is given by: r+oo 1 r+/2 SA = dt' h(t')o(tl) = - dt' So(t') (C.2) -oo T J-/2 where h has been replaced by integration over,. the characteristic time of the receiver filter. The differential solid angle, dQ, can be replaced by the corresponding differential area of the earth's surface according to: dQ = [ (t') * 1 dA. (C.3) L s(t) Ja Exchanging the time and space integration, we can rewrite (C.1) as: TA (PA) = TA(SA) 1 fdA -/2 dt {(So(t), ())[ (t) ]TB (P ' S E (C.4) where the vector PA points from the center of the earth to the intersection of the boresight with the earth's surface. For simplicity, assume that the brightness temperature upwelling from an area dA is both independent of direction and, on average, independent of time. This assumption constitutes a limitation in the analysis only if the TB at a particular dA changes significantly over the range of directions or times at which the overlapping measurements view it. Nevertheless, for practical purposes, the functional dependence of TB is reduced to: TB(p (t'), t') = TB(p-) (C.5) Then TB comes out of the time integral in (C.4) and we can define an effective antenna gain function: PA = 1 / dt' g(o(t) S (t')) [ -(t' (C.6) P - J-.-/2 S2

207 Combining (C.4), (C.5), and (C.6) yields this expression for the measured antenna temperatures: TA(PA) = fJ dA(PA )TB(p. (C.7) E Note that (C.7) is the forward equation for TA(PA) when PA is the geolocation of an SSM/I sensor measurement. We seek a way to use the information in (C.7) to optimize the interpolation of TA to an arbitrary p within the swath. Without a loss of generality, we can constrain the set of desired interpolation loci, Pd, to fall along arcs offset from the scans and along lines parallel to the ground track as shown in Figure C.2. An estimate of the antenna temperature TA(Pd) that would have been measured by the sensor had it been pointed at pd is given by a weighted linear combination of a set of neighboring sensor measurements: N TA(Pd) = E aTA(pAi) (C.8) i=l where {ai} is a set of N weighting coefficients to be determined and {TA(P4i)} is the set of corresponding sensor measurements. Substituting (C.7) into (C.8) yields: N TA(d) = ] dA Z ai(pA, p)TB(p. (C 9) E i=1 Comparison of (C.9) to (C.7) leads to the definition of an interpolated effective antenna gain function as: N 9I(pd, p) = - aG(pAi, ). (C.10) i=1 To determine the {ai}. we would like to minimize the error, e, in the estimate of

208 Scan Angle 12.5 km I I 12.5 km Along Scan (j) Along Track (i) Figure C.2: A section of the SSM/I swath showing parts of five scans. The section is about 27~ past the subsatellite track along the arc of the scans. The circles approximate the 3 dB extent of each measurement's effective field of view. The dashed box delineates the set of 16 measurements used to interpolate to any of the points in the center of the box. For example, an idsealized interpolated antenna pattern for the point marked with an X is shown as a dashed ellipse. The separation between rows and samples is 12.5 km for 85 GHz sampling and 25 km for other SSM/I channels.

209 e = TA(Pd) -TA(Pd)\ =IJJdA [g(pdp) - gI(PdP )TB( (C.ll) E It is clear from (C.11 ) that e is minimized if 9 is closely approximated by the interpolated pattern, I9-which is. at the same time, a desirable condition in that antenna pattern information remains consistent in the resampling process. From [64], the solution for the weighting coefficients, {a;}, is obtained by minimization of Qd Qd = dA [(Pd, P' - ai, (PAiP, J(1Id, P) (C.12) where the function J(pd,p) can be used, for example, to reduce the sidelobe levels in GI but is set to unity in this analysis. Using Backus-Gilbert inversion theory [66], Stogryn [64] gives the {a;} that minimize Qd: u g gu - - [- - - (C 13) where the elements of the matrix g and the vectors ui and v are: g = J dAG(pAp}.(pAj( (C.14) E u= JJdA (PAi p) (C.15) E v = fdA g(4pA), 5)( pd, P) (C.16) E Since this analysis does not take system noise into account, the variance in the resulting interpolated antenna temperatures should be examined to insure that it does not exceed acceptable levels. We assume that the sensor antenna temperature measurement noise is uncorrelated between samples and that the variance of the

210 DMSP Platform F08 F11 Launch date 19 June 1987 28 Nov. 1991 Viewing direction Aft looking Forward looking Ascending equator crossing 06:15 17:04 time (local time) Maximum altitude 882 km 878 km Table C.1: SSM/I operational data for the F08 and F11 satellite platforms [60]. noise is the same for all samples. The resulting variance in the interpolated antenna temperature is given by: (ATA)2 = Ea = (ATA)2 a (C.17) i=1 where E is the covariance matrix of measurement noise for the samples contributing to TA(Pd). In practice, ATA is often less than ATA because many noisy samples are combined in the interpolation process. C.1.1 Implementation Calculation of interpolation coefficients-the {ai} in (C.8)-requires the instrument antenna patterns which were acquired from Hollinger [67]. These antenna patterns were measured from the SSM/I instrument aboard the F08 DMSP platform. Because the 85 GHz channel of the F08 SSM/I failed in 1990, this thesis presents data from the Fll platform. The sensor differences are negligible for the purposes of this work. Table C.1 compares parameters for the F08 and Fll platforms. (Parameters for the individual channels can be found in Table 6.1.) I acquired global NESDIS Level lb Format SSM/I data for the period of REBEX-1 from the Marshall Space Flight Center (MSFC) Distributed Active Archive Center in Huntsville, Alabama. MSFC provided computer programs for extracting data for a particular geographical region (ssmillblatlon.c) and for converting raw sensor counts

211 to antenna temperatures (ssmillbta.f). The interpolation coordinates in the satellite reference frame were based on a four-times increase in scan and sample frequency. That is, the number of samples per scan was increased from 64 to 256 for the low frequency channels (19, 22, and 37 GHz) and from 128 to 512 for the 85 GHz channel. C.1.2 Resampling to a common resolution The above derivation assumes in (C.12) that the optimal antenna pattern is that of the channel that sampled TA(PA ). But the need for multi-channel brightness comparisons suggests that a single pattern is desirable for all polarizations and frequencies. If the chosen pattern is 9d(Pd, P, then the minimized function Qd from (C.12) becomes: Qd= dA[dpdp- ag p JpdP (C.18) and the parameters of (C.13) are exactly the same except for vi, which is now given by: vI= Jf/dA (pi, P-)d(Pd ). (C.19) E In this thesis, the common pattern is that of the 19 GHz V-polarized channel of the SSM/I. The 19 GHz pattern has the largest footprint of the four SSM/I frequencies, and the interpolated patterns of the other channels are able to fit it without high side-lobes. If the 19 GHz pattern were interpolated to a significantly smaller desired pattern-that is, either the 37 GHz or 85 GHz pattern-then the best achievable interpolated pattern would be distorted and have high side-lobe levels. Figure C.3 shows example interpolated antenna patterns for the 37 and 85 GHz channels compared with the optimal 19 GHz pattern. Because samples are more closely distributed near the scan edges, the degree to which the interpolated pattern

212 (a) (b) (c) Figure C.3: Antenna patterns: (a) 19 GHz V-pol. effective antenna pattern, (b) 37 GHz V-pol. interpolated to 19 GHz pattern at center of scan, (c) 37 GHz V-pol. interpolated to 19 GHz pattern between edge and center of scan, (d) 85 GHz V-pol. interpolated to 19 GHz pattern at center of scan, (e) 85 GHz V-pol. interpolated to 19 GHz pattern between edge and center of scan. The contours are at 3, 6, 12, and 24 dB.

213 (d) (e) Figure C.3: (Continued from previous page.) matches the optimal pattern depends upon where in the scan the interpolation occurs. Figure C.3 shows this effect in the 37 GHz patterns. At the scan edge, the pattern closely matches the elliptical shape of the optimal pattern while at the the center of the scan the interpolated pattern is more rectangular. The 85 GHz patterns are distinctly rhomboidal at all scan angles because the antenna pattern at 85 GHz is more than three times smaller than the 19 GHz pattern and the density of samples is four times greater (see Table 6.1). To better approximate a 19 GHz optimal pattern, interpolation could be performed with more than 16 85 GHz samples. This is a subject for future investigation. C.2 Remapping to the Equal-Area SSM/I Earth Grid As discussed above, the interpolation loci are defined in the satellite reference frame relative to the sample points of the SSM/I sensors. This yields a nominal sample spacing of 6.25 km (or 3.125 km for 85 GHz high-resolution resampling). Consequently, nearest neighbor re-registration to fixed coordinates on the earth has

214 a maximum error of about 3.125 km (1.563 km for 85 GHz), which is 7% of the 19 GHz field of view. Fixed earth coordinates are defined by the Equal-Area SSM/I Earth Grid (EASEGrid), used by the National Snow and Ice Data Center (NSIDC) for global SSM/I data archiving [68]. This thesis uses two EASE-Grid projections: Ml (low resolution mercator) and Mh (high resolution mercator, for 85 GHz only). Column number (r) and line number (s) of the projections are based on a cylindrical equal-area map true at 30~ N/S latitude: r= RAcos30+ rO (C.20) C R sin< R= - + so (C.21) C cos 30 where R is the radius of the earth (6371.288 km), C is the nominal cell size, A is the latitude in radians, and < is the longitude. The projection origin (rO,sO) is the point in the projection where A = 0 and 0 = O. Table C.2 lists the parameter values for the Ml and Mh grids. Grid C rO sO Ml 25.067525 691.0 292.5 Mh 12.5337625 1382.0 585.0 Table C.2: Parameters for the Ml and Mh EASE-Grid projections. Figure C.4 shows example images made from EASE-Grid data in both the Ml and Mh formats. The data are displayed in a conical projection but the coordinates of each pixel are at EASE-Grid locations. Note that the four-times increase in sample density of the Mh grid and the smaller EFOV highlight variation not visible in the low resolution image. The displayed data was subset to the region bounded by 41 N latitude, 51 N latitude, 89 W longitude, and 110 W longitude. Note that the edge of the swath is to the east of 110 W.

21.5 -113 85 GHz Horizontal Polarization, Low Res. 1992:399:12:21 UTC -1.05 -97 -89 lll'L 160 180 200 220 240 260 280 300 Rodiobrightness (K) -1,13 85 GHz Horizontal Polarization, High Res. 1992:399 12:21 UTC -1.05 -97 17sl_' c :i i:i 160 180 200 220 240 260 280 300 Radiobrightness (K) Figure C.4: SSM/I 8.5 GHz h-pol. images from February 2, 1992 at 12:21 T'1'( demonstrating resampling to low resolution (top) and high resolution (bottom) EFOV.

BIBLIOGRAPHY [1] Giorgi, F., and Mearns, L. O., "Approaches to the Simulation of Regional Climate Change: A Review," Rev. Geophys., Vol. 29, No. 2, pp. 191-216. 1991. [2] Dickinson, R. E.. Henderson-Sellers, A., Kennedy, P. J., and Wilson. 'M. F.. Biosphere-Atmosphere Transfer Scheme (BATS) for the NCAR Cornulnitnty Climate Model, NCAR Technical Note, NCAR/TN-275+STR, December. 1986. [3] Sellers, P.J., Mintz, Y., Sub Y.C.. and Dalcher, A., "A Simple Biosphere Model (SiB) for use within general circulation models,' J. Atm. Sci., Vol. 43, No. 6. pp. 505-531. 1986. [4] Peixoto. J. P.. and Oort, A. H., Physics of Climate, American Institute of Physics, New York. 1992. [5] Dickinson. R. E., Errico. R. M.. Giorgi, F., and Bates, G. T., "A regional climate model for the Western United States." Climate Change. Vol. 15, pp. 383-422, 1989. [6] C(less. R. D., Potter. G. L., Zhang, M.-H., Blanchet, J.-P., Chalita, S., Colman, R., Dazlich. D. A., Del Genio, A. D., Dvmnikov, V., Galin, V., Jerrett, D., Ieup, E., Lacis, A. A., Le Treut, H., and others "Interpretation of Snow-Climate Feedback as Produced by 17 General Circulation Models," Science. Vol. 253, pp. 888-892. 1991. [7] Edgerton, A. T.. Stogryn, A.. and Poe, G., Microwave radiometric investigations of snowpacks. Final Report, 1285R-4, Aerojet-General Corporation, July. 1971. [8] Schmugge, T., Wilheit, T. T.. Gloerson, P., Meier, M. F., Frank. D., and Dirmhirn, I.. "Microwave signatures of snow and fresh water ice," Presented at the Interdisciplinary Symposium on Advanced Concepts and Techniques in the Study of Snow and Ice Resources, Goddard Space Flight Center. Greenbelt, MD, Nov. 1, 1973. [9] Stiles. W. H., and Ulaby, F. T., "The active and passive microwave response to snow parameter 1. Wetness,"J. Geophys. Res., Vol. 85, No. C2, pp. 1037-1044, 1980. 216

217 [10] Stiles, W. H., and Ulaby, F. T., "The active and passive microwave response to snow parameter 2. Water equivalent of dry snow," J. Geophys. Res., Vol. 85, No. C2, pp. 1045-1049, 1980. [11] Matzler, C., Schanda, E., and Good, W., "Towards the definition of optimum sensor specifications for microwave remote sensing of snow," IEEE Trans. Geosci. Rem. Sens., Vol. GE-20, No. 1, pp. 57-66, 1982. [12] Sturm, M., Grenfell, T. C., and Perovich, D. K., "Passive microwave measurements of tundra and taiga snow covers in Alaska, U.S.A.," Annals Glaciology, Vol. 17, pp. 125-130, 1993. [13] Zwally, H. J., "Microwave emissivity and accumulation rate of polar firn," J. Glaciol., Vol. 18, No. 79, pp. 195-215, 1977. [14] Foster, J. L., Rango, A., Hall, D. K., Chang, A. T. C., Allison, L. J., and Diesen, B. C., "Snowpack monitoring in North America and Eurasia using passive microwave satellite data," Remote Sens. Environ., Vol. 10, pp. 285-298, 1980. [15] Hall, D. K., Sturm, M., Benson, C. S., Chang, A. T. C., Foster, J. L., Garbeil, H., and Chacho E., "Passive microwave remote and in situ measurements of arctic and subarctic snow covers in alaska," Remote Sens. Environ., Vol. 38, pp. 161-172, 1991. [16] Hofer, R., and Good, W., "Snow parameter determination by multichannel microwave radiometry," Remote Sens. Environ., Vol. 8, pp. 211-224, 1979. [17] Hofer, R., and Shanda, E., "Signatures of snow in the 5 to 94 GHz range," Radio Sci., Vol. 13, No. 2, pp. 365-669, 1978. [18] Goodison, B. E.. "Determination of areal snow water equivalent on the Canadian praires using passive microwave satellite data," in Proceedings of IGARSS89, Vancouver, Canada, July 10-14, 1989. [19] Hall, D. K., Foster, J. L., and Chang, A. T. C., "Measurement and modeling of microwave emission from forested snowfields in Michigan," Nord. Hydrol., Vol. 16, pp. 129-138, 1982. [20] Chang, A. T. C., Foster, J. L., and Hall, D. K., "Nimbus-7 SMMR derived global snow cover parameters," Annals Glaciology, Vol. 9, pp. 39-44, 1987. [21] Chang, A. T. C., Foster, J. L., and Rango, A., "Utilization of surface cover composition to improve the microwave determination of snow water equivalent in a mountain basin," Int. J. Remote Sens., Vol. 12, No. 11, pp. 2311-2319, 1991. [22] Rango, A., Chang A. T. C., and Foster, J. L., "The utilization of spaceborne microwave radiometers for monitoring snowpack properties,", Nord. Hydrol., Vol. 10, pp. 25-40, 1979.

218 [23] Wang, J. R., Chang, A. T. C., and Sharma, A. K., "On the estimation of snow depth from microwave radionretric measurements," IEEE Trans. Geosci. Rem. Sens., Vol. 30, No. 4, pp. 785-792, 1992. [24] Kunzi, K. F., Patil, S., and Rott, H., "Snow-cover parameters retrieved from Nimbus-7 Scanning Multichannel Microwave Radiometer (SMMR) Data," IEEE Trans. Geosci. and Remote Sensing, Vol. GE-20, No. 4, pp. 452-467, 1982. [25] Chang, T. C.. Gloersen, P., and Schmugge, T., "Microwave emission from snow and glacier ice,"J. Glaciol., Vol. 16, No. 1, pp. 23-39, 1976. [26] McFarland, M. J., Wilke, G. D., and Harder, P. H. III, "Nimbus 7 SMMR investigation of snowpack properties in the northern great plains for the winter of 1978-1979," IEEE Trans. Geosci. and Remote Sensing, Vol. GE-25, No. 1, pp. 35-46, 1987. [27] Grody, N. C., "Surface identification using satellite microwave radiometers,"IEEE Trans. Geosci. and Remote Sensing, Vol. 26, No. 6, pp. 850-859, 1988. [28] Grody, N. C., "Classification of snow cover and precipitation using the Special Sensor Microwave Imager," J. Geophys. Res., Vol. 96, No. D4, pp. 7423-7435, 1991. [29] Fiore, J. V. jr., and Grody, N. C., "Classification of snow cover and precipitation using SSM/I measurements: case studies,"Int. J. Remote Sens., Vol. 13, No. 17, pp. 3349-3361, 1992. [30] Neale, C. M. U., McFarland, M. J., and Chang, K., "Land-surface-type classification using microwave brightness temperatures from the special sensor microwave/imager," IEEE Trans. Geosci. Rem. Sens., Vol. 28, No. 5, pp. 829-838, 1990. [31] Stiles, W. H., Ulaby, F. T., and Rango, A., "Microwave measurements of snowpack properties,"Nord. Hydrol., Vol. 12, pp. 143-166, 1981. [32] Chang, A. T. C., Foster, J. L., and Rango, A., "The role of passive microwaves in characterizing snow cover in the colorado river basin," GeoJournal, Vol. 36, No. 3, pp. 381-388, 1992. [33] Davis, R. E., Dozier, J., and Perla, R., "Measurement of snow grain properties," in Seasonal Snowcovers: Physics, Chemistry, Hydrology, D. Reidel, p. 63-74, 1987. [34] Shi, J., Davis, J.,. E.Davis, R. E., and Dozier, J., "Stereological determination of dry-snow parameters for dicrete-scatterer microwave modeling," Annals Glaciology, Vol. 17, pp. 295-299, 1993.

219 [35] Brun, E., Martin, E., Simon, V., Gendre, C., and Coleou, C., "An energy and mass model of snow cover suitable for operational avalanche forecasting," J. Glaciol., Vol. 35, No. 121, pp. 333-342, 1989. [36] Loth, B., Graf, H. F., and Oberhuber, J. M., "Snow cover model for global climate simulations," J. Geophys. Res., Vol. 98, No. D6, pp. 10451-10464, 1993. [37] Anderson, E. A., A Point Energy and Mass Balance Model of a Snow Cover, A dissertation submitted to the Dept. of Civil Engineering and the committee on graduate studies of Stanford University, December 1, 1975. [38] Jordan, R., A one-dimensional temperature modelfor a snow cover: Technical documentation for SNTHERM.89, U.S. Army Corps of Engineers, Cold Regions Research and Engineering Laboratory, Special Report 91-16, 1991. [39] England, A. W., "Radiobrightness of Diurnally Heated, Freezing Soil," IEEE Trans. Geosci. Rem. Sens., Vol. 28, No. 4, pp. 464-476, 1990. [40] Nakano, Y., and Brown, S., "Mathematical modeling and validation of the thermal regimes in tundra soils, Barrow, Alaska," Arctic and Alpine Res., Vol. 4, No. 1, pp. 19-38, 1972. [41] Anderson, D. M., Tice, A. R., and McKim, H. L., "The unfrozen water and the apparent specific heat capacity of frozen soil," in Permafrost: North American contributions [to the] Second International Conference, Yakutsk, Siberia, pp. 289 -295, 1973. [42] De Vries, D. A., "Thermal properties of soils", in Physics of Plant Environment, W. R. Van Wijk, ed., North-Holland, Amsterdam, pp. 210-235, 1963. [43] de Vries, D. A., and Philip, J. R., "Soil heat flux, thermal conductivity, and the null-alignment method," Soil Sci. Soc. Am. J., Vol. 50, pp. 12-18, 1986. [44] de Vries, D. A., and Afgan, N. H., Heat and Mass Transfer in the Biosphere 1. Transfer Processes in Plant Environment, John Wiley and Sons, New York, 1975. [45] Brandt, R. E., and Warren, S. G., "Solar-heating rates and temperature profiles in Antarctic snow and ice," J. Glaciol., Vol. 39, No. 131, pp. 99-110, 1993. [46] Panofsky, H. A., and Dutton, J. A., Atmospheric Turbulence, John Wiley, New York, 1984. [47] Smith, E. A., Crosson, W. L., and Tanner, B. D., "Estimation of surface heat and moisture fluxes over a prairie grassland 1. In situ energy budget measurements incorporating a cooled mirror dew point hygrometer," J. Geophys. Res., Vol. 97, No. D17, pp. 18557-18582, 1992.

220 [48] Press, W. H., Flannery, B. P., Teukolsky, S. A., and Vetterling, W. T., Numerical Recipes: The Art of Scientific Computing, Cambridge University Press, Cambridge, 1986. [49] Chandrasekhar, S., Radiative Transfer, Dover Publication, Inc., New York, 1960. [50] Ulaby, F. T., Moore, R. K., and Fung, A. K., Microwave Remote Sensing: Active and Passive, Vol. I, Addison Wesley, 1981. [51] Bohren, C. F., and Huffman, D. R., Absorption and Scattering of Light by Small Particles, John Wiley & Sons, New York, 1983. [52] Ishimaru, A., and Kuga, Y., "Attenuation constant of a coherent field in a dense distribution of particles," J. Opt. Soc. Am., Vol. 72, No. 10, pp. 1317-1320, 1982. [53] Tsang, L., Mandt, C. E., and Ding, K. H., "Monte carlo simulations of the extinction rate of dense media with randomly distributed dielectric spheres based on solution of Maxwells equations," Optics Letters, Vol. 17, No. 5, pp. 314-316, 1992. [54] Kuga, Y., Ulaby, F. T., Haddock, T. F., and DeRoo, R. D., "Millimeter-wave radar scattering from snow: 1. Radiative transfer model," Radio Sci., Vol. 26, No. 2, pp. 329-341, 1991. [55] Tsang, L., Kong, J. A., and Shin, R. T., Theory of Microwave Remote Sensing, John Wiley & Sons, New York, 1985. [56] Hufford, G., "A model for the complex permittivity of ice at frequencies below 1 THz," Int. J. IR and MM Waves, Vol. 12, No. 7, pp. 677-682, 1991. [57] John Kendra, personal communication, 1995. [58] Galantowicz, J.F. and England, A.W., Field data report for the First Radiobrightness Energy Balance Experiment (REBEX-1), October 1992-April 1993, Sioux Falls, South Dakota, Univ. of Michigan, Radiation Lab Technical Report RL-913, February, 1995. [59] Hollinger, J., DMSP Special Sensor Microwave/Imager Calibration/Validation Final Report Volume 1, Naval Research Laboratory, Washington, DC 20375 -5000, 1989. [60] User Guide to Special Sensor Microwave/Imager (SSM/I) Data (NESDIS Level lb Format), Marshall Space Flight Center, Distributed Active Archive Center, Huntsville, Alabama, 1994. [61] Colbeck, S. C., "Vapor-pressure dependence on temperature in models of snow metamorphism," J. Glaciol., Vol. 36, No. 124, pp. 351-353, 1990. [62] Arya, S. P., Introduction to Micrometeorology, Academic Press, San Diego, 1988.

221 [63] Farouki, 0. T., Evaluation of methods for calculating soil thermal conductivity, CRREL Report 82-8, U. S. Army Cold Regions Research and Engineering Laboratory, Hanover, NH, 1982. [64] Stogryn, A., "Estimates of brightness temperatures from scanning radiometer data," IEEE Trans. Antennas and Propagation, Vol. AP-26, No. 5, pp. 720-726, 1978. [65] Poe, G. A., "Optimum interpolation of imaging microwave radiometer data," IEEE Trans. Geosci. Rem. Sens., Vol. 28, No. 5, pp. 800-810, 1990. [66] Backus, G., and Gilbert, F., "Uniqueness in the inversion of inaccurate gross earth data," Phil. Trans. Roy. Soc. London, Vol. A226, pp. 123-192, 1970. [67] J. Hollinger, personal communication. Space Sensing Branch, Naval Research Laboratory, Washington, DC 20375-5000, 1992. [68] README file for NOAA/NASA Pathfinder Program SSM/I Level 3 Brightness Temperatures Equal-Area SSM/I Earth Grid (EASE-Grid) CD, National Snow & Ice Data Center, Boulder, CO, 1995.