WO2016169699A1

WO2016169699A1 - A system, an apparatus and a method for determining mass change in a study area using remote sensing data

Info

Publication number: WO2016169699A1
Application number: PCT/EP2016/055289
Authority: WO
Inventors: Svein SOLBERG
Original assignee: Forest Vision As
Priority date: 2015-04-23
Filing date: 2016-03-11
Publication date: 2016-10-27
Also published as: NO20150499A1; NO340764B1

Abstract

A system, an apparatus and a method for determining mass change in a study area using remote sensing data. First and second data sets of the area are retrieved, which datasets are obtained at respective first and second point of time and comprises respective first and second set of pixel data given by respective first and second wavelength. Each pixel comprised in the respective first and second set of pixel data corresponds to a subarea of the area and a value for each pixel comprised in the respective first and second set is indicative of an amount of mass in the subarea at the respective first and second points of time. By first and second correction functions correcting for differences in penetration between the wavelengths and for technical errors, respectively, determining mass change in a subarea as a corrected difference between values of corresponding pixels of the second and first sets.

Description

A SYSTEM, AN APPARATUS AND A METHOD FOR DETERMINING MASS CHANGE IN A STUDY AREA USING REMOTE SENSING DATA

TECHNICAL FIELD

Embodiments herein relate generally to a system, an apparatus and a method therein for determining mass change in a study area of the earth. In particular they relate to the determination of mass change using remote sensing.

BACKGROUND

In remote sensing, acquisition of information about an object or phenomenon is performed without making physical contact with the object in questions. In modern usage, the term generally refers to the use of aerial sensor technologies to detect and classify objects on Earth, both on the surface, and in the atmosphere and oceans, by means of propagated signals, e.g. electromagnetic radiation. Remote sensing may be split into active remote sensing, wherein a signal is first emitted from for example an aircraft or satellites and then the reflected signal is recorded, and passive remote sensing, wherein no signal is emitted but only reflected radiation from e.g. sunlight is recorded.

Biomass is biological material derived from living, or recently living organisms, or it is used as a term for the mass of living biological organisms in a given area or ecosystem at a given time. With regard to remote sensing, the term "mass" is, as herein, usually not used directly in its physical sense, but as derived from sensing the volume or volume change, and knowing the mass density. However, remote sensing is not limited to electromagnetic waves: For example is the Gravity Recovery And Climate Experiment (GRACE), a joint mission of NASA and the German Aerospace Center, using two satellites and measuring their speed and the distance between them. Changes in Earth's gravity will change the speed and distance, and thus mass change can be detected directly. This can be used for sensing mass changes in the Earth's crust and in the oceans.

Optical satellite data systems and methods are the dominating remote sensing systems and methods of today. However, they have the limitation that the correlation with biomass is weak and tend to saturate at low levels. By the expression "saturate at low levels" is meant that the optical satellite systems and methods are only sensitive to biomass variation at low biomass levels, while the signal is more or less constant from intermediate to high biomass levels. Further, it is weather sensitive and clouds may prevent and/or deteriorate data acquisition in some areas. Another limitation with the optical method is that it merely tracks land cover changes, whereby it is only possible to detect changes from forested land to non-forested land and vice versa. Further, by the optical systems and methods, forest degradation is hardly detectable. This is a drawback with the optical systems and methods, since forest degradation make up a considerable share of cuttings in some areas of the earth, e.g. in the tropics. Finally, the conversion of the annual clear-cut area into changes in forest C stocks is crudely obtained by using fixed emission factors.

Thus, and as described in "Forest biomass change estimated from height change in interferometric SAR height models" by Solberg et al. Carbon Balance and Management 2014, 9:5, there is a need for new satellite remote sensing methods for monitoring tropical forest carbon stocks. Advanced RADAR instruments on board satellites can contribute with novel methods. RADARs can see through clouds, and furthermore, by applying stereo

RADAR imaging forest height and its changes can be measured. Such height changes are related to carbon stock changes in the biomass. Data from the current Tandem-X satellite mission, where two RADAR equipped satellites go in close formation providing stereo imaging, may be applied and combined with similar data acquired with one of the space shuttles in the year 2000, i.e. the so-called Shuttle Radar Topography Mission (SRTM). Height information from a RADAR image pair may be derived using a method called interferometry.

TerraSAR-X add-on for Digital Elevation Measurement (TanDEM-X) is the name of the TerraSAR-X's twin satellite, a German Earth observation satellite using Synthetic Aperture Radar (SAR) that is a modern radar imaging technology. TanDEM-X is also the name of the satellite mission flying the two satellites in a closely controlled formation with typical distances between 250 and 500 m his unique twin satellite constellation generated the WorldDEM, which refers to the seamless global digital elevation models available from 2014.

LANDSAT is a program for acquiring images using multispectral scanning from a number of satellites starting in 1973. Image data from LANDSAT has been used to create a global dataset characterizing global forest extent and change from 2000 through 2013. The data set is in this disclosure referred to as "Hansen data", the full reference is: Hansen, M. C, P. V. Potapov, R. Moore, M. Hancher, S. A. Turubanova, A. Tyukavina, D. Thau, S. V. Stehman, S. J. Goetz, T. R. Loveland, A. Kommareddy, A. Egorov, L. Chini, C. O. Justice, and J. R. G. Townshend. 2013. "High-Resolution Global Maps of 21 st-century Forest Cover Change." Science 342 (15 November): 850-53. Data available on-line from:

http://earthenginepartners.appspot.com/science-2013-global-forest. The version that is used in this disclosure is called version 1 .1 . The forest coverage is assessed and encoded as a percentage for tree canopy cover and as bit 1/0 for loss/no loss or gain/no gain, using computer algorithms. Thus a forest cover data set exists that is based on two-dimensional (2D) optical data. Other data set could be used in the same way, e.g. dataset compiled from aerial or satellite remote sensing, such as laser scanning, from which canopy cover may be derived.

Remote sensing methods that provide three-dimensional (3D) data have a considerable advantage in comparison with 2D optical data. The 3D data provide measurements of forest height, and its changes, which is crucial for forest biomass changes. Airborne 3D remote sensing, i.e. Airborne Laser Scanning (ALS) or stereo photogrammetry, is a more accurate tool than optical satellite data and it can detect also forest degradation and growth. However, for most countries, applications with complete ALS coverage are cost- prohibitive. The feasible application of ALS would be strip sampling, which may provide accurate estimates of C stocks changes compared to other methods.

Therefore, there is i.a. a need for improvements in the remote sensing systems, apparatus and methods for monitoring mass changes.

SUMMARY

An object of embodiments herein is to provide a way of improving the determination of mass change in a study area of the earth.

According to a first aspect of embodiments herein, the object is achieved by a method in an apparatus for determining mass change in a study area of the earth using remote sensing data.

The apparatus retrieves a first three-dimensional, 3D, dataset of the study area, which first 3D dataset is obtained by means of remote sensing at a first point of time and comprises a first set of pixel data given by a first wavelength of the remote sensing. Each pixel comprised in the first set of pixel data corresponds to a subarea of the study area and a data value for each pixel comprised in the first set of pixel data is indicative of an amount of mass in the subarea at the first point of time.

Further, the apparatus retrieves a second 3D dataset of the study area, which second 3D dataset is obtained by means of remote sensing at a second point of time and comprises a second set of pixel data given by a second wavelength of the remote sensing. Each pixel comprised in the second set of pixel data corresponds to a subarea of the study area and a data value for each pixel comprised in the second set of pixel data is indicative of an amount of mass in the subarea at the second point of time.

By means of a first correction function, Diffp_ene. configured to correct for the differences in penetration between the first and second wavelengths, and by means of a second correction function, Difft_ech. configured to correct for technical errors, the apparatus determines mass change in one or more subareas of the study area as a corrected difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas.

According to a second aspect of embodiments herein, the object is achieved by an apparatus for determining mass change in a study area of the earth using remote sensing data.

The apparatus is configured to retrieve a first three-dimensional, 3D, dataset of the study area, which first 3D dataset is obtained by means of remote sensing at a first point of time and comprises a first set of pixel data given by a first wavelength of the remote sensing. Each pixel comprised in the first set of pixel data corresponds to a subarea of the study area and a data value for each pixel comprised in the first set of pixel data is indicative of an amount of mass in the subarea at the first point of time.

Further, the apparatus retrieves a second 3D dataset of the study area, which second

3D dataset is obtained by means of remote sensing at a second point of time and comprises a second set of pixel data given by a second wavelength of the remote sensing. Each pixel comprised in the second set of pixel data corresponds to a subarea of the study area and a data value for each pixel comprised in the second set of pixel data is indicative of an amount of mass in the subarea at the second point of time.

By means of a first correction function, Diffpene. configured to correct for the differences in penetration between the first and second wavelengths, and by means of a second correction function, Difftech. configured to correct for technical errors, the apparatus is configured to determine mass change in one or more subareas of the study area as a corrected difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas.

According to a third aspect of embodiments herein, the object is achieved by a remote sensing system for determining mass change in a study area of the earth using remote sensing.

The remote sensing system comprises at least one transport means comprising at least one sensor and configured to transport the at least one sensor at a distance over the study area.

The at least one sensor is configured to acquire remote sensing data by, at a certain point of time, emit an electromagnetic radiation of a predefined wavelength towards the study area of the earth and receive backscattered electromagnetic radiation from the study area. Thereby, a 3D dataset is obtained for the certain point of time, which 3D data set comprises a set of pixel data for the predefined wavelength of the emitted electromagnetic radiation, wherein each pixel comprised in the set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the set of pixel data is indicative of an amount of mass in the subarea at the certain point of time.

Further, the remote sensing system comprises a data storage means configured to store the 3D dataset acquired by the at least one sensor; and the apparatus for determining mass change in a study area of the earth using remote sensing data.

According to a fourth aspect of embodiments herein, the object is achieved by a computer program, comprising instructions which, when executed on at least one processor, causes the at least one processor to carry out the method in the arrangement for determining mass change in a study area of the earth using remote sensing data.

Since the mass change in one or more subareas of the study area is determined as a corrected difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas, datasets acquired by different remote sensing technologies may be used for the mass change determination. For example, by means of embodiments herein, the first set of pixel data may comprises from a Shuttle Radar Topography Mission (SRTM) Digital Surface Model (DSM), sometimes also referred to as Digital Elevation Model (DEM), acquired at a first point of time and the second set of pixel data may comprise a Tandem-X DSM (e.g. the WorldDEM™) data acquired at a second point of time. This results in an improved determination of mass change in a study area of the earth.

BRIEF DESCRIPTION OF DRAWINGS

Examples of embodiments herein are described in more detail with reference to attached drawings in which:

Figure 1 schematically illustrates embodiments of a remote sensing system;

Figure 2 schematically illustrates the difference in penetration between electromagnetic radiations of different wavelengths;

Figure 3 is a flowchart depicting embodiments of a method in an apparatus; and

Figure 4 is a schematic block diagram illustrating embodiments of an apparatus. DETAILED DESCRIPTION

As part of developing embodiments herein, some problems will first be identified and discussed. Management of forest Carbon (C) stocks is increasingly addressed due to its impact on the global greenhouse gas cycle and climate. Deforestation contributes to a significant fraction of the total anthropogenic C emissions. The suite of methods for mapping, monitoring and estimating parts of the forest C cycle is expanding and they relate to field inventory, modelling and remote sensing.

Deforestation in the tropics is of particular significance due to its rapid speed, and the

Reducing Emissions from Deforestation and Forest Degradation in Developing Countries (REDD) initiative aims at reducing the C losses through performance-based economic compensation by comparison of performance against a reference level, such as a business- as-usual reference emission level. In addition to deforestation, forest degradation and enhancement of carbon stock through forest growth are other components of the forest C changes that should influence the REDD credits. In order to realize this payment-for- ecosystem-service, the tropical countries need to document their annual changes in forest C stocks, and satellite remote sensing is likely to be a major data provider for this. Such data would also enable detection of logging areas for possible counteractions. However, there is a need for new satellite remote sensing methods for REDD. Firstly, there is a need for methods that overcome the limitations of today's methods, which methods are weather sensitive, e.g. sensitive to clouds, have a small areal coverage and fail to detect C stock changes other than deforestation. Secondly, there is a need for historical data on forest changes for the business-as-usual emission level.

Methods based on active remote sensing such as satellite Synthetic Aperture RADAR

(SAR) methods may resolve the limitations with optical data methods. As previously mentioned, the SAR method and system is an imaging RADAR method and system. The cloud problem mentioned above in conjunction with the optical satellite system described is non-existent with SAR, because the longer wavelengths, commonly 3-70 cm, penetrate clouds.

Logging leads to a reduction in height, while forest growth leads to an increase in height. Heights may be derived from SAR data in two or more ways, mainly by phase differences, e.g. by means of Interferometric SAR (InSAR), or by parallaxes e.g. by means of radargrammetry, in a SAR image pair.

Further, from changes in surface height obtained from InSAR, changes in biomass and carbon stocks may be retrieved. As mentioned above, with the satellite based InSAR method the heights are derived from phase differences between two SAR images taken from different positions in space. This may provide accurate height measurements. However, the accuracy depends on various acquisition properties. In particular for short wavelength SAR, e.g. microwaves in the X-band, over forested areas, the SAR imaging needs to be carried out by two satellites going together in a close formation. This is called a bi-static, or a single pass, acquisition, where one satellite is submitting microwave RADAR pulses and both satellites receive the same echoes from Earth's surface. A bi-static acquisition removes temporal de-correlation, i.e. phase noise caused by differences in the position of branches and in moisture. Phase noise is also influenced by the distance between the satellites, i.e. the baseline, which should neither be too large nor too small. With increasing baseline there is an increasing volume de- correlation caused by an increasing difference in the look angle, herein sometimes also referred to as the local incidence angle, into the canopy volume. Contrary to this, when the baseline is very small the noise increases because of quantization errors, i.e. a given height correspond to a tiny fraction of a 2π cycle of phase difference. In order to compress data onboard the satellite the data is typically compressed to 5 or 8 bit, and a tiny fraction of such a number will correspond to a crude height measurement. Finally, random errors and phase noise will be relatively large when the backscatter signal is weak, i.e. a low signal-to-noise ratio, depending on the incidence angle, the polarization and the topography.

It has been demonstrated that forest biomass, or the equivalent stem volume, is strongly related to the so-called InSAR height. The InSAR height is the height above ground of the center of the SAR echo. This relationship may vary with stand structure and in particular with tree number density. In addition, loggings, i.e. clear-cutting and forest degradation, may be detected as temporal changes in surface height.

In order to compare a Shuttle Radar Topography Mission (SRTM) Digital Surface Model (DSM), sometimes also referred to as Digital Elevation Model (DEM), and a Tandem- X DSM (e.g. the WorldDEM™), and for deriving changes in forest height and furthermore changes in forest carbon stocks, there is a need to perform one or more corrections, e.g. two corrections, of the acquired data which will be described below. The first correction function Diff_pen_e provides height correction for differences in penetration into a vegetation volume or any other land cover with semi-permeable property or porous property. The second correction function Diff_teCh provides height correction for technical errors typically stemming from phase-unwrapping errors.

Thus, firstly, a correction based on the different wavelengths of the electromagnetic radiation, such as microwaves used may be performed. The near-global SRTM DSM used radar microwaves in the C-band having a wavelength of 5.6 cm, while Tandem-X uses radar microwaves in the X-band having a wavelength of 3.1 cm. The penetration of the radar microwaves into vegetation increases with wavelength, and hence, a C-band DSM is deeper down in the forest canopy than the X-band DSM. However, in areas of the earth without forest, i.e. in no-forested areas, the two models are identical for most cases. However, one exception is dry sand. The height difference between a C-band DSM and an X-band DSM in a forest depends on the forest cover. In some embodiments, and in order to be able to compare the C-band DSM and the X-band DSM, a correction value that increases monotonically with increasing forest cover is used. Further, in some embodiments, the correction value increases with approximately 1-3 cm per %-unit increase in forest cover, e.g. approximately 2 cm per %-unit increase in forest cover.

Secondly, large-scale errors, such as errors stemming from problems in steep terrain and phase unwrapping, may exist in both the SRTM DSM and the Tandem-X DSM. Such errors may be considerably larger in the SRTM DSM than in the Tandem-X DSM since the Tandem-X data have a higher spatial resolution, more accurate orbit data and a more accurate processing scheme including radargrammetry as a first step and Interferometric processing (InSAR) as a second step. Thus, it is meaningful to do the corrections on the SRTM data. The SRTM DSM data were acquired during February 12-20, 2000, and the Tandem-X DSM data were acquired in both 2000 and 2012. For the correction, pixels having no forest cover in 2000 and no change in forest cover during the years 2000 - 2012, as derived from the Hansen data were identified. These pixels should have no height change.

The Hansen data is divided into 10x10 degree tiles, consisting of seven files per tile. All files contain unsigned 8-bit values and have a spatial resolution of 1 arc-second per pixel, or approximately 30 meters per pixel at the equator. The seven files includes files for tree cover in the year 2000, defined as canopy closure for all vegetation taller than 5m in height and encoded as a percentage per output grid cell, and files for loss, gain and loss per year.

Reference composite imagery are median observations from a set of quality assessed growing season observations in four spectral bands. Normalized top-of-atmosphere (TOA) reflectance values (p) have been scaled to an 8-bit data range. Below, embodiments herein will be illustrated in more detail by a number of exemplary embodiments. It should be noted that these embodiments are not mutually exclusive.

Components from one embodiment may be tacitly assumed to be present in another embodiment and it will be obvious to a person skilled in the art how those components may be used in the other exemplary embodiments.

Further, the term "wavelength" is sometimes herein used as a synonym to "a physical property of a remote measurement". Most remote sensing is based on direct sensing and recording of electromagnetic or other waves, such as light or radar. For remote sensing using gravity, such as GRACE, the sensing is not directly using a wavelength - and a first and second wavelength as used in the present disclosure is then equal to using two different set of satellites with different physical properties.

Furthermore, the term "pixel" is sometimes herein used for a single scalar element of a multi-component representation of the sensed image, and is then the smallest single component of that digital image. It is usually represented as a string of bits stored in a memory position. The bits may represent several sensed physical properties such as for light, the colour or the wavelength and the intensity together with the address of the pixel within the image, corresponding to the sensed location. For a radar-based image, there is of course no colour, rather the pixels are a complex phasor (also called complex amplitude) representation of the coherent backscatter from the resolution element on the ground and the propagation phase delay. When SAR with two sensors are used, the coherent backscatter will be the same, but the propagation phase delay will be different. The representation of the coherent backscatter is then an intensity that can be used for deducing the mass of the resolution element, i.e the object that causes the backscatter.

Yet further, even if some embodiments herein are described with reference to forest carbon stocks it should be understood that some embodiments may relate to other applications wherein a mass change in a study area of the earth using remote sensing is to be determined. The mass to study may be biomass, ice mass, water mass, land mass, mass or a part of the earth crust, and/or tectonic mass just to give some examples. Further, depending on the mass to study different wavelengths of the electromagnetic radiation of the remote sensing may be suitable or other differences as for the GRACE data measuring changes in gravity. Thus, the wavelengths mentioned herein and the data sets acquired by those wavelengths are only given as examples.

Figure 1 schematically illustrates embodiments of a remote sensing system 100 for determining mass change in a study area of the earth using remote sensing. As previously mentioned, the mass to study may be biomass, ice mass, water mass, land mass, mass or a part of the earth crust, and/or tectonic mass just to give some examples.

The remote sensing system 100 comprises at least one transport means

102, 102a, 102b such as an airplane, a satellite or a space shuttle. The at least one transport means 102, 102a, 102b is configured to move, e.g. fly, at a distance from the earth surface, e.g. at a distance from the study area of the earth. In some embodiments, the transport means 102, 102a, 102b is configured to orbit around the earth or part thereof while collecting data by remote sensing. At least one sensor 104,106 is comprised in or attached to the at least one transport means 102, 102a, 102b. The at least one sensor 104,106 is or comprises an active remote sensing sensor such as radar, e.g. radar imaging or radargrammetry, and lidar, and/or a passive remote sensing sensor base on e.g. photogrammetry and stereophotogrammetry.

In some embodiments, the at least one sensor 104,106 is configured to acquire remote sensing data by, at a certain point of time, emit an electromagnetic radiation E of a predefined wavelength towards the study area of the earth and receive backscattered electromagnetic radiation B from the study area. Thereby, a 3D dataset is obtained for the certain point of time, which 3D data set comprises a set of pixel data for the predefined wavelength of the emitted electromagnetic radiation E. Each pixel comprised in the set of pixel data corresponds to a subarea of the study area and a data value for each pixel comprised in the set of pixel data is indicative of an amount of mass in the subarea at the certain point of time.

The remote sensing system 100 comprises further a data storage means 108 configured to store remote sensing data, e.g. one or more 3D datasets, acquired by the at least one sensor 104,106. The data storage means 108 may be comprised in the transport means 102, 102a, 102b or it may be arranged externally of the transport means 102,

102a, 102b. In some embodiments, wherein the data storage means 108 is an external storage, the transport means 102, 102a, 102b and/or the at least one sensor 104,106 may use a wireless communication path for communicating with, e.g. storing remote sensing data in, the data storage means 108.

The data storage means 108 may comprise one or more memory units. Further, the data storage means 108 may be a computer data storage or a semiconductor memory such as a computer memory, a read-only memory, a volatile memory or a non-volatile memory.

The remote sensing system 100 comprises further an apparatus 200 for determining mass change in a study area of the earth using remote sensing data. The apparatus will be described in more detail with reference to Figures 3 and 4.

In some embodiments, the remote sensing system 100 comprises a first transport means 102a and a second transport means 102b, such as two satellites 102a, 102b. The first and second transport means 102a, 102b comprise a first sensor 104 and a second sensor 106, respectively. One example of such a remote sensing system 100 is the TanDEM-X, wherein the two sensors 104,106 comprised in or attached to a respective satellite

102a, 102b, are flying in a closely controlled formation with typical distances between 250 and 500 m. This twin satellite constellation will allow the generation of the remote sensing data referred to as WorldDEM™. Herein this kind of remote sensing data is sometimes referred to as X-band TanDEM-X data. In embodiments having a twin satellite constellation, one of the sensors 104,106 transmits the electromagnetic radiation having a desired wavelength in the X-band and both of the sensors simultaneously register the backscattered electromagnetic radiation. Thereby, the sensors 104,106 may acquire an interferometric data pair without any time difference.

Figure 2 schematically illustrates the difference in penetration between

electromagnetic radiations of different wavelengths. As schematically illustrated in Figure 2, the penetration of electromagnetic radiation into the vegetation increases with increasing wavelength of the radiation. Further, in Figure 2, the reference numeral E indicates emitted electromagnetic radiation, i.e. electromagnetic radiation emitted from a radiation source, and the reference numeral B indicates backscattered electromagnetic radiation, i.e. reflected electromagnetic radiation. The wavelengths of the illustrated electromagnetic radiation will now be exemplified. The sun in Figure 2 illustrates visible light that has a wavelength in the range of 400 nanometres to 700 nanometres. Electromagnetic radiation in the so-called X- Band has a wavelength in the range of 2.5 - 3.75 centimetres. The wavelength of electromagnetic radiation in the C-Band is in the range of 3.75 - 7.5 centimetres.

Electromagnetic radiation in the L-band has a wavelength in the range of 15-30 centimetres, and electromagnetic radiation in the P-band has a wavelength in the range of 30-100 centimetres.

A method in an apparatus 200 for determining mass change in a study area of the earth using remote sensing will now be described with reference to a flow chart depicted in Figure 3. The method comprises one or more of the following actions. It should be understood that actions may be taken in another suitable order and that actions may be combined.

Action 301

The apparatus 200 retrieves a first three-dimensional (3D) dataset DS1 of the study area. This may also be expressed as the first 3D dataset DS1 is associated with the study area. The apparatus 200 may be configured to retrieve the first 3D dataset DS1 from a data storage means, such as the data storage means 108 of the remote sensing system 100. The first 3D dataset DS1 is obtained by means of remote sensing, such as by means of sensors 102,104 of the remote sensing system 100, at a first point of time T-| .

The first 3D dataset DS1 comprises a first set of pixel data given by a first wavelength of the remote sensing. Further, each pixel comprised in the first set of pixel data corresponds to a subarea of the study area and a data value for each pixel comprised in the first set of pixel data is indicative of an amount of mass in the subarea at the first point of time TV

For example, the first set of pixel data may comprise C-band SRTM data acquired at the first point of time T In some embodiments, the first set of pixel data comprises C-band SRTM data acquired in the year 2000.

Action 302

The apparatus 200 retrieves a second 3D dataset DS2 of the study area. This may also be expressed as the second 3D dataset DS2 is associated with the study area. The apparatus 200 may be configured to retrieve the second 3D dataset DS2 from a data storage means, such as the data storage means 108 of the remote sensing system 100. The second 3D dataset DS2 is obtained by means of remote sensing, such as by means of the remote sensing system 100, at a second point of time T2. T2 is any point of time subsequent the first point of time Ti , and Ti and T2 may be any future points of time.

The second 3D dataset DS2 comprises a second set of pixel data given by a second wavelength of the remote sensing. Further, each pixel comprised in the second set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the second set of pixel data is indicative of an amount of mass in the subarea at the second point of time T2.

For example, the second set of pixel data may comprise X-band TanDEM-X data acquired at the second point of time T2. In some embodiments, the second set of pixel data comprises X-band TanDEM-X data acquired in the year 2012.

Action 303

In some embodiments, the apparatus 200 determines a third 3D dataset Diffl as a difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas of the study area. The determined third 3D dataset Diffl corresponds to an uncorrected estimate of the mass change in the one or more subareas between the first point of time Ύ and the second point of time T₂.

Action 304

In some embodiments, the apparatus 200 determines a first correction function Diffpene- The first correction function Diffp_ene corrects for differences in penetration into a volume of mass, e.g. a vegetation volume, due to differences in wavelengths of the electromagnetic radiation used for collecting the remote sensing data. The first correction function Diffpene comprises a third set of pixel data, wherein each pixel comprised in the third set of pixel data corresponds to a subarea of the study area, and wherein a data value for each pixel comprised in the third set of pixel data is indicative of a difference in penetration between the first wavelength and the second wavelength.

In some embodiments relating to forest carbon stocks, the mass is biomass, the remote sensing is radar imaging, and the first and second 3D datasets are respective first and second radar images. In such embodiments, the apparatus 200 determines the first correction function Diffp_ene by selecting an image area of the first radar image taken at the first point of time ΤΊ , determining a sixth 3D dataset Diffx-c and by determining the first correction function Diff_pene by correcting the sixth 3D dataset Diffx.c for technical errors using the second correction function, Difftech-

The sixth 3D dataset Diffx.c may comprise the difference between data values of corresponding one or more pixels of a seventh set of pixel data and the first set of pixel data, wherein the seventh set of pixel data is given by the first wavelength of the radar imaging at the first point of time T-| , and wherein corresponding one or more pixels correspond to the image area.

For example, the seventh set of pixel data may comprise X-band SRTM data acquired at the first point of time T-| . In some embodiments, the seventh set of pixel data comprises X-band SRTM data acquired in the year 2000.

As previously described in Action 301 above, the first set of pixel data may comprise C-band SRTM data acquired at the first point of time T-| . In some embodiments, the first set of pixel data comprises C-band SRTM data acquired in the year 2000.

Thus, in some embodiments, the sixth 3D dataset Diffx.c comprises the difference between X-band SRTM data acquired at the first point of time Ti and C-band SRTM data acquired at the first point of time T-| .

In some embodiments relating to forest carbon stocks, the apparatus 200 determines the first correction function Diffp_ene as approximately 1 -3 cm per percentage increase in forest cover, especially as approximately 2 cm per percentage increase in forest cover. In such embodiments, the forest cover is given by a fifth set of pixel data comprising data acquired at the first point of time T-| by means of passive remote sensing and sometimes herein referred to as the Hansen data. The fifth set of pixel data will be described in Action 305 below. In some embodiments relating to ice mass/water mass/land mass/mass/a part of the earth crust, and/or tectonic mass, the apparatus 200 determines the first correction function Diffpene based on relevant physical parameters for the reflection of the used wavelengths. E.g. for ice moisture and temperature are such data from a a radar, such as that in the SRTM data, , will have different reflections from new, dry snow, where the radar penetrates, and from wet snow that reflects the radar. Thus, meteorological data with temperature and moisture for the study area can be used as basis for the correction function.

Action 305

In some embodiments, the apparatus 200 determines a second correction function

Difftech- The second correction function Difftech corrects for technical errors such as errors stemming from for example phase-unwrapping errors. The second correction function Difftech comprises a fourth set of pixel data, wherein each pixel comprised in the fourth set of pixel data corresponds to a subarea of the study area, and wherein a data value for each pixel comprised in the fourth set of pixel data is indicative of a technical error.

In some embodiments relating to forest carbon stocks, the mass is biomass, the remote sensing is radar imaging, and the first and second 3D datasets are respective first and second radar images. In such embodiments, the apparatus 200 determines the second correction function Difftech by retrieving a fourth 3D dataset DS4 of the study area. The fourth 3D dataset DS4 may be obtained by means of passive remote sensing at the first point of time Ti and comprises a fifth set of pixel data. The passive remote sensing may be as in action 304. In some embodiments, the fifth set of pixel data comprises the so-called Hansen data acquired in the year 2000. Further, each pixel comprised in the fifth set of pixel data corresponds to a subarea of the study area. Furthermore, a data value for each pixel comprised in the fifth set of pixel data is indicative of forest cover in the subarea at the first point of time T-| .

Yet further, the apparatus 200 may retrieve a fifth 3D dataset DS5 of the study area, which fifth 3D dataset is obtained by means of passive remote sensing at the second point of time T₂ and comprises a sixth set of pixel data, wherein each pixel comprised in the sixth set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the sixth set of pixel data is indicative of forest cover in the subarea at the second point of time. In some embodiments, the sixth set of pixel data comprises the so- called Hansen data acquired in the year 2012.

Furthermore, the apparatus 200, by means of the fifth and sixth sets of pixel data, determines one or more pixels in the first radar image which one or more pixels comprise less than a predetermined forest cover and which one or more pixels have no change in forest cover between the first and second points of time. These pixels represent vegetation free pixels during the period between the first point of time Ti and the second point of time T2, e.g. between the years 2000 - 2012. Further, these pixels should ideally have identical heights, and any difference here is due to technical errors.

When one or more data values indicative of an artefact are comprised in the third 3D dataset Diffi described under Action 303 above, the apparatus 200 removes the one or more values indicative of the artefact from the third 3D dataset Diffi and determines the second correction function Difftech as an interpolation of the third 3D dataset Diffi . In some embodiments, data values smaller than -20m and larger than +20m are indicative of one or more artefacts and the apparatus 200 removes such data values from the third 3D dataset

Diffi

Thus, the height differences in the pixels were calculated, extreme values removed and then the height difference from the pixels were interpolated for the entire subarea. The interpolation may be performed by using the Minimum Curvature algorithm, splines or kriging. The interpolated height difference raster represents the second correction function Difftech to be applied on for example the SRTM DSM data in order to correct for the technical errors between the X-band Tandem-X and the C-band SRTM DSMs. In some embodiments relating to ice mass/water mass/land mass/mass/a part of the earth crust, and/or tectonic mass, the apparatus 200 determines the second correction function Difftech in the same way as for another embodiment, such as the one for forest coverage. The correction function can be similar when the data is from the same apparatus, such as from the same pair of satellites.

Action 306

By means of the first correction function Diff_pene and the second correction function Difftech. the apparatus 200 determines mass change in one or more subareas of the study area as a corrected difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas.

As previously mentioned, in some embodiments, the first set of pixel data may comprise C-band SRTM data acquired at the first point of time Ti e.g. acquired in the year 2000. In some embodiments, the second set of pixel data may comprise X-band TanDEM-X data acquired at the second point of time J 2, e.g. acquired in the year 2012.

Further and as previously mentioned, the first correction function Diff_pene is configured to correct for the differences in penetration between the first and second wavelengths, and the second correction function Difftec is configured to correct for technical errors.

In some embodiments, the apparatus 200 determines mass change in one or more subareas by further correcting the data values of the determined third 3D dataset Diffi by means of the first correction function Diff_pene and the second correction function Difftech. whereby the corrected data values of the determined third 3D dataset Diff-i are indicative of the mass change in one or more subareas of the study area.

As previously mentioned in Action 303 above, the third 3D dataset Diffi comprises the difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data.

Action 307

In some embodiments, the apparatus 200 displays the determined mass change in the one or more subareas of the study area. The apparatus 200 may for example display an image of the determined mass change in the one or more subareas of the study area.

Action 308

In some embodiments, the apparatus 200 stores the determined mass change. For example, the apparatus 200 may store the determined mass change in a Tagged Image File Format (TIFF), such as a GeoTIFF format. The apparatus 200 may also store additional information, in the same or in different files, such as where and how much corrections have been made from correction functions, such as Diffpene and Difftech .

To perform the method for determining mass change in a study area of the earth using remote sensing data, the apparatus 200 may comprise an arrangement depicted i Figure 4.

In some embodiments, the apparatus 200 comprises an input and/or output interface 400. The input and/or output interface 400 may be configured to receive input data, instructions, commands, etc. and/or to transmit output data or information such as result data or information relating to processed data. The apparatus 200 comprises further means, such as e.g. a retrieving module 401 , adapted to retrieve datasets, such as first and second 3D datasets DS1 ,DS2 of or associated with the study area. The retrieving module 401 may be implemented as a processor 407 of the apparatus 200. The processor 407 will be described in more detail below.

The retrieving module 401 may be configured to retrieve the datasets from a data storage means, such as the data storage means 108 of the remote sensing system 100, by means of wireless communication.

As previously mentioned, the first 3D dataset is obtained by means of remote sensing at a first point of time Ti and comprises a first set of pixel data given by a first wavelength of the remote sensing, wherein each pixel comprised in the first set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the first set of pixel data is indicative of an amount of mass in the subarea at the first point of time T-i. As previously mentioned, in some embodiments, the first set of pixel data comprises C-band SRTM data acquired in the year 2000. Further, as also previously mentioned, the second 3D dataset is obtained by means of remote sensing at a second point of time T2 and comprises a second set of pixel data given by a second wavelength of the remote sensing, wherein each pixel comprised in the second set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the second set of pixel data is indicative of an amount of mass in the subarea at the second point of time T2. As previously mentioned, in some embodiments, the second set of pixel data comprises X-band TanDEM-X data acquired in the year 2012.

In some exemplifying embodiments, wherein the mass is biomass, the remote sensing is radar imaging, the first and second 3D datasets are respective first and second radar images, the retrieving module 401 is configured to retrieve a fourth 3D dataset of the study area, which fourth 3D dataset is obtained by means of passive remote sensing at the first point of time Ti and comprises a fifth set of pixel data, wherein each pixel comprised in the fifth set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the fifth set of pixel data is indicative of forest cover in the subarea at the first point of time T-| . In some embodiments, the fifth set of pixel data comprises the so-called Hansen data acquired in the year 2000.

Further, in such embodiments the retrieving module 401 is configured to retrieve a fifth 3D dataset of the study area, which fifth 3D dataset is obtained by means of passive remote sensing at the second point of time T2 and comprises a sixth set of pixel data, wherein each pixel comprised in the sixth set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the sixth set of pixel data is indicative of forest cover in the subarea at the second point of time T2. In some

embodiments, the sixth set of pixel data comprises the so-called Hansen data acquired in the year 2012.

The apparatus 200 comprises further means, such as e.g. a determining module 402, adapted to determine mass change in one or more subareas of the study area. The determining module 402 may be implemented as the processor 407 of the apparatus 200.

The determining module 402 is configured to, by means of the first correction function Diffpene and the second correction function Difftech. determine mass change in one or more subareas of the study area as a corrected difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas.

As previously mentioned the first correction function Diff_pene is configured to correct for the differences in penetration between the first and second wavelengths, and the second correction function Difftech is configured to correct for technical errors.

In some embodiments, the determining module 402 is configured to determine the first correction function Diffpene. which first correction function Diff_pene comprises the third set of pixel data. Each pixel comprised in the third set of pixel data corresponds to a subarea of the study area, and a data value for each pixel comprised in the third set of pixel data is indicative of the difference in penetration between the first wavelength and the second wavelength.

Further, in some embodiments, the determining module 402 is configured determine the second correction function Difftech. which second correction function, Difftech. comprises the fourth set of pixel data. Each pixel comprised in the fourth set of pixel data corresponds to a subarea of the study area, and wherein a data value for each pixel comprised in the fourth set of pixel data is indicative of a technical error.

In some embodiments, the determining module 401 may further be configured to determine a third 3D dataset Diffi as a difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which

corresponding one or more pixels correspond to the same one or more subareas. In such embodiments and in order to determine mass change in one or more subareas of the study area as the corrected difference, the determining module 401 is further configured to correct the data values of the determined third 3D dataset Diffi by means of the first correction function Diffpene and the second correction function Difft_ech. whereby the corrected data values of the determined third 3D dataset Diffi are indicative of the mass change in the one or more subareas of the study area.

In some exemplifying embodiments, wherein the mass is biomass, the remote sensing is radar imaging, wherein the mass is biomass, the remote sensing is radar imaging, the first and second 3D datasets are respective first and second radar images, the determining module 402 is further configured to, by means of the fifth and sixth sets of pixel data, determine one or more pixels in the first radar image which one or more pixels comprise less than a predetermined forest cover and which one or more pixels have no change in forest cover between the first and second points of time. Further, when one or more data values indicative of an artefact are comprised in the third 3D dataset, Diffi , the determining module is configured to remove the one or more values indicative of the artefact from the third 3D dataset Diffi and to determine the second correction function Difftech as an interpolation of the third 3D dataset Diffi _.

In order to determine the first correction function Diffp_ene. the determining module 402 may further be configured to select an image area of the first radar image taken at the first point of time Ti . The image area may be part of the image of the study area but it could also be part of an image of another area having similar properties, e .g. forest type, as the study area. The determining module 402 may further be configured to determine a sixth 3D dataset, Diffx-c, comprising the difference between data values of corresponding one or more pixels of a seventh set of pixel data and the first set of pixel data, wherein the seventh set of pixel data is given by the first wavelength of the radar imaging at the first point of time

Ti , and wherein corresponding one or more pixels correspond to the image area.

As mentioned above, the seventh set of pixel data may comprise X-band SRTM data acquired at the first point of time Ti . In some embodiments, the seventh set of pixel data comprises X-band SRTM data acquired in the year 2000.

The determining module 402 may then determine the first correction function Diff_pene by correcting the sixth 3D dataset Diffx-c for technical errors using the second correction function Difftech- In some embodiments relating to forest carbon stocks and in order to determine the first correction function Diffp_ene. the determining module 402 may further be configured to determine the first correction function Diff_pene as approximately 1 -3 cm per percentage increase in forest cover, especially as approximately 2 cm per percentage increase in forest cover, wherein the forest cover is given by the fifth set of pixel data.

In some embodiments, the apparatus 200 comprises means, such as e.g. a displaying module 403, adapted to display a determined mass change, e.g. on a display 406 comprised in or associated with the apparatus 200. The displaying module 403 may be implemented as the processor 407 of the apparatus 200.

The displaying module 403 may be configured to display the determined mass change in the one or more subareas of the study area. In some embodiments, the displaying module 403 displays, e.g. on the display 406, an image of the determined mass change in the one or more subareas of the study area.

In some embodiments, the apparatus 200 comprises means, such as e.g. a storing module 404, adapted to store a determined mass change e.g. on a memory 405 comprised in or associated with the apparatus 200. The storing module 404 may be implemented as the processor 407 of the apparatus 200.

The storing module 404 may be configured to store the determined mass change in a Tagged Image File Format (TIFF), such as a GeoTIFF format. The apparatus 200 may also comprise means for storing data such as data comprised in or related to the datasets and/or the correction functions described herein or such as the data comprised in or related to parts of the datasets and/or parts of the correction functions described herein. As mentioned above, in some embodiments, the apparatus 200 comprises the memory 405 configured to store the data. The data may be processed or non-processed data and/or information relating thereto. The memory 405 may comprise one or more memory units. Further, the memory 405 may be a computer data storage or a semiconductor memory such as a computer memory, a read-only memory, a volatile memory or a non-volatile memory. The memory is arranged to be used to store obtained information, data, configurations, and applications etc. to perform the methods herein when being executed in the apparatus. The data may be stored in a Tagged Image File Format (TIFF), such as GeoTIFF.

TIFF is a computer file format for storing raster graphics images. The TIFF format is widely supported by image-manipulation applications, by publishing and page layout applications, and by scanning, faxing, word processing, optical character recognition and other applications. The format has tags, that is fields that may be coded with metadata, and include "Private" tags that can be used by any organization, like the GeoTIFF tags. GeoTIFF is a public domain metadata standard, using such private tags, which allows georeferencing information to be embedded within a TIFF file. The potential additional information includes map projection, coordinate systems, ellipsoids, datums, and everything else necessary to establish the exact spatial reference for the file. Thus the TIFF format can also be used with similar, proprietary additions, or additions can be made within the GeoTIFF format using "Private User Implementations" keys, identifying such data.

An alternative is the ^*.tfw World File sidecar file format where the data is stored in a separate file which may sit in the same folder as the regular TIFF file to provide a subset of the functionality of the standard GeoTIFF described here.

Embodiments herein for determining mass change in a study area of the earth using remote sensing may be implemented through one or more processors, such as the processor 407 in the arrangement depicted in Fig. 4, together with computer program code for performing the functions and/or method actions of embodiments herein. The program code mentioned above may also be provided as a computer program product, for instance in the form of a data carrier carrying computer program code for performing the embodiments herein when being loaded into the in the apparatus 200. One such carrier may be in the form of an electronic signal, optical signal, radio signal or a non-transitory computer-readable medium such as a computer readable storage medium. The computer readable storage medium may be a CD ROM disc or a memory stick.

The computer program code may furthermore be provided as pure program code on a server and downloaded to the apparatus 200.

Those skilled in the art will also appreciate that the retrieving module 401 , the determining module 402, the displaying module 403, and the storing module 404 above may refer to a combination of analog and digital circuits, and/or one or more processors configured with software and/or firmware, e.g. stored in the memory, that when executed by the one or more processors such as the processors in the apparatus 200 perform as described above. One or more of these processors, as well as the other digital hardware, may be included in a single application-specific integrated circuitry (ASIC), or several processors and various digital hardware may be distributed among several separate components, whether individually packaged or assembled into a system-on-a-chip (SoC).

When using the word "comprise" or "comprising" it shall be interpreted as non-limiting, i.e. meaning "consist at least of". The embodiments herein are not limited to the above described preferred

embodiments. Various alternatives, modifications and equivalents may be used. Therefore, the above embodiments should not be taken as limiting the scope of the invention, which is defined by the appending claims.

EMBODIMENT FOR THE STUDY AREA UGANDA

The embodiment demonstrates how a correction function, Diffp_ene. configured to correct for the differences in penetration between a first and a second wavelengths, and a second correction function, Difft_ech. configured to correct for technical errors, are used in a study area being the country of Uganda. The mass change that was determined was in the field of forest carbon monitoring. Forest height changes over time were used to estimate change in forest carbon stock. The height changes were derived from two data sets from two different points in time covering the study area. The forest heights were derived from interferometric SAR, i.e. an advanced radar sensor system used to take radar images. Both forest height decreases from logging and forest height increases from growth were seen in the Uganda study area. On average the height decreased by 2.8 cm per year during 2000- 2012. This decrease corresponds to an annual C02 emission having a mass of 20.7 million tonnes.

The apparatus, 200, used to determine the correction functions, calculate the mass changes and display it, was operated by NIBIO (the Norwegian Institute Bioeconomy). The retrieving module 401 was connected using an internal system computer bus of the SATA Express type to memory 405 of the SSD (Solid State Drive) type. Using SSD as memory rather than an ordinary magnetic hard disk, greatly increases the processing speed of the apparatus, 200. The same bus was used to interconnect the 403 Displaying module that was implemented using both a video card with a separate Graphics Processing Unit (GPU) in addition to an Integrated Graphics Processor unit (IGP). The Processor 407 was thus implemented with a GPU, so that it connected with the system bus to an Intel processor of the Core i7 type, which includes an IGP. The 404 Storing module was implemented using Random-Access Memory (RAM), whereas the 402 Determining module and 401 Retrieving module were implemented as part of the Intel processor. As the Intel Core i7 processor has several cores, the processing load was distributed between the cores, the IGP and the GPU. The display 406 was an ordinary computer display with resolution 1920 ^* 1200 pixels. The Input/Output interface 400 used to download the data sets into the apparatus 200 was an 1 Gbit/s Ethernet interface using the ordinary Internet protocol TCP/IP.

One data set was acquired around the year 2012 (Tandem-X) and another in the year

2000 (SRTM). The data were obtained with radar sensors, or Synthetic Aperture Radar (SAR), creating radar images; and the Digital Surface Models, (DSMs) were produced using interferometric processing (InSAR). These core data were supplemented with auxilliary data on land cover types, Landsat-based forest cover estimates and their changes over time, water-masks as well as external models obtained from other studies on the conversion between variables.

Both the X- and C-band InSAR DSMs from SRTM were used, and acquired in February 2000 from the space shuttle, 102, 102a, 102b. The SRTM mission acquired bi-static data by mounting SAR sensors both in the space shuttle body, 104, and on the tip of a 60 m extended boom, 106. The C-band data were acquired with a scanSAR system of the Spaceborne Imaging Radar-C (SIR-C) hardware, while the X-band data were acquired with a stripmap system. The C-band data had full areal coverage, while the X-band provides partly coverage (-30%) only, making up an X-shaped coverage pattern over Uganda. The data were in geographic (lat/lon) projection in meters above sea level with a spatial resolution of 1 arc sec. Data were downloaded as Digital Elevation Model (DEM) tiles from the DIVA-GIS repository (diva-gis.org), 108. The C-band data were downloaded from United States

Geological Survey (USGS), 108, and the X-band from the German Aerospace Centre (DLR), 108. The tiles were mosaicked together, using actions 301 and 302, to one C-band and one X-band DSM file over Uganda, and elevation data converted to ellipsoidal heights using ENVI / Sarscape software.

A particular issue with obtaining height changes from SRTM C-band to Tandem-X X- band is that they used different wavelengths, and the penetration depth into vegetation is influenced by wavelength. The wavelength is 3.1 cm for X-band and 5.6 cm for C-band. The longer the wavelength the deeper is the penetration into a vegetation volume, which means that the radar echo from C-band will be somewhat deeper down in the canopy than the X- band echo. Hence, before comparing the two data sets, and extract height changes over time, the C-band DSM were corrected, i.e. lifted upwards to a simulated X-band DSM, using a correction function as in action 304, so that no correction shall occur in treeless areas, while it should be a maximum correction in dense forests.

From the SRTM C band DEM a slope map was derived to discard extreme values. InSAR data is sensitive to steep slope angles that cause random extreme values due to irregularities in the signal reflection. The slope data set was employed to support the exclusion of false information in further processing steps.

USGS provided an open source raster product derived from MODIS data between 2001 and 2010. Since the data product had a global coverage with 17 land cover types, the spatial resolution was with 15 arc seconds respectively low. Further processing steps required equal extent and resolution as found for the TanDEM-X DSM and SRTM C-band data sets. Therefor the land cover type data set (LCT) was resampled to 1 arc second, using 'nearest neighborhood' and clipped to the rectangular extent of the SAR products. Like the prior data sets, also the LCT data was downloaded in latitude/longitude WGS84 geographic projection.

The 17 land cover types, 10 of which were found in Uganda, were further reduced to 8 types by combining land cover types of small extent.

The country boundaries, that were the boundaries of the study area, as well as the water mask were provided by the open source spatial GIS data platform Diva-GIS.org. Both were vector data sets in geographic projection WGS84.

The University of Maryland provided the Hansen Global Forest Change data set, which was a raster product derived from a Landsat 7 ETM+ time series between 2000 and 2012. Each tile covered an area of 10 by 10 degrees on a global scale. The extent of Uganda was covered with four tiles. Three different data products were downloaded, containing information about the tree canopy cover in 2000 from 0 to 100 % for all trees with a height more than 5m, the forest gain and forest lost from 2000 to 2012 as a 8-bit tiff file where 0 represents no gain (or no loss respectively) and 1 represents gain (or loss). All tiles were provided in geographic projection (WGS84) and could be downloaded in the 'download_v1.0' version with 1 arc second spatial resolution.

In this embodiment for the study area Uganda, actions 301 to 306 were performed in the following manner:

The first part of the processing in apparatus 200 was to extract the forest height difference over time. The difference between the Tandem-X DEM and the SRTM C-band DEM contains three different sources:

a) Real height changes over time (logging, forest growth)

b) Errors, mainly or only errors in the SRTM C-band DEM

c) The penetration difference of X- and C-band into the forest canopy

The processing steps were as follows:

Forest height change from year 2000 to year 2012 derived as ΔΗ1 = DEMTDX

- DEMSRTM

Processing step 1 , Correcting errors (artifacts) of the SRTM C-band DEM, obtaining AHcorrl

Checking for remaining bias and ramp: No correction needed.

Processing step 2, Correcting for C-band to X-band penetration difference, including removal of errors (artifacts) of the SRTM X-band DEM, obtaining AHcorr2

Processing step 3, Filling of voids and illegal, extreme values, obtaining AHcorr3

The second part was to recalculate the height changes to biomass and carbon stock changes. Processing step 4: Recalculating AHcorr2 to Above Ground Biomass (AGB) change, AAGBcorr2

Processing step 5: Recalculating AAGBcorr2 to Δ0Ο2 emissions Processing step 6: Uncertainty description and estimation

All data sets were available with geographic coordinates, referring to WGS 84 ellipsoid, with 1 arc second spatial resolution and pixel alignment, - a few of them after resampling from other coordinate systems or from orthometric heights. The resampling was carried out with ENVI / Sarscape and ArcMap software. The data sets were clipped to a rectangular extent covering the study area Uganda. All raster data sets were transformed to tables, and statistical analyses were carried out using SAS (Statistical Analyses Software).

The errors occurring as line and band artifacts were clearly visible in the difference between the Tandem-X DEM and the SRTM C-band DEM, and were assumed to be errors in the SRTM C-band DEM. These artifacts were successfully removed using an Analysis of Variance (ANOVA) model to make a tailored filter, and by subtracting this a difference image was obtained representing real height changes. The artifacts occurred as lines and belts having fixed angles. Every pixel was assigned to four variables L1 ; L2; B1 and B2 representing these artifacts: L1 was a line number with a 330.2° Southwest - Northeast direction; L2 was a line number with a 29.6° Southeast - Northwest direction; and B1 was a belt number perpendicular to L1 and B2 was a belt number perpendicular to L2. It is notable that these artifacts were mixed up with real changes, and hence, pixels were excluded having apparent real height changes from the ANOVA: Pixels having height changes larger than 10 m, as well as pixels classified as 'gain' or 'loss' in the Global Forest Cover change dataset based on Landsat data were excluded. This left about 400 million pixels which were entered into the ANOVA model:

(1) ΔΗ^ βο + Ll + L2 + B1 + B2 + e where β₀ was an intercept and e was the residual, which equals AH_COIXi, i.e. the corrected height difference between the Tandem-X DEM and the SRTM C-band DEM. Action 305 was performed such that Diff_tech = β,, + L1 + L2 + B1 + B2= ΔΗι - e of equation (1 ).

Pixels were selected having no forest cover and no forest cover change evenly around Uganda, and the mean and slope of these pixels were calculated. This was to possibly make a final large-scale correction of the AH_conl. A model was developed for a height correction between the C-band and the X-band DEMs, representing the difference in penetration depth into the forest canopy. First, it was assumed that the artifacts estimated in step 1 , above, were errors in the C-band DEM, and then they were used to correct the C-band DEM. The height difference between the X-band DEM and the corrected C-band DEM has two components, - penetration differences and artifacts of the X-band DEM. It was crucial to estimate these two sources of difference simultaneously with one model. The model type used was based on the assumption that the penetration difference increase from zero in non-forest pixels linearly with the density of the forest canopy, and forest cover from the GFC data set as a proxy for canopy density was used. Altogether, the X-band artifact removal and penetration difference were estimated with a GLM model:

(2) AH_corrl = Po + Px_Li + FC*Pj where β₀ was an intercept, βχπ was fixed effect of line segments, XL, crossing the belts of S TM X-band data, and + FC*P_j was the penetration difference depending on land cover type j and forest cover (%). Specific correction models were made for 8 land cover types, i.e. Evergreen Broadleaf Forest; Woody Savanna; Savannas; Grasslands; Permanent Wetlands; Croplands; Urban and Built-Up; and Cropland/Natural Vegetation Mosaic. With these models the corrections varied from zero in no-forest pixels up to some decimeters in densely forested pixels.

Action 304 was performed such that Diff_pene = FC*P_j in equation (2), and this equals ΔΗ_∞ΓΓ1 - βο + βχϋ-

In addition to the actions previously described, an action was performed to replace void areas in the SRTM C-band data, as well some extreme values outside the range of plausible values. Such pixel values were discarded and replaced them with plausible values. First, the pixels to be discarded were identified in steep areas. The occurrence of such pixels were found mainly in steep terrain and near water bodies indicating inaccuracy in the outline of water bodies. A careful and visual examination of these occurrences revealed that pixels in terrain steeper than about 30° had unreliable height change values. Pixels above the 99% quantile of slopes were discarded, i.e. pixels having a slope > 29.6 degrees. Furthermore, pixels were discarded for height change values exceeding +-25 m, which made up 0.56 %o of the observations. The +-25 threshold was found after careful examination of extreme values. Large and real changes were seen as clear-cut areas, fitting well with the 'loss' category in GFC changes, and as well in some cases forest growth was evident e.g. following exactly the outline of protected areas. The maximum height changes found in these cases were ±20, while a few cases of height changes up to ±25 m were found.

These pixels were replaced in the following way: They were assigned with the mean height change of pixels in the same land cover type and the same Global Forest Cover (GFC) change category. This means for example that a given void pixel occurring in a woody savannah and tagged as 'forest loss' in GFC was assigned the mean value of all woody savannah pixels tagged as 'forest loss' in GFC.

(3) Δ¾_0Ιτ2 = Po + GL*LCT + e,

where βο was an intercept, GL^*LCT was the interaction between the gain-loss category (-1 , 0, 1 ) and the land cover type. In this model pixels were discarded having a slope steeper than 30° or having a height change more than 25 m.

Having decided the correction functions Diff_pene and Diff_teCh,for this study area of Uganda, the apparatus 200 can then be used to determine mass changes for other study areas, with similar characteristics for landscape types and similar coverage of remote sensing, such as other study areas in Africa. If a study area differs, the correction functions best suited for that area can be calculated and the apparatus 200 of this example can be used with those correction functions.

Claims

A method for determining mass change in a study area of the earth using remote sensing data, the method comprises:

- retrieving (301) a first three-dimensional, 3D, dataset of the study area, which first 3D dataset is obtained by means of remote sensing at a first point of time and comprises a first set of pixel data given by a first wavelength of the remote sensing, wherein each pixel comprised in the first set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the first set of pixel data is indicative of an amount of mass in the subarea at the first point of time;

- retrieving (302) a second 3D dataset of the study area, which second 3D dataset is obtained by means of remote sensing at a second point of time and comprises a second set of pixel data given by a second wavelength of the remote sensing, wherein each pixel comprised in the second set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the second set of pixel data is indicative of an amount of mass in the subarea at the second point of time; and

- by means of a first correction function, Diffp_ene_> configured to correct for the differences in penetration between the first and second wavelengths, and by means of a second correction function, Difftech. configured to correct for technical errors, determining (306) mass change in one or more subareas of the study area as a corrected difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas. 2. The method of claim 1 , further comprising:

- determining (304) the first correction function, Diffp_ene. which first correction function, Diffp_ene. comprises a third set of pixel data, wherein each pixel comprised in the third set of pixel data corresponds to a subarea of the study area, and wherein a data value for each pixel comprised in the third set of pixel data is indicative of a difference in penetration between the first wavelength and the second wavelength; and

- determining (305) the second correction function, Difftech. which the second correction function, Difftech. comprises a fourth set of pixel data, wherein each pixel comprised in the fourth set of pixel data corresponds to a subarea of the study area, and wherein a data value for each pixel comprised in the fourth set of pixel data is indicative of a technical error.

The method of claim 1 or 2, further comprises:

- determining (303) a third 3D dataset, Diffi, as a difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas; and

wherein determining (306) mass change in one or more subareas of the study area as the corrected difference further comprises:

- correcting the data values of the determined third 3D dataset, Diffi , by means of the first correction function, Diffp_ene_> and the second correction function, Difftech. whereby the corrected data values of the determined third 3D dataset, Diffi , are indicative of the mass change in one or more subareas of the study area.

The method of claim 3, wherein the mass is biomass, the remote sensing is radar imaging, the first and second 3D datasets are respective first and second radar images, and wherein determining (305) the second correction function, Difftech. further comprises:

- retrieving a fourth 3D dataset of the study area, which fourth 3D dataset is obtained by means of passive remote sensing at the first point of time and comprises a fifth set of pixel data, wherein each pixel comprised in the fifth set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the fifth set of pixel data is indicative of forest cover in the subarea at the first point of time;

- retrieving a fifth 3D dataset of the study area, which fifth 3D dataset is obtained by means of passive remote sensing at the second point of time and comprises a sixth set of pixel data, wherein each pixel comprised in the sixth set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the sixth set of pixel data is indicative of forest cover in the subarea at the second point of time;

- by means of the fifth and sixth sets of pixel data, determining one or more pixels in the first radar image which one or more pixels comprise less than a predetermined forest cover and which one or more pixels have no change in forest cover between the first and second points of time;

- when one or more data values indicative of an artefact are comprised in the third 3D dataset, Diffi , removing the one or more values indicative of the artefact from the third 3D dataset, Diffi ;and

- determining the second correction function, Difftech. as an interpolation of the third

3D dataset, Diffi

The method of claim 4, wherein determining (304) the first correction function, Diffp_ene_> further comprises:

- selecting an image area of the first radar image taken at the first point of time;

- determining a sixth 3D dataset, Diffx.c, comprising the difference between data values of corresponding one or more pixels of a seventh set of pixel data and the first set of pixel data, wherein the seventh set of pixel data is given by the first wavelength of the radar imaging at the first point of time, and wherein corresponding one or more pixels correspond to the image area; and

- determining the first correction function, Diffpene. by correcting the sixth 3D dataset, Diffx-c, for technical errors using the second correction function, Difftech-

The method of claim 4 or 5, wherein determining (304) the first correction function, Diffpene. further comprises:

- determining the first correction function, Diff_pene. as approximately 1 -3 cm per percentage increase in forest cover, especially as approximately 2 cm per percentage increase in forest cover, wherein the forest cover is given by the fifth set of pixel data.

The method of any one of claims 1 -6, further comprising:

- displaying (307) an image of the determined mass change in the one or more subareas of the study area.

The method of any one of claims 1 -7, further comprising:

- storing (308) the determined mass change in a Tagged Image File Format, TIFF, such as a GeoTIFF format.

An apparatus (200) for determining mass change in a study area of the earth using remote sensing data, wherein the apparatus is configured to:

- retrieve a first three-dimensional, 3D, dataset of the study area, which first 3D dataset is obtained by means of remote sensing at a first point of time and comprises a first set of pixel data given by a first wavelength of the remote sensing, wherein each pixel comprised in the first set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the first set of pixel data is indicative of an amount of mass in the subarea at the first point of time; - retrieve a second 3D dataset of the study area, which second 3D dataset is obtained by means of remote sensing at a second point of time and comprises a second set of pixel data given by a second wavelength of the remote sensing, wherein each pixel comprised in the second set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the second set of pixel data is indicative of an amount of mass in the subarea at the second point of time; and

- by means of a first correction function, Diffpene. configured to correct for the differences in penetration between the first and second wavelengths, and by means of a second correction function, Difft_ech. configured to correct for technical errors, determine mass change in one or more subareas of the study area as a corrected difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas.

10. The apparatus (200) of claim 9, further being configured to:

- determine the first correction function, Diffpene_* which first correction function,

Diffpene. comprises a third set of pixel data, wherein each pixel comprised in the third set of pixel data corresponds to a subarea of the study area, and wherein a data value for each pixel comprised in the third set of pixel data is indicative of a difference in penetration between the first wavelength and the second wavelength; and

- determine the second correction function, Difftech. which the second correction function, Difftech. comprises a fourth set of pixel data, wherein each pixel comprised in the fourth set of pixel data corresponds to a subarea of the study area, and wherein a data value for each pixel comprised in the fourth set of pixel data is indicative of a technical error.

1 1. The apparatus (200) of claim 9 or 10, further being configured to:

- determine a third 3D dataset, Diff-i , as a difference between data values of corresponding one or more pixels of the second set of pixel data and the first set of pixel data, which corresponding one or more pixels correspond to the same one or more subareas; and

to determine mass change in one or more subareas of the study area as the corrected difference the apparatus is further configured to:

- correct the data values of the determined third 3D dataset, Diff-i , by means of the first correction function, Diffp_ene_> and the second correction function, Difft_ech_> whereby the corrected data values of the determined third 3D dataset, Diff-i , are indicative of the mass change in one or more subareas of the study area.

12. The apparatus (200) of claim 1 1 , wherein the mass is biomass, the remote sensing is radar imaging, the first and second 3D datasets are respective first and second radar images, and wherein the apparatus in order to determine the second correction function, Difftech. further is configured to:

- retrieve a fourth 3D dataset of the study area, which fourth 3D dataset is obtained by means of passive remote sensing at the first point of time and comprises a fifth set of pixel data, wherein each pixel comprised in the fifth set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the fifth set of pixel data is indicative of forest cover in the subarea at the first point of time;

- retrieve a fifth 3D dataset of the study area, which fifth 3D dataset is obtained by means of passive remote sensing at the second point of time and comprises a sixth set of pixel data, wherein each pixel comprised in the sixth set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the sixth set of pixel data is indicative of forest cover in the subarea at the second point of time;

- by means of the fifth and sixth sets of pixel data, determine one or more pixels in the first radar image which one or more pixels comprise less than a predetermined forest cover and which one or more pixels have no change in forest cover between the first and second points of time;

- when one or more data values indicative of an artefact are comprised in the third 3D dataset, Diff _: remove the one or more values indicative of the artefact from the third 3D dataset, Diff^and

- determine the second correction function, Difftech. as an interpolation of the third 3D dataset, Difft

13. The apparatus (200) of claim 12, wherein the apparatus in order to determine the first correction function, Diff_pene. further is configured to:

- select an image area of the first radar image taken at the first point of time;

- determine a sixth 3D dataset, Diffx.c, comprising the difference between data values of corresponding one or more pixels of a seventh set of pixel data and the first set of pixel data, wherein the seventh set of pixel data is given by the first wavelength of the radar imaging at the first point of time, and wherein corresponding one or more pixels correspond to the image area; and

- determine the first correction function, Diffp_ene_> by correcting the sixth 3D dataset, Diffx-c, for technical errors using the second correction function, Difftech-

14. The apparatus (200) of claim 12 or 13, wherein the apparatus in order to determine the first correction function, Diffpene. further is configured to:

- determine the first correction function, Diffpene. as approximately 1 -3 cm per percentage increase in forest cover, especially as approximately 2 cm per percentage increase in forest cover, wherein the forest cover is given by the fifth set of pixel data.

15. The apparatus (200) of any one of claims 9-14, further being configured to:

- display an image of the determined mass change in the one or more subareas of the study area.

16. The apparatus (200) of any one of claims 9-15, further being configured to:

- store the determined mass change in a Tagged Image File Format, TIFF, such as a GeoTIFF format.

17. A remote sensing system (100) for determining mass change in a study area of the earth using remote sensing data, wherein the remote sensing system (100) comprises:

at least one transport means (102,102a, 102b) comprising at least one sensor (104,106) and configured to transport the at least one sensor (104,106) at a distance over the study area;

- wherein the at least one sensor (104,106) is configured to acquire remote sensing data by, at a certain point of time, emit an electromagnetic radiation (E) of a predefined wavelength towards the study area of the earth and receive backscattered

electromagnetic radiation (B) from the study area, whereby a 3D dataset is obtained for the certain point of time, which 3D data set comprises a set of pixel data for the predefined wavelength of the emitted electromagnetic radiation (E), wherein each pixel comprised in the set of pixel data corresponds to a subarea of the study area and wherein a data value for each pixel comprised in the set of pixel data is indicative of an amount of mass in the subarea at the certain point of time; wherein the remote sensing system (100) further comprises:

- a data storage means (108) configured to store the 3D dataset acquired by the at least one sensor (104,106); and - an apparatus according to any one of claims 9-16.

18. A computer program, comprising instructions which, when executed on at least one processor, causes the at least one processor to carry out the method according to any one of claims 1 -8.

19. A carrier comprising the computer program of claim 18, wherein the carrier is one of an electronic signal, an optical signal, a radio signal or a computer readable storage medium.