Methods and systems for detecting heat emission sources with nighttime satellite data

Abstract

The disclosure relates to a computer-implemented method of detecting heat emitting sources in nighttime satellite data. It assumes an independent and identical Gaussian distribution on the background data and looks for parts of the images with abnormally high values. For this, the mean and variance of the background model is estimated. Then, a region growing algorithm is used to find connected regions with high values. Finally, a statistical test is used to decide whether the sum of the values of each region is significantly higher than expected on the background model. Only regions detected in the two short-wave infrared bands are validated. The test level is selected in order to control the number of false detections. Compared to classical pixel-based methods, the method allows the detection of hotspots with lower radiance while keeping a low commission error rate.

Description

[0001] Methods and systems for detecting heat emission sources with nighttime satellite data

[0002] TECHNICAL FIELD

[0003] The field of the invention is that of Earth remote sensing methods and systems for detecting heat emission sources such as active fires and heat-emitting industrial facilities.

[0004] DESCRIPTION OF RELATED ART

[0005] Satellite data has been used to detect and monitor heat emission sources (also known as hotspots) since the early 1980s. One prominent application is the mapping of active fires, notably using spectral bands specifically designed to identify fires and estimate their emitted thermal energy. Many other applications rely on hotspot detection, including the monitoring of volcanic activity, the monitoring of industrial activity and the detection of gas flares.

[0006] The identification and characterization of hot sources is greatly facilitated in the absence of sunlight. Without sunlight, emitting objects on the ground can indeed be detected when their radiance exceeds the noise floor of the sensor. Hence, most existing methods use a simple threshold to detect hotspots at night. A value typically used is 1 W m-2sr1pm-1.

[0007] Most methods add complementary tests to remove commission errors. However, it is worth noting that this radiance level is the minimum detectable using these methods.

[0008] One exception is the HOTMAP method from Murphy et al., "HOTMAP: Global hot target detection at moderate spatial resolution" Remote Sensing of Environment, vol. 177, pp. 78-88, May 2016. For the detection of hotspots in nighttime data, this method relies on two thresholds applied to the SWIR band B7. The first threshold, at 1 W m“2sr-1pm-1, is used to detect "obviously hot pixels". A second threshold, at 0.015 m-2sr1pm-1, is applied to detect "candidate pixels." Those candidate pixels are validated only if they are juxtaposed to another candidate or obviously hot pixel. BRIEF DESCRIPTION OF THE INVENTION

[0009] The invention aims at enabling the detection of hotspots having lower radiance, hence increasing the number of detected hotspots, while maintaining a similarly high precision.

[0010] To this purpose, the invention proposes a computer-implemented method of detecting heat emitting sources in satellite data, comprising the steps of: obtaining a first night-time satellite image of a region of interest acquired in a first spectral band; obtaining a second night-time satellite image of the region of interest acquired in a second spectral band different from the first spectral band; for each image among the first and second night-time satellite image, determining a hot sources detection mask by: o determining a distribution of a background model of the image; o based on the determined distribution, separating pixels of the image between dark pixels and bright pixels; o identifying regions of contiguous bright pixels within the image; o for each identified region:

[0011] ■ judging whether or not the identified region belongs to the background based on comparing the identified region and the determined distribution; and

[0012] ■ adding the identified region to the hot sources detection mask when it is judged that the identified region does not belong to the background; determining heat emitting regions as a result of intersecting the hot sources detection mask determined for the first night-time satellite image with the hot sources detection mask determined for the second night-time satellite image.

[0013] Certain preferred, but non-limiting aspects of the method are as follows: it further comprises determining a source location of a determined heat emitting region; determining a source location of a determined heat emitting region comprises determining a point of the determined heat emitting region having the highest cumulative brightness over the first and second night-time satellite image; it further comprises outputting a map of the region of interest having a visual indicator at each source location of a determined heat emitting region; intersecting the hot sources detection mask determined for the first night-time satellite image with the hot sources detection mask determined for the second night-time satellite image produces a final detection mask and the method further comprises filtering the final detection mask so that each heat emitting region in the filtered final detection mask comes from a different identified region in one of the hot sources detection masks; filtering the final detection mask comprises validating, among heat emitting regions resulting from the intersection of the hot sources detection mask determined for the first night-time satellite image with a single identified region in the hot sources detection mask determined for the second night-time satellite image, only the largest region; identifying regions of contiguous bright pixels within the image implements a greedy region algorithm; the determined distribution comprises a central intensity, the dark pixels being pixels having an intensity lowerthan the central intensity and the bright pixels being pixels having an intensity higher than the central intensity; for each identified region, judging whether or not the identified region belongs to the background based on comparing the identified region to the determined distribution comprises: o calculating a sum of bright pixel values of the bright pixels of the identified region; o determining a probability that a sum of background pixel values is higher than the calculated sum of bright pixel values, wherein the sum of background pixel values is derived from the determined distribution considering as many background pixels as the number of bright pixels in the identified region; o comparing the determined probability to a threshold; o when the determined probability is lower to the threshold, judging that the identified region does not belong to the background. the threshold is set based on a size of the identified region and on a size of the image; determining a distribution of a background model of the image comprises determining a background standard deviation; determining a background standard deviation comprises calculating an image biweight midvariance; the first and second spectral bands are shortwave infrared bands; the first night-time satellite image and the second night-time satellite image are acquired by a same imager onboard a satellite.

[0014] BRIEF DESCRIPTION OF THE DRAWINGS

[0015] Other aspects, aims, advantages and features of the invention will better appear upon reading the following detailed description of preferred embodiments thereof, provided as a non-limiting example, and done in reference to the appended drawings, in which:

[0016] - Figure 1 shows steps of a method according to a possible embodiment of the invention;

[0017] - Figure 2 shows example detection on a Landsat 9 product with (a) an image acquired in band B6, (b) the hot sources detection mask determined based on the image acquired in band B6, (c) an image acquired in band B7, (d) the hot sources detection mask determined based on the image acquired in band B7 and (e) the final detection resulting from the intersection of the two hot sources detection masks and including the filtering of false positives. DETAILED DESCRIPTION OF THE INVENTION

[0018] The present disclosure relates to a computer-implemented method of detecting heat emitting sources (also referred to as hotspots or hot sources) in satellite data. The heat emission sources may be active fires (such as gas flares, wildfires, agricultural burning or campfires) or heat-emitting industrial facilities (such as cement plants, steel mills, power plants or Liquified Natural Gas plants).

[0019] With reference to figure 1, this method comprises a step OBT1 of obtaining a first night-time satellite image of a region of interest acquired in a first spectral band and a step OBT2 of obtaining a second night-time satellite image of the region of interest acquired in a second spectral band different from the first spectral band. An example of such first and second night-time satellite images are shown in (a) and (c) on figure 2.

[0020] The region of interest may be an oil and gas-producing region, the heat emitting sources being gas flares located in the oil and gas-producing region. In another embodiment, the region of interest may be an industrial region, the heat emitting sources being heat-emitting industrial facilities located in the industrial region.

[0021] The first and second spectral bands may be shortwave infrared bands. The first nighttime satellite image and the second night-time satellite image may be acquired by a same imager onboard an Earth observing satellite. For instance, the first night-time satellite image and the second night-time satellite image are satellite images acquired by an Operational Land Imager (OLI) aboard a Landsat 8 or Landsat 9 satellite. The first spectral band may be Band 6 (SWIR 1) which has a spectral bandpass of 1.57-1.65 pm and the second spectral band may be Band 7 (SWIR 2) which has a spectral bandpass of 2.11-2.29 pm.

[0022] The disclosure makes the assumption that in the absence of hot sources, the pixel values in each of the first and second night-time satellite image are Independent and Identically Distributed (HD), following a statistical distribution such as a Gaussian distribution. While this is approximately true for most part of the SWIR bands acquired at night, this assumption does not hold for pixels corresponding to hot sources, whose value depart from the distribution of the background and neighbours are correlated. This is precisely how hotspots are detected: by identifying statistically significant outliers in the background model.

[0023] Still with reference to figure 1, for each image u: R2-> R of size X x Y among the first and second night-time satellite images, the method comprises a step HSM1, HSM2 of determining a hot sources detection mask which identifies groups of pixels within the image which correspond to hotspots. The hot sources detection mask is of size X x Y and may have white pixels for the hotspots and black pixels for the background pixels. Examples of hot spot detection masks determined for the first and the second image are shown on figure 2 in (b) and (d) respectively.

[0024] The hot sources detection mask determination is performed by first determining a distribution of a background model of the image, then extracting candidate regions and finally validating the candidate regions according to a statistical test.

[0025] Background estimation

[0026] The background pixels follow a stochastic model in which their values are HD random variables with a statistical distribution. Taking for instance the assumption of a Gaussian distribution, then determining a distribution of a background model of the image comprises estimating two parameters, a (background) central value (such as the mean / z) and a (background) spreading value (such as the standard deviation <J).

[0027] In the context of the disclosure, robust estimators may be required, able to provide reliable estimations even in the presence of outliers, namely the hotspots. In a possible embodiment, the median of the values of u is used as an estimation of / z and the biweight midvariance is used to estimate <J .

[0028] Be m the median of the image pixel values: m = medianfiz .

[0029] Then, the median absolute deviation (MAD) is defined as

[0030] MAD = medianfliij — m|}.

[0031] The biweight midvariance is defined as: where n is the total number of samples, here XY. The biweight midvariance is a robust estimator of the variance and V / ca nbe used as an estimator for <7. Candidate regions

[0032] Based on the determined distribution of the background model, the pixels of the image can be separated between dark pixels and bright pixels. In a possible embodiment, dark pixels are those pixels of the image which value is lower than the background central value (the image median m for instance) and bright pixels are those pixels of the image which value is larger than the background central value (the image median m for instance).

[0033] Then, extracting candidate regions from the image can comprise identifying regions of contiguous bright pixels within the image.

[0034] Hot sources can take any shape on the image. The smallest one appears as just one pixel, while large ones (for example wildfires) can spread into any shape. In an embodiment, identifying regions of contiguous bright pixels within the image implements a greedy region algorithm.

[0035] The greedy region growing algorithm may focus on groups of pixels connected under 4-connectivity, in which a pixel is connected to the pixels at coordinates (x ± 1, y) and (x, y ± 1). The shapes connected under 4-connectivity correspond to the figures called polyominoes. An exhaustive evaluation of all polyomino configurations on the image is not possible given its colossal number. This is why instead, the greedy region growing algorithm is used.

[0036] A possible implementation of this algorithm is as follows. First, the pixels are sorted by pixel value. They are considered in decreasing order, so the brightest is considered first. If the pixel value is largerthan the image median, a new region is started, initially consisting of the seed pixel. Then, the four 4-connected neighbors of each pixel already in the region are considered, and the largest among them is selected provided that three conditions are satisfied: a) the pixel is not already in the region; b) it is a bright pixel (i.e., its value is larger than the image median); c) its value is smaller or equal to the last pixel added to the region. The last condition is a simple way of stopping the region growing without spreading on the background or growing into a neighboring hot source. Indeed, a hot source is roughly represented as a cone, with a brightest pixel and the pixels around gradually decreasing values until they fuse to the background. At each step of the region growing process, the partially obtained region of pixels is evaluated using the statistical test described below. If the test validates the hot source detection, then the region is added to an output map (i.e. the hot sources detection mask), and the process continues. This can lead to a hot source being detected several times with increasingly larger regions; this is merged in the output map as the largest region will contain all the others.

[0037] Finally, when no new pixel is found satisfying the conditions, the region growing process stops. At this point, the pixel values of all the pixels included in the region can be set to a small value (for example the median value) to prevent them from being used again as seed pixels or being included in any further region. This results in a greedy and fast algorithm. Then, the process starts again from the next pixel in the list, provided that its value is larger than the image median.

[0038] Statistical test

[0039] For each region r in the image identified as a region of |r| contiguous bright pixels, it is then judged whether or not the identified region belongs to the background based on comparing the identified region and the determined distribution and adding the identified region to the hot sources detection mask when it is judged that the identified region does not belong to the background.

[0040] Given a region r of size |r|, it is first decided whether the pixel values in r belong to the background or to a hotspot. To this purpose, in a possible embodiment, first a sum srof the values of the bright pixels of the identified region is calculated: sr= iErui- Then, a sum of values of background pixels is considered, for as many background pixels as the number |r| of bright pixels in the identified region r. This sum of values of background pixels is a random variable Srcorresponding to the sum of |r| Gaussian random variables of parameters / z and <J. If all the pixels in the region r follow the background model, then srshould follow the same distribution as Sr. Given the 11 D assumption in the background model, Srcan be determined from the distribution of the background model. For instance, it also follows a Gaussian distribution of parameters

[0041] Then, a probability P(Sr> sr) that the sum of background pixel values is higher than the sum of bright pixel values is determined. This probability can be computed using the "I * I error function. For instance, P(Sr> sr) = - erfc J> with erfc the complementary error function, erfc

[0042] The determined probability P(Sr> sr) is compared to a threshold and when the determined probability is lower to the threshold, it is judged that the region r does not belong to the background but belongs to a hotspot.

[0043] In a possible embodiment, the threshold is set based on a size |r| of the identified region r and the total number X x Y of pixels in the image, in order that an expected number of false detections is controlled.

[0044] According to the a-contrario methodology, the number of false detections depends on the total number of tests performed. Following this methodology, the Number of False Alarms (NFA) can be defined as NFA(sr, r) = Nr■ P(Sr> sr), provided that £reR— < 1, where R is the set of all the regions tested and Nris a weight allowing to distribute the risk of false detection among all the tests.

[0045] The a-contrario methodology prescribes accepting as valid detections the candidates with NFA < e for a predefined value e. It can be shown that under the background model J£o, the expected number of tests with NFA < is bounded by £. As a result, £ gives an a priori estimate of the mean number of false detections under J£o. In many practical applications, including the present one, the value £ = 1 can be adopted. Indeed, it allows for less than one false detection on an image, which is usually quite tolerable.

[0046] It is natural to think that R would correspond to the actual set of candidate regions evaluated. But doing this would complicate the formulation. The procedure used to select the candidates (the region growing process) implies that, even if it were applied to data following the background model exactly, the resulting regions would not. Indeed, the procedure selects large values and only values larger than the median; the Gaussian model would surely be violated.

[0047] A safer alternative is to consider a theoretical set R, taking into account all potential regions theoretically possible. That is, all the polyomino regions in the image. The exact number bkof polyomino configurations of given size k is not known in general, but there are good approximations of this number. In the present case, it is enough to use an estimate of the order of magnitude, so the following approximate formula is sufficient. It reads bk« a • (3k / k, where a « 0.316915 and / ? « 4.062570.

[0048] It has to be considered that each particular polyomino may be placed at any position in the image and also the different sizes of the connected regions. To consider these factors, it is possible to set

[0049] This weight is less strict on regions of small size. A simple computation shows that it satisfies

[0050] All in all, for each region r, the Number of False Alarms is computed as NFA(sr, r) = XY ■ b\r\ ■ 2|r|• P(Sr> sr) and a hot source is detected when NFA(sr, r) < 1, i.e., when

[0051] Bi-band detection

[0052] We have described so far the detection of groups of bright pixels in the first image and the detection of groups of bright pixels in the second image, resulting in a hot sources detection mask determined for the first image and in a the hot sources mask detection determined for the second image.

[0053] With reference to figure 1, the method then comprises a step ITS of determining heat emitting regions as a result of intersecting the hot sources detection mask determined for the first satellite image with the hot sources detection mask determined for the second satellite image. A heat emitting region is therefore composed of pixels each determined as being part of a region of contiguous bright pixels in both the first and the second images. A final detection mask is obtained at step ITS by intersecting the two hot sources detection mask determined at steps HSM1 and HSM2.

[0054] In other words, hotspots are detected independently in the two images and consistent detections are kept. This makes the method even more restrictive, as only detections that are made in both bands are validated. Thus, the expected number of false detections in remains inferior to e.

[0055] In an embodiment, the method further comprises determining a source location of a determined heat emitting region. Determining a source location of a determined heat emitting region can comprise determining a point of the determined heat emitting region having the highest cumulative brightness over the first and second satellite image.

[0056] Then the method may comprise outputting a map of the region of interest having a visual indicator at each source location of a determined heat emitting region.

[0057] The hot source is in general much smaller than the detected connected component. To identify it in the group of detected pixels, the brightest pixel in the sum of radiances of the two spectral bands can be looked for, and this pixel's center is used as location for the hotspot.

[0058] In addition, strong hotspots are typically surrounded by a bright halo. This can lead to false detections as some nearby but separated connected components can be identified as hotspots due to the halo contamination. In order to remove the false detections, the method may further comprise a step FILT of filtering the final detection mask so that each heat emitting region in the filtered final detection mask comes from a different identified region in one of the hot sources detection masks (such as the B7 mask). This filtering may consist in validating, among heat emitting regions resulting from the intersection of the hot sources detection mask determined for the first night-time satellite image (such as the B6 mask) with a single identified region in the hot sources detection mask determined for the second night-time satellite image (such as the B7 mask), only the largest region.

[0059] Indeed, the detected polyominoes are generally larger in B7 than in B6. Hence, the false positives in the final detection mask located in a hotspot's halo are generally part of a unique B7 connected component. A method of removing the false detections in the final detection mask therefore consists in validating only one of the connected components (a heat emitting region) when several originate from the same connected component (identified region) in the B7 detection mask. The validated connected component may be the largest one in the final detection mask. With this procedure, each detection in the filtered final mask comes from a different B7 connected component.

[0060] Figure 2 illustrates the steps of the proposed method with (a) a Landsat 9 image acquired in band B6, (b) the hot sources detection mask determined based on the Landsat 9 image acquired in band B6, (c) a Landsat 9 image acquired in band B7, (d) the hot sources detection mask determined based on the Landsat 9 image acquired in band B7 and (e) the final detection mask resulting from the intersection of the two hot sources detection masks and including the filtering of false positives. In (e), references Cl, C2 and C3 each designates a source location of a determined heat emitting region.

[0061] Experiments

[0062] Nighttime data acquired over the Permian Basin, a major oil and gas producing area in the United States with numerous gas flares, is considered. Landsat 8 and 9 both acquired data in this region from March to August 2023. Products on row 129 and path 206 are used. The temporal consistency of the proposed method is evaluated using this time series of 22 products, and compared to Schroeder et al., "Active fire detection using Landsat-8 / OLI data," Remote Sensing of Environment, vol. 185, pp. 210-220, nov 2016.

[0063] Given the low probability of a false detection occurring at the same location in two different acquisitions, it is assumed that any detection present in more than one acquisition is a true positive. Conversely, a detection found in a single acquisition is categorized as a false positive. With this definition, some detections classified as false positives can nonetheless correspond to actual valid detections, for example in the case of intermittent phenomena or because of clouds hiding the ground in several acquisitions. Hence, any hotspot found at the same location in at least one other product of the time-series is validated. The eight neighboring pixels were included to account for the slight inaccuracy in the localization of the hotspots, the maximum pixel position in the polyominoes being indeed not perfectly stable, notably because of the digital number folding.

[0064] The below table displays the evaluation results. method True Positives False Positives precision

[0065] Schroeder et a / . 5675 746 88.4%

[0066] Proposed method 18843 2598 87.9%

[0067] The proposed method detects over three times more hotspots than Schroeder et al.'s method while maintaining a similarly high precision level.

[0068] The primary source of false positives is in bright regions surrounding powerful flares, especially in the presence of clouds. Indeed, they reflect the light from the flares. These reflections sometimes appear disconnected from their source, so the removal of false detections in the hotspot halos may fail. Given the model, the detection from reflected light in the clouds are reasonable commission errors. To verify that the precision score increases in the absence of clouds, five acquisitions with the lowest cloud coverage (none is perfectly free of clouds) were manually selected and considered their detections only: their average precision is 92.0%. This confirms that clouds cause false positives.

[0069] The present disclosure presents a statistically based method for detecting hotspots, for instance in short-wave infrared data acquired at nighttime by Landsat 8 and 9 or Sentinel-2 (when nighttime data comes available). Thresholds are set automatically in order to control the expected number of false detections. Experiments show that, compared to classical techniques, the proposed approach greatly increases the number of hotspots that can be detected in nighttime images, while maintaining a similarly high precision.

[0070] In particular, the proposed method differs in three notable ways from the HOTMAP method from Murphy et al.. First, the parameters of the background model are estimated from the image itself, which makes the method more robust to changes in the parameters of the noise distribution between satellite products. Second, the proposed method considers regions of all sizes, which allows to use even lower thresholds as soon as the considered region is larger than two pixels. Third, the proposed method relies on two spectral bands (such as SWIR bands B6 and B7) and validates detections only if they are consistent between the two bands. This makes the proposed method robust to extreme noise events, which are not uncommon, especially in the South-Atlantic Anomaly region where the sensor is more exposed to very high energy particles. Furthermore, the proposed method has a single parameter with a clear signification: the expected number of false alarms one is ready to tolerate in an image.

[0071] The invention is not limited to the method disclosed above, but also relates to a system comprising a processor configured to perform said method, to a computer program product comprising instructions which, when the program is executed by a computer, cause the computer to carry out said method and to a non-transitory computer-readable medium storing a program including instructions that, when executed by a processor of a computing device, cause the computing device to carry out said method.

Claims

CLAIMS1. A computer-implemented method of detecting heat emitting sources in satellite data, comprising the steps of: obtaining (OBT1) a first night-time satellite image of a region of interest acquired in a first spectral band; obtaining (OBT2) a second night-time satellite image of the region of interest acquired in a second spectral band different from the first spectral band; for each image among the first and second night-time satellite image, determining a hot sources detection mask (HSM1, HSM2) by: o determining a distribution of a background model of the image; o based on the determined distribution, separating pixels of the image between dark pixels and bright pixels; o identifying regions of contiguous bright pixels within the image; o for each identified region:■ judging whether or not the identified region belongs to the background based on comparing the identified region and the determined distribution; and■ adding the identified region to the hot sources detection mask when it is judged that the identified region does not belong to the background; determining heat emitting regions as a result of intersecting (ITS) the hot sources detection mask determined for the first night-time satellite image with the hot sources detection mask determined for the second night-time satellite image.

2. The method of claim 1, further comprising determining a source location of a determined heat emitting region.

3. The method of claim 2, wherein determining a source location of a determined heat emitting region comprises determining a point of the determined heat emitting region having the highest cumulative brightness over the first and second satellite image.

4. The method of any one of claims 1 and 3, further comprising outputting a map of the region of interest having a visual indicator at each source location of a determined heat emitting region.

5. The method of any one of claims 1 to 4, wherein intersecting (ITS) the hot sources detection mask determined for the first night-time satellite image with the hot sources detection mask determined for the second night-time satellite image produces a final detection mask and further comprising a step of filtering (FILT) the final detection mask so that each heat emitting region in the filtered final detection mask comes from a different identified region in one of the hot sources detection masks.

6. The method of claim 5, wherein filtering the final detection mask comprises validating, among heat emitting regions resulting from the intersection of the hot sources detection mask determined for the first night-time satellite image with a single identified region in the hot sources detection mask determined for the second night-time satellite image, only the largest region.

7. The method of any one of claims 1 to 6, wherein identifying regions of contiguous bright pixels within the image implements a greedy region algorithm.

8. The method of any one of claims 1 to 7, wherein the determined distribution comprises a central intensity, the dark pixels being pixels having an intensity lowerthan the central intensity and the bright pixels being pixels having an intensity higher than the central intensity.

9. The method of any one of claims I to 8, wherein for each identified region, judging whether or not the identified region belongs to the background based on comparing the identified region to the determined distribution comprises: calculating a sum of bright pixel values of the bright pixels of the identified region; determining a probability that a sum of background pixel values is higher than the calculated sum of bright pixel values, wherein the sum of background pixel values is derived from the determined distribution considering as many background pixels as the number of bright pixels in the identified region; comparing the determined probability to a threshold; when the determined probability is lower to the threshold, judging that the identified region does not belong to the background.

10. The method of claim 9, wherein the threshold is set based on a size of the identified region and on a size of the image.

11. The method of any one of claims 1 to 10, wherein determining a distribution of a background model of the image comprises determining a background standard deviation.

12. The method of claim 11, wherein determining a background standard deviation comprises calculating an image biweight midvariance.

13. The method of any one of claims 1 to 12, wherein the first and second spectral bands are shortwave infrared bands.

14. The method of claim 13, wherein the first night-time satellite image and the second night-time satellite image are acquired by a same imager onboard a satellite.

15. A computer program product comprising instructions which, when the program is executed by a computer, cause the computer to carry out the method of any one of claims 1 to 14.

16. A system, comprising a processor configured to perform the method of any one of claims 1 to 14.

Images