The terrain effect correction : how it works

=>


Caution, this post contains formulas.

 

The topographic (or terrain) effects on the observed reflectances are due to several phenomena, illustrated below :

  • the closer the surface is perpendicular to sun direction, the more energy it receives per surface unit (we talk about irradiance). If the surface is parallel to sun direction, it does not receive direct sunlight. We can model it this way :
    • For an horizontal surface :  E_h= E_0.T_{dir}^\downarrow.cos(\theta_s)

      Definition of sun zenith angle and sun incidence angle

    • For a sloped surface  E_i=E_0.T_{dir}^\downarrow.cos(\theta_i)
    •  E_0 is the Top of Atmosphere irradiance, and  T_{dir}^\downarrow is the downward direct transmission, i.e. the proportion of the light that reaches directly the surface without being scattered by the atmosphere.
    • Assuming that all the irradiance is direct, the measured reflectance if the surface was horizontal is calculated from the following formula:  \rho_h=\rho_i \frac{cos(\theta_s)}{cos(\theta_i)} . However, the above assumption is not true and this formula tends to over correct terrain effects
  • The surfaces also receive a diffuse sun irradiance scattered by the atmosphere. If the surface is not horizontal, a part of the sky is obscured by the slope reducing the diffuse irradiance. Moreover, the diffuse irradiance depends on the amount of aerosols (and clouds) in the atmosphere. In addition, the surrounding terrain can also hide a part of the sky, but we do not take this effect into account here in our modelling. We use the following approximation, which is equivalent to assuming that the slope is alone in a horizontal region.
    • If surface is horizontal, the visible sky fraction is 1, if it is vertical, this fraction is 1/2
    •  \displaystyle F_{sky}= \frac{1+cos(slope)}{2}
  • Finally, the slope can receive light from surrounding surfaces, which become directly visible. In our simplified model, we always assume the entire environment of our slope is flat and,for instance, we do not take the effect of the opposite side in a valley into account :
    • If surface is horizontal, the visible ground fraction is 0, if it is vertical, it is 1/2.
    •  \displaystyle F_{fround}= \frac{1-cos(slope)}{2}

 

Finally, we use the following formula  to compute the reflectance that would be observed if the surface was horizontal   \rho_{h} , as a function of the slope (inclined) reflectance   \rho_{i} :

  \rho_{h}=\displaystyle \rho_{i}.\frac{T^{\downarrow}}{T_{dir}^{\downarrow}.\frac{cos(\theta_i)}{cos(\theta_s)} + T_{dif}^{\downarrow} F_{sky} + T^{\downarrow} F_{ground} \rho_{env}}

 T^{\downarrow} is the downward transmission, sum of direct and diffuse irradiances :  T^{\downarrow}= T_{dir}^{\downarrow}+ T_{dif}^{\downarrow} , and \rho_{env} is the average reflectance of the neighbourhood.

 

Finally, we can also account for bidirectional reflectance effects, but this correction is tricky since directional effects depend on the surface cover type. See for instance : Dymond, J.R.; Shepherd, J.D. 1999: Correction of the topographic effect in remote sensing. IEEE Trans. Geosci. Remote Sens. 37(5): 2618-2620.

 

It is very difficult to validate the correction of directional effects : , we could compare the correction results for satellite overpasses at different times in the day. But all the satellite optical imagers have nearly the same overpass time. A qualitative way of estimating the accuracy is to check that similar land covers on opposite slopes in a valley ( a meadow, a forest) have a similar reflectance after correction. The most suitable points are North-South valleys. Here are some examples of terrain effects correction results.

Formosat-2 image in the Alps, before (left) and after (right) terrain effect correction. The image on the right seems to have been flattened, and the opposite slope reflectances seem much more alike after correction.

Finally, an essential part of the method's accuracy is the availability of a highly accurate digital  elevation model (DEM),  up to now, only the SRTM DEM is available globally, and it only has a 90 meter resolution. Its accuracy is somewhat inadequate and sometimes leaves artefacts if the slope changes are poorly located.

 

 

Just back from Living Planet Symposium

I had the opportunity to participate to the ESA Living Planet Symposium, in Edinburgh, last week, together with 1800 other remote sensing data user, mostly from Europe. It is incredible to see that after more than 20 years in this domain, I know less than 10% of the audience !

There were 10 parallel sessions each day, and in the biggest room, 2 sessions were organised for each Sentinel mission, and these sessions have been filmed. These films are on the web, accessible from this page. As a reader of this blog, you might be interested by the sessions about Sentinel-2, available here :

  • Sentinel-2 mission, to have a complete update about the mission and its ground segment (official launch dates : 09/2014 (S2-A) and 09/2015 (S2-B)).
  • Sentinel-2 potential applications and services, in which I presented the SPOT(Take5) experiment (the third presentation, from 46' to 1h06'). Please be lenient, I was a little nervous in front of 150-200 people !

I also wrote (quickly) a seven page paper about the SPOT4 (Take5) experiment.

Feedback on the irrigation scheduling experiment using remote sensing images

=>

CESBIO contributes to an international joint laboratory in Morocco, called TREMA, "Télédétection et Ressources en Eau en Méditerranée semi-Aride", which means "Remote Sensing and Water Resources in Semi-Arid Mediterranean". This year, this laboratory has embarked on an ambitious experiment of irrigation scheduling by satellite imagery, on a wheat plot near Marrakech. This experiment was already described in March, and it gave very promising results.


The main objective of the experiment was to see if  the logistics of irrigation scheduling by water balance model were feasible in real conditions. For this, a farmer accepted to play our game on two four hectares plots of wheat: Irrigation of the reference plot was driven by the farmer in the usual way. The test plot irrigations was driven by our tool SAMIR (FAO-56 model forced by satellite imagery).

 

Since the sowing late December to the harvest in early June, a weather station installed on a reference culture has given us the daily reference evapo-transpiration measurements. On the other hand, to control a posteriori the quality of our estimates of water requirements for irrigation, two flux measuring stations were set up. We also acquired a series of images SPOT5 early in the season to compensate for the slightly late start of SPOT4 experience (TAKE5) which began in February.

 

In addition to a clear weather throughout the season, we were able to benefit from the excellent work of the SPOT4 (TAKE5) team which provided us with the georeferenced images very quickly. The NDVI evolutions were thus available in a relatively short time. As an end user, the Office of Agricultural Haouz allowed us to perform the irrigation of the test plot in the best conditions while being subjected to the constraints of the canal system.

On the ground, everything did not work as well as we planned. Following a misunderstanding with the farmer, we completely missed the second irrigation and the fertilizer application was not timely. Indeed, the study plot is installed on a heavy clay soil that forms a crust. We were not aware that, a few days after sowing, a specific irrigation is needed to ease the emergence of plants. On the other hand, the farmer applied nitrogen fertilizer on two plots just after irrigation of the reference parcel and relatively far from the irrigation of the test plot. Under these conditions the nitrogen is relatively less soluble, and our test plot lacked fertilizers.

Our experiment has been seriously hampered by the misunderstandings with the farmer. But despite the bad start, the experiment was pursued to its end.

 

This plot shows the changes throughout the course of the experiment of the water supply from rainfall and irrigation, the evapo-transpiration ETobs measured in the field and the Evapo-Transpiration ET estimated by SAMIR model, using the vegetation status from SPOT4 (Take5) images. On this plot, the dates of irrigation were suggested by the model.

 

To our surprise, the results are extremely promising. Indeed, despite a 20% lower biomass compared to the plot driven by the farmer, we got a equivalent performance in grain yield. This can be explained by the fact that, although the average number of wheat blades was much lower on the test plot, it is very likely that the reference plot, irrigated by the traditional method, has suffered water stress in late March limiting the filling of grain.

 

 

This full-scale experiment finally turned out to be very instructive. First,  imaging/weather/irrigation logistics worked great : the weather data transmission, the reception and the geometric and radiometric correction of images, the model runs and  irrigation decision were largely automated. The SPOT4 (Take5) data, that prefigure those of Sentinel-2, proved perfectly suited to this application. Unfortunately, the clay crust has severely limited the emergence of culture. Yet this phenomenon, well-known to our farmer, taught us to cultivate humility ;-) , and we will consider the introduction of the risk in a decision support system. Finally, the functional constraints of the gravity irrigation system have taught us that our tool should be more flexible to recommend an irrigation period instead of a single date, and that we should link the service to weather forecasts.

 

Following this experiment, we started developing a Web service (SAT-IRR) that should shortly provide the essential functions of an irrigation decision support with a simplified interface.

 

A SPOT4(Take5) meeting at CNES, on October the 2nd.

=>

CNES is organizing a one day meeting about SPOT4(Take5) on October the 2nd, in its premises in Toulouse.

 

This meeting is aimed at explaining how the experiment turned out and at presenting the available data set. These presentations wil take place during the morning and the afternoon will be dedicated to presentations of users projects and eventually their first results. Questions/answers sessions will also be organized.

 

Because of CNES security rules, this meeting is only open to persons having a passport from one of the European Union countries. Otherwise a 2 month delay is necessary, which is no possible here.If you are interested, you will need to send an email before September the 17th, to sylvia.sylvander @ cnes.fr (mentioning your name, nationality and affiliation). If you would like to present your project with SPOT4(Take5) data, please also mention the title of your presentation.

 

A new version of the SPOT4(Take5) products is available.

Here are the thumbnails from the China(2) site, for which several dates were missing on the version 1.0. Please note that on the server, you may download all the dates at once by clicking on the 1C or 2A buttons.

=>

The CNES teams of the THEIA Land Data Center have reprocessed the SPOT4 (Take5) data, in order to take into account a large number of images that were not processed in the first place, because some data had not been yet received or because their processing had failed due to a few little bugs.

 

The same processors and parameters were used and the only difference is the increased number of available dates, but as the L2A methods are multi-temporal and recurrent, when we add an image, the results on the subsequent images are also changed. It is thus advisable that you download again all the products of the sites you are interested in, from the following address : http://spirit.cnes.fr/take5

 

On this prototype ground segment, our management of product versions is basic, and only takes the processors into account. As the processors are unchanged, the new version 1.1 products are still identified as level 1.0 products in the Metadata. We are sorry for this inconvenience, you will need to pay attention not to mix them with the older version.

 

A few missing images

=>

I just took a work break in the middle of my holidays, but as I was away, we received a few feedbacks from users, and CNES PTSC teams, with Mireille's help at CESBIO verified the data sets released on July the 16th, in quite a rush...

 

They found out that a few scenes were missing. For some of them, it was due to the late arrival of some images (just as for planes at the airport). These images have already been added to the server.

And there were a couple of bugs that mostly affected the sites made of several SPOT images (CNES and NASA sites), and ESA Chinese site. These glitches have been corrected and the reprocessing started. The whole data set will be updated before end of August, which will constitute the version 1.1 of the SPOT4(Take5) data set.

 

Keep posted on this blog, we will update it as soon as the data are available. Meanwhile, version 1 is still accessible here, and the format described there.

 

The adjacency effects, how they work.

As explained in the post about atmospheric effects, the scattering of light by molecules and aerosols in the atmosphere brings about several effects : scattering adds some haze on the images (the atmospheric reflectance), lessens the signal from the surface (the atmospheric transmission), and blurs the images (the adjacency effects). This post is about the adjacency effects, the other aspects have already been quickly explained in the above post.

 

The figure on the right shows the types of paths that light can follow before getting to the satellite. Path 1 corresponds to the atmospheric reflectance, path 2 is path that interacts with the target, it is the one which is useful to determine the surface reflectance, paths 3 and 4 contribute to the total reflectance but interact with the surface away from the target. These paths are thus the cause of adjacency effects and they blur the images.

 

 

If not corrected, adjacency effects may cause large errors. Let's take the case of a fully developed irrigated field surrounded by bare soil. For such a case, the second figure on the right shows the relative percentage of errors for reflectances and NDVI as a function of aerosol optical thickness, if adjacency effect is not corrected.

 

 

 

An approximate correction can be applied, but it thus requires to know the aerosol optical thickness. In our MACCS processor, here is how it works :

 

  1. We first correct the images under the assumption that the Landscape is uniform. We obtain a surface reflectance under uniform absorption which is noted  \rho_{s,unif} .
  2. We compute the neighbourhood reflectance (  \rho_{s,adj} ) using a convolution filter with a 2km radius, that computes the average neighborhood reflectance weighted by the distance to the target. To be fully rigorous, this filter should depend on the optical thickness and on the viewing and sun angle (The less aerosols, the larger radius), but as we did not work on an accurate model, we used a constant radius.
  3. We correct for the contribution of paths 3 and 4 using :

 \rho_{s}=\frac{\rho_{s,unif}.T^{\uparrow}.\frac{1-\rho_{s,unif}.s}{1-\rho_{s,adj}.s}-\rho_{s,adj}. T_{dif}^{\uparrow}}{T_{dir}^{\uparrow}}

  • where  T^{\uparrow}=T_{dif}^{\uparrow}+T_{dir}^{\uparrow} is the total upward transmission, sum of diffuse and direct upward transmissions, and s is the atmosphere spheric albedo. These quantities depend on the wavelength, on the aerosol model and on the AOT. They are computed using Look up Tables based on radiative transfer calculations.

 

As this processing uses convolution with a large radius, it takes quite a large part of the atmospheric processing time.

 

Result Exemples

The images below show 3 stages of the atmospheric processing, for 2 Formosat-2 images obtained over Montreal (Canada) with a 2 days interval. The first image was acquired on a hazy day (aerosol optical thickness (AOT) of 0.47 according to MACCS estimate); and the second one on a clear day (AOT=0.1).

  • The first line corresponds to the Top Of Atmosphere images, without atmospheric correction. The left image is obviously blurred compared to the right image.
  • The second line corresponds to the atmospheric correction under uniform landscape assumption (as in step 1). The left image is still obviously blurred compared to the right image.
  • the third line show the same images after adjacency effect correction. In that case, the left image is not blurred any more, it is even maybe a little over corrected as it seems somewhat sharper that the right image.

TOA Images (On the left, the hazy image)


Surface reflectance under uniform landscape assumption (on the left, the hazy image)

 

Surface reflectance after adjacency effect correction (on the left, the hazy image)

 

The pixel wise comparison of reflectances is also a way to show the enhancement due to the adjacency effect correction. The plot below compares the images of both dates corrected under the uniform landscape assumption (on the left), and after adjacency effect correction (on the right). You may observe that the dots are closer the the black diagonal on the right. On the hazy image (May 27th), the high reflectances are a little too low, while the low reflectances are a little too high, which is the symptom of a loss of contrast.

Release of the first version of SPOT4 (Take5) data

 

Phew ! These last days were quite intense, but we did it (almost) ! We are very proud to announce that the first version of SPOT4 (Take5) will be released on Monday the 15th of July at the following site : https://www.ptsc.fr/en/products/spot4-take5. The data format is explained here.


It will be possible to download the whole time series by clicking on 1C and 2A buttons. However downloading a whole series might take some time. Several small issues have been observed but not yet corrected (see this paragraph), but we still think it is useful that you start working on our data despite the few bugs that affect them. Downloading a data set requires the user to accept a licence. With the current version of the web site, you will need to accept this license each time you download a file.

In our initial schedule, we had announced a distribution of data during the month of June, but we are just a bit late, and happy to have done it, accounting for the numerous tasks we had to fulfil.

 

  • We had to obtain the experiment decision, which was done by the 11th of December 2012, and we had to implement it. The last data were acquired less than one month ago, and they will be available with the whole time series.
  • We had to modify and tune the processing chains, and the geometry  kept us busy for several weeks (and still does over equatorial forests). The level 2A processor worked quite well from the beginning, thanks to Mireille Huc (CESBIO) !
  • It was the first production of the French Land Data Center. The ground segment was built and implemented simultaneously to the data acquisition, its development team did a very nice work (namely, Dominique Clesse from CAP GEMINI, Hassan Makhmara and Joelle Donadieu at CNES). The data processing in this changing environment was perfectly managed by the exploitation team, Nicolas Prugent, Karl Rodriguez (Steria), Eric Faucher (CNES)...
  • A simple but very nice data server was set up in a very short time by Jerome Gasperi and Bernard Specht (CNES). The result is very simple and convenient (it helped us finding the last bugs...)

Issues and bugs

The products on the data server should be considered as a preliminary version of the processing. Other versions will be distributed, because the data sets still have a few defects :

 

  • Because of a little bug, the overlapping sites (Midi-Pyrénées East and West, and Britanny) were badly processed. These sites will be processed again and released in about 10 days. For the same reason, Maricopa site (USA) could not be processed, it will be released probably in September.
  • One image at the south of Languedoc site is always missing. We do not know why yet, but there is no doubt it will be corrected shortly.
  • Some Level 1A  products were not provided by Spot Image (25 out of 1600). We will need to reprocess the affected sites.
  • For all sites, the reference image for ortho-rectification were taken from LANDSAT  5 or 7. From our point of view, their location accuracy is a little insufficient. On the next version, we will replace LANDSAT data in FRANCE by the GeoSud high resolution cover from 2009. Out of France, we will replace them by LANDSAT 8 data whenever possible.
  • Data ortho-rectification often fails on some uniform equatorial forest sites (in particular for Congo(2), Gabon, Borneo, Sumatra). As  a result several dates are missing and those available are not perfectly registered.
  • On the level 2A, atmospheric corrections were performed with a constant aerosol model globally. We will use different aerosol models depending on the location in a future version.

 

Finally, we want to draw your attention to the fact that SPOT4 data are coded on 8 bits and can be saturated. SPOT system uses a data base of histograms to determine the gains to use depending on the location and date. This database is far from perfect and saturated data are quite common. For sites resulting from the merger of several images, it is also possible that the saturation thresholds differ between the left and the right half of the site. We provide you, within Level 1C and level 2A products with a mask saturations that you really should use.

 

Despite these few defects, we wish that this data set will be useful for you, and that you will obtain good results for your experiment and for your preparation for the arrival of Sentinel-2.

 

Many of the SPOT4(Take5) team will leave for vacations in the coming weeks, this blog will b much less active as I leave today !

 

 

 

 

 

L1C registration performances for SPOT4(Take5) V1 products

Now that all SPOT4(Take5) images have been processed (pfew !), we can make an appraisal of the performances. Let's start by the geometry, which caused us a lot of trouble :

  • SPOT4 has a location accuracy around 400 mètres, but during the experiment, it went through a fifteen day period when the location errors could reach 1500 m.
  • We seek a multi-temporal registration performance of 0.3 pixel RMS. This performance is difficult to measure because the measurement technique itself (correlation image matching) is not perfectly accurate.
  • We provide as a criterion the maximum registration error observed for the 50% best results or for the 80% best results. It is likely that the last criterion includes less inaccurate measurements.

 

Here are the observed performances for 3 very different sites :

  • CMaroc site, which is an arid site with a green period in march, a lot of blue sky, and high mountains (the Atlas). performances are excellent, with errors lower than 0.3 pixels for 50% of the measurementsl.

 

  • CBretagneLoireE site, which is a rather flat coastal area with large tides, and is often very cloudy. In that case, performances are still better than 0.5 pixels. The worse dates correspond to images with a large cloud cover, for which it is not easy to automatically collect accurate ground control points.

 

  • JSumatra site is a very flat area, covered with very uniform equatorial forest, and a large river whose limits change with time. In that case, the performance is really bad, with registration errors up to 10 pixels. This uniform site does not enable to find good control points, and the ones that are found are often along the river whose contour changes with the water level.

 

Conclusions

We have obtained very good results for most sites, with registration errors below 0.5 pixels (10m) even when the initial location error reaches 1500m. However, 4 sites are resisting to this processing. These 4 sites correspond to flat forest sites covered by equatorial forest : JSumatra, JBorneo, EGabon, ECongo. The ECongo site is even so uniform that it is not possible to measure its registration performance.

These sites will be distributed with the others in a few days with the first version of the products, but you should use them cautiously.

Finally, if the registration of 95% of images is good, the location performance is inherited from our reference images, ie LANDSAT (5 et 7). The next versions will be based on Geosud (IGN) images in France and on LANDSAT 8 data elsewhere. Performances should be enhanced in the next versions.

All the quicklooks

Before the distribution starts at the land data center (before mid July, or even sooner), I have updated all the quicklooks of all the images taken during the experiment. I have checked in the catalog to see if a few images had not been forgotten. I found about 20 images (on a total of 1600). These images will be processed soon.

 

You may find all the SPOT4 (Take5)  quicklooks following the links below, or via the SPOT4 (Take5) menu.