Melting trends over the Greenland ice sheet ( 1958 – 2009 ) from spaceborne microwave data and regional climate models

Introduction Conclusions References


Introduction
Melting over the Greenland ice sheet (GrIS) has been accelerating over the 2000's, as suggested by recent satellite-based observations (Hall et al., 2008;Tedesco et al., 2008;Wouters Correspondence to: X. Fettweis (xavier.fettweis@ulg.ac.be) et al., 2008) and model simulations (Van den Broeke et al., 2009;Tedesco et al., 2011).This can be partially attributed to recent observed warming of the Arctic, generally attributed to increases in atmospheric greenhouse gas loading (Hanna et al., 2008).However, supporting this recent melt acceleration over a long enough period (at least 30 yr) is needed to confirm that the changes currently occurring over the GrIS are significant.
Satellite-based mass estimates began in the early 1990's with the European Remote-sensing (ERS) satellites (Zwally et al., 2005), followed by the NASA Gravity Recovery and Climate Experiment (GRACE) (Luthcke et al., 2006;Wouters et al., 2008) but a 15 yr time series cannot capture long-term mass variations of the GrIS, given the large inter-annual variability observed in the mass balance and surface melt (Tedesco et al., 2011).Measurements from coastal weather stations managed by the Danish Meteorological Institute (DMI) have been available for more than a century.However, in this case, data have been collected at local spatial scales and therefore, cannot be assumed to be representative of other places (e.g., interior of the GrIS).Nevertheless, these data can be used to quantify the inter-annual variability (Fettweis et al., 2008;Wake et al., 2009) at selected locations.Results from regional climate models (RCMs) can provide estimates of GrIS surface changes over a long period at high resolution, but measurements to validate these models are limited.
With more than 30 yr of data, the spaceborne passive microwave brightness temperature data set offers a unique opportunity to study changes of near surface melting over the GrIS (Abdalati and Steffen, 2001;Mote, 2007;Tedesco, 2007).When snow melts, the presence of liquid water within the snowpack increases the snow microwave emissivity (Ulaby and Stiles, 1980).Consequently, microwave Published by Copernicus Publications on behalf of the European Geosciences Union.
X. Fettweis et al.: Melting trends over the Greenland ice sheet  brightness temperatures (e.g., the product of snow temperature and emissivity) are considerably higher in the case of wet snow than those in the case of dry snow.Several algorithms have been developed to detect melting from spaceborne passive microwave data (Abdalati and Steffen, 1997;Picard and Fily, 2006;Mote, 2007;Tedesco, 2007).
In this study, we analyze results from two Regional Climate Models (RCMs) in conjunction with those obtained from spaceborne microwave brightness temperatures to study surface and near-surface melt variability over the GrIS since 1979 and extend the analysis back to 1958.The RCMs used in this study are MAR (for Modèle Atmosphérique Régional), described by Gallée and Schayes (1994) and Fettweis (2007) and RACMO2, described by Ettema et al. (2009).Both RCMs are fully coupled with an energy balance-based snow model allowing feedbacks between the surface and the atmosphere.The models have been forced at the boundaries by the re-analyses of the European Centre for Medium-Range Weather Forecasts (ECMWF) since September 1957.Both models have demonstrated their ability to reliably simulate the GrIS SMB since 1958 (Ettema et al., 2009;Fettweis et al., 2011).Moreover, the MAR model was already satisfactorily compared with passive microwavederived observations (Fettweis et al., 2005(Fettweis et al., , 2007) ) and was used to improve the microwave-based algorithm developed by Abdalati and Steffen (2001) for retrieving melt extent (Fettweis et al., 2006).Combining results from modeling and measurement tools (e.g., satellite data) would provide an ideal framework to overcome some of the intrinsic limitations/uncertainties of each tool considered separately.A major outcome of this study is a first assessment of the potential assimilation of melt detection from spaceborne microwave data into the RCMs, with the ultimate scope of reducing uncertainties in model-based Surface Mass Balance (SMB) estimates.
This study complements a previously published study of Fettweis et al. (2007) by using two different RCMs outputs rather than only MAR, hence allowing a more robust analysis, especially in those cases when the results from models and remote sensing algorithms data are not consistent.The data used are described in Sect. 2. In Sect.3, two spaceborne passive microwave melt retrieval techniques are discussed and their results compared with in-situ measurements.In Sect.4, RCMs outputs are used to assess the applicability of these algorithms over areas and periods where no measurements are available.Finally, in Sect.5, the melting variability is estimated before beginning satellite observations using the two models.The aims of this paper are triple: (i) to select melt retrieval algorithms that can be easily compared with RCM outputs, (ii) to inter-validate both RCMs and the satellite-retrieved melt extent and, (iii) to show the unprecedented behavior of the melt extent of the last 50 yr.

Data
Satellite data used here consist of spaceborne passive microwave brightness temperatures recorded by the Scanning Multichannel Microwave Radiometer (SMMR) satellite (1979)(1980)(1981)(1982)(1983)(1984)(1985)(1986)(1987), the Special Sensor Microwave/Imager (SSM/I) flying on the F-08 (1987-1991), F-11 (1992-1994) and F-13 (1995-2009) satellites and by the Special Sensor Microwave Imager/Sounder (SSMIS) on the F-17 satellites (2009).Microwave data are available at four frequencies: 18.7 GHz (respectively 19.35 GHz for SSM/I(S) data) at both vertical (V) and horizontal (H) polarizations, 22.2 (V), 37.0 (V, H) and, 85.5 GHz (V, H).Data are reprojected on a 25 km × 25 km grid (EASE-Grid) and distributed by the National Snow and Ice Data Center (NSIDC, Boulder, Colorado) (Armstrong et al., 1994;Knowles et al., 2002).We use the near-real time SSMI/S brightness temperatures after May 2009 (Brodzik and Armstrong, 2008).Data acquired during both ascending and descending passes are averaged to generate daily values of brightness temperature, as in Abdalati andSteffen (1997, 2001).In addition, data gaps for periods shorter than five days are filled through linear interpolation.Most of the gaps occur in the SMMR data set (knowing that there is data only each 2 days) or in the end of live of the SSM/I satellites (see Annex A in the Supplement).For the gap larger than 5 days, we consider that both satellite and RCM do not detect melt.This represents about 1% of the whole sample.
MAR outputs are available at 25 km horizontal resolution  using the same model set-up described by Fettweis (2007) but with some improvements in the surface albedo scheme according to Fettweis et al. (2011).Preliminary comparisons between the modeled and satellite-derived melt extent using the results of the MAR version from Fettweis et al. (2011) highlighted some biases in the MAR surface albedo scheme which are now corrected for the results presented here.For RACMO2, we used the 11 km-results from Ettema et al. (2009Ettema et al. ( , 2010)).The atmospheric modules of MAR and RACMO2 (atmospheric dynamics, microphysics, turbulence scheme, etc.) are of similar complexity.Both RCMS are fully coupled with a multi-layered energy balance snow model taking into account the meltwater percolation, retention and refreezing, the bare ice presence, the increase in snow density due to refreezing of liquid water and the packing of dry snow.Only the surface albedo scheme is rather different.In RACMO2, albedo is a simple function of the snow density and cloudiness (Greuell and Konzelmann, 1994) while in MAR it depends on the shape and size of the snow grains (following the CROCUS snow albedo parametrization developed by Brun et al., 1992), on the cloudiness and on the zenithal angle (Lefebre et al., 2003).
Results from RACMO2 and satellite data are re-projected on the 25 km MAR grid using an inverse distance weighted interpolation.Ice sheet masks for the three different data sets are saved and the analysis is performed only over those pixels where the ice sheet is assumed to be present on all of the three different ice masks.For the spaceborne passive microwave data set, we use the ice mask from the EASE-Grid LOCI mask (derived from Boston University Version of Global 1 km Land Cover from MODIS (Moderate Resolution Imaging Spectroradiometer) 2001, Version 4).Finally, ice pixels along the coast are excluded from our analysis in order to reduce the uncertainty related to the sea ice/ocean mixed pixel effect.

Satellite-based algorithms for melt detection
Several approaches have been proposed for mapping wet snow over the Greenland and Antarctica ice sheets from spaceborne passive microwave measurements.Most of them use data collected at K-band (18.7 GHz for SMMR and 19.35 GHz for SSM/I(S) data) horizontal polarization (T19H) because of the largest observed difference of this combination of frequency and polarization between dry and wet snow conditions (e.g.Liu et al., 2006;Tedesco, 2007).
The two main methods for wet snow detection are -the use of a threshold value of T19H as follows Equation ( 1) is also often used with a combination of brightness temperatures using multiple frequencies and polarizations in replacement of only T19H to take advantage of their differing responses to the liquid water content (LWC) increase inside the snow pack.
The most well-known is the approach of Abdalati and Steffen (2001), which is based on a combination of the K-band horizontal polarized brightness temperature (T19H) and the Ka-band (36.5 GHz) vertical polarized brightness temperature (T37V) (called cross-polarized gradient ratio, XPGR) to detect liquid water: -the tracking of steep rises and drops in the brightness temperature time series (corresponding to the onset and respectively the end of a melt event) for delimiting the melting periods.With this last approach, a pixel is detected as melting for those days included between a successive upward and downward edge pair in the time series (Joshi et al., 2001;Liu et al., 2006).
In this paper, however, we will specifically focus on two algorithms using a constant T19H and XPGR threshold value in Eq. ( 1) over time and space to retrieve the melt extent from the brightness temperatures.Other more complex algorithms using time and/or space varying melt threshold found in the literature are discussed in the Annex B of the Supplement.Knowing that each method has intrinsic limitations and that models compare already well with the simplest approach using a constant melt threshold, the one used here will suffice for a first attempt.In the future, we plan to extend the analysis to other more complex algorithms as well.

Discussion about the use of a constant melt threshold over time and space
Using a spatially and temporally fixed threshold value of T19H thsd in Eq. ( 1) over the entire GrIS to detect the presence of meltwater within the snowpack is the most elementary approach.The threshold value can be estimated either from theoretical considerations or from the comparison between satellite data and ground observations.However, it is clear that the melt threshold in Eq. ( 1) depends in fact on one hand on the melting snowpack behaviors (density, grain size,...) and on the other, on different satellites passing times and sensors used in the SMMR-SSM/I data base (Picard and Fily, 2006).
If the melting snowpack behaviors were known, the best way for determining a snowpack-dependent T19H threshold should be to use a microwave-emission model like Mote (2007).However, even if the microwave scattering coefficients of the snowpack can be derived from the T19H at the beginning of the melting season (Mote, 2007), the grain size/density values used in the Mote's model should be updated as melting and refreezing cycles continue, in view of the constructive metamorphism.Therefore, using a microwave model asks to empirically set several properties of the snowpack (temperature, density and grain size) which induces uncertainties.
like here where we use a constant threshold for each pixel.The limitation of the use of a snowpack independent threshold will be discussed in Sect. 4.
Due to the dependence on which sensor is used, the use of daily mean brightness temperature values might affect retrieval performances for those areas where refreezing occurs at night: if the ascending and descending passes occur around local midnight and noon (as in the case of SMMR), then the mean daily brightness temperature could be lower than that computed using data collected around 06:00 a.m. and 06:00 p.m. (as for the SSM/I) when brightness temperature reaches its maximum (Picard and Fily, 2006).The SSM/I sensors (F08, F11, F13, F17, we will refer here as SSM/I to both SSM/I and SSMIS sensors) indeed collect data over the GrIS in the early morning during descending orbit and in the late afternoon during ascending orbit (there is a shift of 1 or 2 h between the different SSM/I sensors).The difference between the SSM/I and SMMR passing time can be a relevant issue, especially at the beginning and at the end of the melting season, when melting may only last a few hours during the mid to late afternoon.www.the-cryosphere.net/5/359/2011/The Cryosphere, 5, 359-375, 2011 X. Fettweis et al.: Melting trends over the Greenland ice sheet  To reduce the impact of using a constant threshold in Eq. (1) over the whole SMMR-SSM/I covered period, a cross-calibration of the five sensors is carried out here (for example by computing a linear regression of the brightness temperatures over a common period covered by any two sensors, Abdalati et al., 1995).At least one month of common data is available for the several sensors (1987,1991,1995,2009).The five sensors are cross-calibrated to the F8 baseline as in Abdalati and Steffen (2001) and in Liu et al. (2006) and only the GrIS pixels are retained here for the inter-calibration.Adoption of an adaptive T19H-threshold re-computed every year (Zwally and Fiegles, 1994;Torinesi et al., 2003;Picard and Fily, 2006;Tedesco, 2009) or a passing time homogenization of the SMMR-SSM/I data set Picard and Fily (2006) are other more complex solutions to the problem of uncertainties linked to the sensors intercalibration.They are discussed in the Annex B of the Supplement.

T19H-based melt detection
In the following, the simplest algorithm using a constant T19H threshold is applied to the cross-calibrated data.This algorithm will be named T19H melt hereafter.The value of the T19H threshold is calculated by using near-surface temperature measurements performed by AWS of GC-Net between 1995 and 2005 (Steffen and Box, 2001).We assume that melt occurs if the daily mean near-surface temperature is above 0 • C. A different and more sensitive threshold for detecting near-surface melt could be calculated using the maximum 2m temperature above 0 • C, as used by Torinesi et al. (2003).However, melt occurs only if the surface temperature reaches the melting point for a period long enough, whereas a positive maximum 2-m temperature does not necessarily induce surface melt in case of surface inversion or radiative heating of the sensor.A positive daily mean near-surface temperature suggests that the temperature was positive in average over the day, indicating that near-surface melt is likely to occur.
Figure 1 shows time series of near-surface temperature simulated by the two RCMs (MAR and RACMO2) together with GC-Net-measured air temperature for the CP-1, Dye-2, ETH-Camp and JAR1 GC-Net stations for selected years.For each station, the temporal trend of the K-band horizontally polarized brightness temperature is also presented.In addition, we also show the temporal trend of the XPGR, together with a binary plot in which melting was derived using different approaches is indicated by the value 1 and nomelting is indicated by the value 0. The comparison between the near-surface temperature simulated by the RCMs and the measured one at GC-Net stations shows that, over the period 1995-2005, the biases in modeled temperature are very low: MAR (respectively RACMO2) underestimates (resp.over-estimates) the measured near-surface temperature by 0.2 • C (resp.1.3 • C) (see Table 1).
We use two criteria for detecting the melt events in the RCM modeled time series (a list of the algorithms used to retrieve melt can be found in Table 2).The first one, (RCM temp ), uses a modeled positive daily mean 2 mtemperature as melt threshold and compares favorably with the GC-Net temperature-based time series melt detection (see Table 1).The second approach, (RCM melt ), uses a threshold in modeled daily meltwater production as melt events detection and compares better on average with the T19H-derived melt time series, as we will see in the next section.The best agreement between the melt/no melt time series using observed daily mean temperature over the 20 GC-NET AWS listed by Steffen and Box (2001) and the melt time series derived from T19H by using Eq. ( 1) is obtained with T19H thsd =227.5 ± 2.5 K (see Table 3), with more than 95% of melt events detected.Figure 1 shows melt and temperature time series for four GC-Net AWS where considerable melting is recorded.We have selected in Fig. 1 four summers for four different AWS's where no gap is present in the data.
Obviously, it would be ideal to compare melt events detected with the T19H time series with measurements of snow properties, in particular liquid water content (LWC), because the snowpack response to positive air temperatures is complex and melting also depends on other factors, such as albedo, for example.Unfortunately, LWC is not recorded by the GC-Net AWS.Therefore, we use the results from RCMs coupled with a snow model to refine our analysis.
The daily Eq. ( 1)-derived melt signal has been compared for the period 1979-2009 with corresponding time series simulated by the models and applying different melt detection criteria based on a threshold using -surface and near-surface temperature, -internal snowpack temperature, -LWC in the snowpack, -presence of bare ice at the surface (surface snow density higher than 900 kg m −3 ), -meltwater run-off, -meltwater production.
Fig. 2. Time series of the maximum melt extent (in percentage of the GrIS area) simulated by MAR for different melt thresholds: daily meltwater production higher than 1, 2, 5, 7.5 and 10 mm WE and daily runoff higher than 0 mm WE.Finally, the linear trend over 1979-2009 is given by a dashed line.
see hereafter, this suggests that the T19H-based algorithm might not detect melting at the end of the ablation season, when the surface is refreezing but liquid water may still be present in the snowpack or no-melting bare ice may appear on the surface.Unfortunately, no in-situ measurements are available to confirm this hypothesis.We used the same value of daily produced meltwater threshold for both RCMs, which reduced the uncertainty on the relation between T19H thsd and model-simulated daily meltwater production.Figure 2 shows that the choice of the value (listed in the caption) of the melt threshold used in the RCMs could considerably influence the estimated melt extent area while the linear trends over 1979-2009 are similar.The comparison between model-simulated and satellite-derived melt extent is discussed in more detail in the next section.

XPGR-based melt detection
The second algorithm used here to retrieve the melt extent from the passive microwave brightness temperature is based on Abdalati and Steffen (2001) (Abdalati and Steffen, 1997).Bare ice (melting or not) in the ablation zone is also detected as melting by XPGR.
The XPGR threshold varies with the spaceborne passive microwave sensors following Abdalati and Steffen (2001).According to Abdalati and Steffen (2001), the XPGR algorithm is sensitive to (sub-)surface melting and to the presence of liquid water in the snowpack when the snow is refreezing at the surface at the end of the summer.While this approach uses a combination of frequencies compared to T19H melt , Fettweis et al. (2005Fettweis et al. ( , 2006Fettweis et al. ( , 2007) ) showed that XPGR does not detect melt during rainfall events due to the use of the 37 GHz channel.Consequently, Fettweis et al. (2006) proposed an improved version of XPGR (called ImpXPGR), imposing the continuity of the melt season to remove gaps shorter than three days between two melting days.These gaps are found to be associated with dense clouds, mostly causing liquid precipitation on the ice sheet.However, we have not enough in situ data to affirm that the melt season (presence of liquid meltwater) is continuous in time and therefore, this hypothesis needs to be further assessed.
However, knowing that 19 GHz channel is less sensitive to the dense clouds and atmospheric components than the 37 GHz, T19H melt can provide improved performance during these events.That is why we propose here a new approach for improving XPGR by merging the XPGR algorithm of Abdalati and Steffen (2001) and T19H melt to give: In this case, melting is considered to be occurring if either XPGR or T19H is higher than the corresponding threshold values.We point out that, as in Fettweis et al. (2006), 227.5 K is close to the mean T19H temperature (= 225.3 K) plus half of the standard deviation (= 3.1 K) over 1979-2009 and over all pixels of the whole melt area where XPGR detects melting.This last addition to the XPGR algorithm corresponds to Correction #3 of ImpXPGR (see Fettweis et al., 2006 for more details).This improved version of XPGR will be called hereafter ExtXPGR for Extended XPGR.
In order to compare the XPGR-like algorithms with RCM outputs, we use a modeled LWC higher than 1.1% in the top meter of snow as melt threshold value above which the model indicates melting in conjunction with the hypothesis of a surface snow density greater than 900 kg m −3 for detecting bare ice.The value of 1.1% was chosen (instead of 1% as suggested by Abdalati and Steffen, 1997), because it compares better with the satellite-derived data set (see Table 5).It is likely that having a better comparison with 1.1% as melt threshold instead of 1.0% is an artifact from the models.Although this does not affect the comparison significantly, (the 1.0% value of melt threshold is included in the error bars as inferior boundary), this shows the sensitivity of the models to a difference of 0.1% in the LWC.

Melt extent analysis
Results from both models and satellites in Fig. 3 indicate that, on the average, the melting season begins in mid-May, culminates in mid-July (when the highest melt extent occurs on average), and lasts until the end of September (see Fig. 3), when melting is small or negligible.The average maximum melt extent is between 15 and 20% of the GrIS area (∼1.61 × 10 6 km 2 ).The melt extent is larger if retrieved by (Ext)XPGR than by T19H melt because (Ext)XPGR is also sensitive to the presence of LWC into a no-melting snowpack as well as the presence of no-melting bare ice at the surface.The mean correlation between the daily melt extent evolution simulated by the models and that derived from spaceborne microwave remote sensing observations is ∼0.93 with the mean RMSE being ∼3% of the ice sheet surface.
In general, the agreement between model-and satellitederived results is better during the years with low melting (the SMMR years, 1992, 1996 with RMSE ∼2-3%) than during those years with high melting (2000's with RMSE ∼3-4%), whereas the RCMs-based time series can be considered as quasi-homogeneous.The RCMs are forced by the ERA-40 reanalysis (1977ERA-40 reanalysis ( -2002) ) and after that by the ERA-Interim reanalysis (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009) from the ECMWF.A potential issue could be that the passive microwave data are used in the reanalysis data (notably to prescribe the sea ice cover) which drives the RCM that are compared to the passive microwave data afterwards.However, the passive microwave data drives the ECMWF reanalysis only in surface and the surface conditions simulated by the RCM over the GrIS are quasi-independent of the surface conditions The Cryosphere, 5, 359-375, 2011 www.the-cryosphere.net/5/359/2011/2. The statistics (i.e. the mean correlation coefficient as well as the mean Root Mean Square Error over 1979Error over -2009) ) between MAR (in blue)/RACMO2 (in red) and the T19Hmelt (resp.ExtXPGR) algorithms are also listed.2).The error bars show a positive/negative change of the threshold used to derive the melt extent as listed in Table 2.The statistics (i.e. the mean correlation coefficient as well as the mean Root Mean Square Error over 1979Error over -2009) ) between MAR (in blue)/RACMO2 (in red) and the T19H melt (resp.ExtXPGR) algorithms are also listed.
coming from the reanalysis over the ice sheet.According to Hanna et al. (2009), the impacts of SST changes are for example negligible in MAR over the GrIS.Therefore, the use (or not) of passive microwave data in reanalysis does not significantly influence the comparison of the RCMs with the satellite-derived data.
For the reader's convenience, comparisons between satellite and RCMs for each summer between 1979 and 2009 are reported in the Annex C of the Supplement.We can notably observe a large difference between the melt extent estimated by T19H melt and by (Ext)XPGR in September 1982September , 1983September , 1985September , 1987September , 1990September , 1994September , 1999September , 2004September and 2006 when T19H melt seems, as assumed earlier, not to be sensitive to the presence of non-melting bare ice at the surface or internal liquid water into the snow pack, differently from (Ext)XPGR, the latter being in agreement with the results from the models.
Looking at the number of melting days/summer, it appears that areas with more than 100 melting days/summer in average lie below 70 • N along the ice sheet margin and at the northeast of the ice sheet (Fig. 4).The melt season lasts generally longer in the Southwest with more than four months of melt days detected by ExtXPGR, on average.Over the whole ice sheet, the average difference between the total mean num-ber of melt days by summer simulated by the models and retrieved from satellite data with T19H melt (resp.ExtXPGR) is 2.5 days (resp.4.5 days).In 95% of the cases, there is agreement among pixel-by-pixel and day-by-day cases between RCMs and satellite-derived melt detected events.
The disagreement between satellite and models can be due to (i) limitations in the algorithms used to derive the melt extent from the microwave data set but also to (ii) biases in the models.The error bars show that using a melt threshold slightly different on both models and on the satellite algorithm does not considerably affect the comparison.However, because of the lack of in-situ measurements, it is difficult to clearly identify (i) where and when the satellite estimates are wrong and (ii) when and why the models fail, knowing that the biases may depend on several errors or feedback.

Some limitations in the remote sensing algorithms
A well-known limitation of spaceborne microwave data is the large field of view of the sensors (being of the order of several tens of km, depending on the frequency).This may be responsible for biases along the ice sheet margin where the microwave signal could be affected by the tundra or sea surrounding the ice or when the ablation zone is very narrow.
www.the-cryosphere.net/5/359/2011/The Cryosphere, 5, 359-375, 2011 Fig. 4. Top -Annual mean total number of melt days derived from the spaceborne passive microwave data and simulated by the RCMs using algorithms described in Table 2. Below -The difference between models and the T19Hmelt and ExtXPGR algorithms.The mean number of GrIS pixels when RCM and the algorithms detect melt (RCM = SAT), when RCM detects melt but the retrieving algorithms do not (RCM > SAT) and when RCM does not detect melt while the algorithms do (RCM < SAT) is also listed as a percentage of the number of GrIS pixels × summer days.

Fig. 4.
Top -annual mean total number of melt days derived from the passive microwave data and simulated by the RCMs using algorithms described in Table 2. Below -the difference between models and the T19H melt and ExtXPGR algorithms.The mean number of GrIS pixels when RCM and the algorithms detect melt (RCM = SAT), when RCM detects melt but the retrieving algorithms do not (RCM > SAT) and when RCM does not detect melt while the algorithms do, (RCM < SAT) is also listed as a percentage of the number of GrIS pixels × summer days.
The increase in microwave emissivity measured by the radiometer when melting occurs may be significantly reduced if a pixel combines both melting and no-melting areas.As considerable melting occurs in these areas, we decided not to remove these pixels from the satellite-models comparison.
We only removed pixels near open water (ocean, fjord) as explained in Sect. 2. The XPGR algorithm detects fewer melt days (about 10%) compared to ExtXPGR (see Fig. 4).This difference is clearly visible along the ice sheet margin where the probability that XPGR is biased by rainfall events or by the presence of clouds with liquid water is high.Figure 1 and Table 3 show that XPGR is less adequate than T19H melt or ExtXPGR in detecting melt in the low elevation part of the ablation zone where the probability of rainfall events is the highest and where the probability of melt is the highest, too.Finally, Fig. 3 shows that the differences between XPGR and Ex-tXPGR occur mainly before August and decrease after mid-August when most precipitation falls as snow.
Both RCMs simulate fewer melt days (about 20 days) along the eastern and southeastern mountainous regions of the ice sheet than the microwave-derived estimates (see Fig. 4).Possibly, microwave brightness temperatures could be biased by rock outcrops found in these regions in the Antarctic Peninsula and therefore, the remote sensing algorithms may overestimate the occurrence of melt events over this area, as suggested by Torinesi et al. (2003).In addition, this apparent RCM overestimation in respect to the microwave-derived one could also be an artifact induced by the high accumulation rates found in this region, knowing that high accumulation rates induce increases of brightness temperature (Bolzan and Jezek, 2000;Winebrenner et al., 2001).And then, as suggested by the RCMs (see Fig. 5), higher T19H thsd than 227.5 K should be used in this area to detect the melt events from T19H.Value of T19Hthsd (given here as anomaly in K in respect to 227.5K) which allows the best comparison between the melt extent derived from T19Hmelt and simulated by the RCM (where a daily meltwater production higher than 8.25mm is used as melt threshold).The GC-Net locations quoted in the text are also shown.

Fig. 5.
Value of T19H thsd (given here as anomaly in K in respect to 227.5 K) which allows the best comparison between the melt extent derived from T19H melt and simulated by the RCM (where a daily meltwater production higher than 8.25 mm is used as melt threshold).The GC-Net locations quoted in the text are also shown.
improvement of the comparison between RCMs and satellite and more generally, the RCMs agree to suggest that differences in accumulation rate, snow density,... influence the T19H thsd value.Therefore, using a constant T19H thsd value over each pixel, as assumed in this study, can explain in part the difference between models and satellite.However, it should be noted that biases/discrepancies in both RCMs impact also the T19H thsd values computed in Fig. 5.For example, we can see this impact along the western margin where MAR and RACMO2 suggest different values of T19H thsd whereas the same daily meltwater production threshold is used in both RCMs for detecting melt events.Unfortunately, most of the GC-net AWSs used to calibrate the T19H thsd ) value are situated in the percolation zone where 227.5 K is the best value for comparing RCM and satellite melt extent.Only one station is located in the ablation zone: JAR1, where the GC-net measurements suggest 220-225 K as the melt threshold (see Table 3) in agreement with the models (see Fig. 5).Therefore, more in-situ observations in the ablation zone and along the southeastern mountainous regions are needed to confirm the T19H thsd values suggested by the RCMs.

Biases in the models
Maximum melt extent area for the melt threshold using the meltwater produced by day (resp.LWC) occurs, on the average, during mid-July.The timing of this maximum is well reproduced by both models (see Fig. 3).However, from May to mid-June (resp. in July), RACMO2 underestimates (resp.overestimates) the melt area with respect to both satellite-based melt retrieval algorithms.The largest bias in MAR occurs in July when MAR underestimates the melt area if the meltwater production is used as the melt threshold in the passive microwave data.As we will see, discrepancies in the modeled surface radiation balance most likely explain these disagreements between both RCMs. Figure 6 shows that MAR, compared to RACMO2, predicts more incoming longwave radiation in winter and less in summer.Compared to the K-Transect measurements (West Greenland), RACMO2 tends to underestimate the longwave radiation mainly in winter (Ettema et al., 2010).This suggests that MAR underestimates even more the longwave radiation in summer.Shortwave incoming fluxes are higher (∼+5-10 W m −2 ) in RACMO2 than in MAR before June and lower afterwards.The differences are however negligible compared to the absolute values of the longwave fluxes driving the seasonal variability of the net fluxes.Finally, www.the-cryosphere.net/5/359/2011/The Cryosphere, 5, 359-375, 2011 the surface snow albedo decreases earlier in MAR than in RACMO2 at the end of spring and reaches the dry snow albedo earlier as suggested by Torinesi et al. (2003) at the beginning of the winter, because the density of the snow surface (which prescribes the albedo in RACMO2) varies more slowly than the surface snow grains morphological characteristics (which determine the albedo in MAR) (Lefebre et al., 2003).However, it should be noted that part of differences between the MAR and RACMO2 values of albedo are enhanced by the melt-albedo positive feedback.Knowing that RACMO2 underestimates the satellitederived melt extent in May-June unlike MAR, we can suggest that at the beginning of the melt season, the snow albedo parametrization used in RACMO2 is not sensitive enough to wet snow conditions compared to the one used in MAR which simulates a lower albedo than RACMO2 during this period (see Fig. 6).In addition to differences in the surface albedo parametrization, the MAR infrared radiative surplus in winter also favors an earlier beginning of the ablation season because the snowpack is warmer in MAR at the end of spring than in RACMO2.
In July, RACMO2 overestimates (note that it is the opposite for MAR) the melt extent compared to the satellites, likely because the RACMO2 albedo is too low as already suggested in Ettema et al. (2010).The downward longwave underestimation in MAR probably explains the melt extent underestimation in July.
Differences in the response time between the two albedo parametrization schemes were shown by Lefebre et al. (2003), concluding that the CROCUS albedo parametrization gives better results in MAR (trough comparison with measurements at ETH-Camp).In addition, the comparison with GC-net measurements in Table 1 shows that MAR underestimates slightly the near-surface temperature in summer and RACMO2 overestimates surface temperature with respect to GC-net measurements as also shown by Ettema et al. (2010).Therefore, it is hard to identify which modules in the models fail.Both albedo parametrization and radiation schemes (based in MAR and RACMO2 on the ECMWF model) should be critically assessed using high quality measurements of the radiative budget of the atmosphere.At the end of the melting season (mid-August to September), RACMO2 underestimates the melt area only if the LWC is used as melt threshold (Fig. 3).The LWC in RACMO2 is limited to 2% of open pore space, whereas this is 6% in the MAR model.Excess liquid water runs directly off in both models.Therefore, at the end of the melting season, there is still liquid water in the MAR snowpack while most of the meltwater is already refrozen in RACMO2.This could be the reason why RACMO2 underestimates the melt extent in August if LWC is used as melt threshold and also why the MAR snowpack is warmer (refreezing of the retained meltwater prevents initial winter cooling of the snow pack).However, a LWC maximum reaching 6% in MAR could be also too high as we can see after an early (June 1979(June , 1988) ) or late (September 1999(September , 2004) melt event (see the Supplement) when MAR overestimates the area where liquid water is present in the snowpack while no more melt is simulated and observed (see T19H melt ) at the surface.
The patterns of differences between the RCMs and T19H melt are quite similar (see Fig. 4).Simulated melt extent is underestimated with respect to T19H melt along the southeastern margin and overestimated over the rest of the melt area.The differences are small compared to the differences with ExtXPGR.The main difference between MAR and RACMO2 occurs along the northwestern coast, where MAR overestimates the melt extent while RACMO2 underestimates it with respect to T19H melt .The differences between the outputs of the two RCMs are further discussed in the following paragraph.
If we use the LWC as the melt threshold, MAR overestimates the melt extent along the northwestern coast and in South Greenland and underestimates it in West and North Greenland with respect to ExtXPGR.The RACMO2 biases compared to ExtXPGR are very different, with an overestimation in the North and in the interior of South Greenland.These discrepancies are mainly due to differences in the snow model and in the snowfall pattern, affecting the number of days when bare ice is present on the surface.
1.As discussed earlier, MAR identifies melting areas that are absent in RACMO2 after Mid-August, likely because the maximum LWC in MAR could be 3 times larger that in RACMO2.This also explains why MAR and RACMO2 outputs differ in South Greenland.
2. In addition to the presence of LWC, a pixel is considered as melting by XPGR if there is bare ice at the surface, even if this pixel is not melting.In the models, a pixel is assumed to be covered by bare ice if the density of the first layer in the snow model is higher than 900 kg m −3 .Figure 7 shows that a large part of the discrepancies between MAR and RACMO2 in Fig. 4 are explained by the difference in the modeled number of pixels with bare ice at the surface.The appearance of bare ice in summer is conditioned mainly by winter snow accumulation (Fig. 7).Indeed, MAR simulates more pixels with bare ice in areas of the ablation zone along the west coast, where MAR simulates less solid precipitation than RACMO2.Summer snowfall www.the-cryosphere.net/5/359/2011/The Cryosphere, 5, 359-375, 2011  .The cumulated melt area is defined as the annual total sum of every daily ice sheet melt area.The linear trend (dashed line) as well as the error bars (see Table 2) are also shown.The uncertainty range of the trend denotes three standard deviations of the trend i.e. a significance of 99%.
also impacts melting by increasing the albedo for a short period after a snowfall event.Differences in snowfall likely explains why RACMO2 underestimates melt extent (compared to T19H melt ) along the northwestern ice sheet margin.
3. In the North, both modeled precipitation distributions are similar, although MAR simulates fewer bare ice pixels.However, the smallest thickness for a new top layer in the snow model is 2 mm water equivalent (WE) in MAR and 6.5 cm WE in RACMO2.Therefore, RACMO2 needs 6.5 cm WE of fresh snow to create a new layer above bare ice, while a few mm of WE is enough in MAR to change a bare ice (and then melting) pixel to a non-melting pixel.This difference in the minimum layer thickness impacts the results mainly in the dry North where only a few snowfall events occur in summer.This explains also, in addition to the differences in the surface albedo parametrization, why the surface albedo reaches dry snow values at the end of the ablation season sooner in MAR than in RACMO2.

Interannual variability
According to previous studies (Fettweis et al., 2007;Mote, 2007;Tedesco et al., 2008), both microwave-derived and model-simulated melt extent show a significant (99% confidence level) positive trend of melt extent (Fig. 8).Since 1979, the annual cumulated melt extent area (i.e. the annual total sum of every daily ice sheet melt area) has doubled.Figure 9 shows that positive trends in melt extent occur everywhere in both the ablation and percolation zones.The slope of the trend in the case of the XPGR-based melt area following Abdalati and Steffen (2001) is not significant (see Fig. 8).This is most likely because rainfall and clouds with liquid water (which perturbs the melt detection by XPGR) also increased in summer over the period 1979-2009(see Fettweis et al., 2006 and Annex D in the Supplement).
According to both models and remote sensing algorithms, the maximum cumulated melt area occurred in 2007, followed by 1998.Low melt extent values occurred in 1983 and 1992 after, respectively, the El Chichon and Mont Pinatubo eruptions.The agreement between RCMs and microwavederived melt over the full spaceborne passive microwave period (see the Supplement) suggests that (i) using a linear  regression to cross-calibrate the five sensors is sufficient to apply the same threshold throughout the whole period in Eq. ( 1), (ii) the modeled melt variability is good and then, (iii) that the RCMs can be used to estimate the melt extent prior to 1979.The RCMs show in Fig. 8  5. there is a periodicity of 2-3 yr (about 30 months) in the time series, i.e. a summer with a high cumulated melt extent is mostly followed by a summer with less melt.A similar variability is found in the North Atlantic Oscillation (Nicolay et al., 2008), which is known to influence the summer temperature over the GrIS (Fettweis, 2007).

Discussion and conclusion
The results from two RCMs have been compared with microwave brightness temperature-based melt extent estimates with the objective of studying the evolution of melting over the GrIS since 1958.We selected two algorithms to retrieve melt extent from the spaceborne passive microwave data set.The first one, using a fixed threshold on the T19H brightness temperature (see Eq. 1), is especially sensitive to the production of surface and sub-surface meltwater by day.The second one, based on the XPGR algorithm from Abdalati and Steffen (1997) (see Eq. 3), detects melt when liquid water is present in the snow pack.The time evolution of melt area over the GrIS through summer and the number of melt days compare very well with the results from both models.The retrieving algorithms as well as the RCMs show a significant (p-value = 0.01) increase of the melt area over the GrIS since 1979.The models show that the recent cumulated melt extents are without precedent in the 50 last years and are two times larger that the ones from the 1980's.
The close agreement between remote sensing results and model outputs provides robustness to the melt retrieving algorithms and melt threshold used.This intercomparison also allowed us to highlight some biases in the models: 1. Differences in response time to surface snow changes in the albedo parametrizations used by RACMO2 explains likely why RACMO2 underestimates the melt extent at the beginning of the melt season and overestimates it in July.
2. The longwave radiation should be increased in MAR in summer to correct the melt extent underestimation occurring in July.
X. Fettweis et al.: Melting trends over the Greenland ice sheet  3. The LWC in the RACMO2 snow model is currently limited to 2% (compared to 6% in the MAR model) which prevents melt detection at the end of the ablation season when freezing conditions occur at the surface and there is normally still liquid meltwater in the snow pack.In addition, a maximum LWC value of 6% in MAR likely overestimates the melt extent after high melt events.High quality snow observations at several locations on the ice sheet, helped to select the maximum LWC value in the snow pack.
4. The minimum thickness of a new layer in the RACMO2 snow model is 6.5 cm WE, compared to 2 mm WE in the MAR model.This overestimates in RACMO2 the pixels with bare ice in the ablation zones where few snowfall events occur in summer.But in MAR, knowing that it underestimates the melt extent at the north of the ice sheet, the minimum layer thickness threshold should be increased to enhance the agreement with remote sensing data.
5. The comparison with the LWC-sensitive algorithms shows that there are some biases in both simulated precipitation patterns.The winter snowfall accumulation impacts the appearance of bare ice at the surface in summer which is detected as wet by the XPGR-based algorithms.
However, part of the discrepancies between satellite and models comes from the microwave data itself which is biased by the tundra or rock pixels nearby the ice sheet or by use of a constant T19H melt threshold over each pixel, whereas the RCMs suggest that snowpack behaviors may significantly influence the T19H thsd value.More in-situ measurements are needed to support these hypothesises.
The close agreement between the RCM results and the microwave-derived melt extent suggests that the RCMs could be run in an assimilation mode, constrained by the spaceborne passive microwave data.Indeed, the snowpack properties of the pixels where RCM and satellite disagree could be adjusted.For example, if the RCM has to detect melt in respect to the satellite observations, but does not simulate it, the snowpack temperature could be increased to reach conditions more favorable to melt while the water mass is conserved.As a future development, we could associate a meltwater threshold as in Table 4 to each brightness temperature and assimilate it in the model to improve the meltwater production simulation.This could have significant impact on the modeled meltwater run-off.
Moreover, this study has shown that large uncertainties remain in the modeled bare ice appearance at the surface while most of the meltwater run-off occurs on this area.Therefore, the MODIS-based albedo could be used for evaluating the bare ice zones in the models and afterwards, for forcing the models in an assimilation mode.Knowing that the GrIS run-off variability is driven to a large part by the bare ice area variability, this assimilation should help to improve the matching with other satellite data sets (GRACE,...), with the objective of reducing the uncertainties of the SMB modelbased estimates.

16
Fig. 1.(a) Top -Daily mean near-surface temperature (in black) observed in summer 2002 at Crawford Point 1 (69.9 • N, 47• W, 2022 m), 3 m-temperature simulated by the MAR model (in blue), 2 m-temperature simulated by the RACMO2 model (in red) and, the T19H brightness temperature (in green on the left axe) for the pixel nearest CP1.Middle -the XPGR value as defined byAbdalati and Steffen (2001).The dotted line shows the XPGR threshold value used byAbdalati and Steffen (2001) for melt detection.Below -The Melt/no melt time series derived from observation (daily mean temperature ≥0 • C), simulated by the RCMs (see Table2), derived from the T19H temperature by using Eq. 1 and using the XPGR algorithm.(b) The same as (a) but for DYE-2 (66.5 • N, 46.3 • W, 2165 m) in 2005.(c) The same as (a) but for ETH-Camp (69.6 • N, 49.2 • W, 1149 m) in 1998.(d) The same as (a) but for JAR-1 (69.5 • N, 49.6 • W, 962 m) in 1999.

Fig. 1 .
Fig. 1.(a) Top -daily mean near-surface temperature (in black) observed in summer 2002 at Crawford Point 1 (69.9 • N, 47• W, 2022 m), 3 m-temperature simulated by the MAR model (in blue), 2 m-temperature simulated by the RACMO2 model (in red) and, the T19H brightness temperature (in green on the left axe) for the pixel nearest CP1.Middle -the XPGR value as defined byAbdalati and Steffen (2001).The dotted line shows the XPGR threshold value used byAbdalati and Steffen (2001) for melt detection.Below -the Melt/no melt time series derived from observation (daily mean temperature ≥0 • C), simulated by the RCMs (see Table2), derived from the T19H temperature by using Eq. 1 and using the XPGR algorithm.(b) The same as (a) but for DYE-2 (66.5 • N, 46.3 • W, 2165 m) in 2005.(c) The same as (a) but for ETH-Camp (69.6 • N, 49.2 • W, 1149 m) in 1998.(d) The same as (a) but for JAR-1 (69.5 • N, 49.6 • W, 962 m) in 1999.

Fig. 3 .
Fig.3.Mean seasonal cycle of melt area (in % of GrIS area) simulated by MAR and RACMO2 as well as retrieved from the spaceborne passive microwave data set with the different techniques discussed in the text (see 2).The error bars show a positive/negative change of the threshold used to derive the melt extent as listed in Table2.The statistics (i.e. the mean correlation coefficient as well as the mean Root Mean SquareError over  1979Error over   -2009) )  between MAR (in blue)/RACMO2 (in red) and the T19Hmelt (resp.ExtXPGR) algorithms are also listed.

Fig. 3 .
Fig. 3. Mean seasonal cycle of melt area (in % of GrIS area) simulated by MAR and RACMO2 as well as retrieved from the spaceborne passive microwave data set with the different techniques discussed in the text (see Table2).The error bars show a positive/negative change of the threshold used to derive the melt extent as listed in Table2.The statistics (i.e. the mean correlation coefficient as well as the mean Root Mean SquareError over 1979Error over  -2009) )  between MAR (in blue)/RACMO2 (in red) and the T19H melt (resp.ExtXPGR) algorithms are also listed.

Figure 5
Figure5suggests also that in the ablation zone where bare ice appears in summer (i.e.area with high snow density), lower values of T19H thsd than 227.5 K will allow the

Fig. 6 .
Fig.6.Top -time series of the mean shortwave (solid line), longwave (dotted line), albedo (dashed line on right y-axis) and net radiation (dashed line) averaged over the melting zone of the GrIS (i.e.surface height below 2500 m) and over the period 1979-2009 simulated by MAR (in blue) and by RACMO2 (in red).Below) Time series of the mean difference in surface temperature (resp.the temperature of the top meter of snow) between MAR and RACMO2.

Fig. 7 .
Fig. 7. Left -The mean difference of the number of days (per year) with bare ice at the surface simulated by MAR and by RACMO2.Right -The same but for the annual snowfall in mmWE.

Fig. 7 .
Fig. 7. Left -the mean difference in the number of days (per year) with bare ice at the surface simulated by MAR and by RACMO2.Right -the same but for the annual snowfall in mm WE.

Fig. 8 .
Fig. 8. Time evolution of the annual cumulated GrIS melt area simulated by the RCMs and retrieved from the spaceborne passive microwave data set with the different techniques discussed in the text (see 2).The cumulated melt area is defined as the annual total sum of every daily ice sheet melt area.The linear trend (in dashed line) as well as the error bars (see Table2) are also shown.The uncertainty range of the trend denotes three standard deviations of the trend i.e. a significance of 99%.

Fig. 8 .
Fig.8.Time evolution of the annual cumulated GrIS melt area simulated by the RCMs and retrieved from the spaceborne passive microwave data set with the different techniques discussed in the text (see 2).The cumulated melt area is defined as the annual total sum of every daily ice sheet melt area.The linear trend (dashed line) as well as the error bars (see Table2) are also shown.The uncertainty range of the trend denotes three standard deviations of the trend i.e. a significance of 99%.

Fig. 9 .
Fig. 9. Linear melt trend (in ablation days year −1 ) detected by T19Hmelt and simulated by the models for the period 1979-2009.

Fig. 9 .
Fig. 9. Linear melt trend (in ablation days yr −1 ) detected by T19H melt and simulated by the models for the period 1979-2009.

Table 1 .
Mean (over the 1995Mean (over the  -2005 summers)  summers)percentage of number of melt days at 10 GC-Net AWS's (where melt is observed) when the RCM and the observations are in agreement, when the RCM fails to detect melt and when RCM detects melt while the observations do not.Two approaches, RCM-temp and RCM-melt (explained in Table2), are used for detecting melt in the RCM modeled time series.The location of the GC-Net AWS's is also listed.Finally, the mean modeled daily temperature bias (in • C) with respect to the GC-net measurements over the 1995-2005 summers is also given.

Table 2 .
Abbreviation of the different techniques discussed in the text to detect the melt.The errors bar range is also listed.

Table 3 .
The same as Table1but with the T19H-based algorithms for different values of T19H thsd .The best agreement is written in bold.

Table 4 .
Correspondence between Eq. (1)-based and RCM-based melt detection.The mean number of GrIS pixels when both RCM and T19H-based algorithms detect melt (RCM = SAT), when RCM detects melt and not the retrieving algorithms (RCM > SAT) and when RCM does not detect melt and well the T19H-based algorithms (RCM < SAT) is also listed in percent over the number of GrIS pixels × summer days.The remaining percent correspond to days where both RCM and T19H-based algorithms do not detect melt.The statistics (i.e. the mean correlation coefficient (Corr.)aswellas the mean Root Mean Square Error(RMSE) over 1979(RMSE) over  -2009) )between MAR/RACMO2 and the T19H-based GrIS melt extent time series are also listed.