Apparatus and Method for Cloud Masking

1. Field of the Invention

The presently disclosed invention relates to an apparatus and method for collecting satellite images of bodies of water such as the oceans and removing interfering clouds from the satellite images to improve the accuracy of the images. More particularly, the presently disclosed invention relates to the use of an apparatus and method to extract cloud edge information through local segmentation of the image and a procedure to discriminate between cloud free and cloud contaminated pixels. The cloud contaminated pixels are then removed from daytime and nighttime 1-km advanced very high resolution radiometer (AVHRR) satellite image.

1. Description of the Related Art

The quality of sea surface temperature (SST) measurements derived from spaceborn instruments depends largely on the accuracy of the cloud removal technique to eliminate cloud-contaminated pixels from the analysis. Clouds present a complicated problem. Generally, clouds are characterized by high albedos (i.e. the proportion of light or radiation reflected by a surface, especially of a planet or moon, for example, the whiteness of a reflection, The Concise Oxford Dictionary of Current English, Eighth Edition, R. E. Allen Ed., Clarendon Press, Oxford (1990)). In the daytime data collected, clouds have a greater whiteness or albedo than the surrounding non-cloud image. Clouds are also characterized by lower temperatures than the surrounding image. When clouds do not fill the instantaneous field of view (IFOV) of the instrument sensor, the cloud signatures become difficult to separate from the surrounding targets. Cirrus, low stratus and small cumulus are examples of clouds whose signatures can introduce large errors to the SST measurements.

Methods of cloud removal from satellite images vary from simple radiance thresholds to rather complicated procedures involving artificial intelligence techniques. Because a single technique does not appear to be sufficient in eliminating time and geographic dependency inherent within most techniques, the tendency is to combine a series of techniques. For instance, the operational cloud masking routine developed for the National Environmental Satellite Data and Information Service of the National Oceanic and Atmospheric Administration (NOAA NESDIS) uses a series of tests (i.e. gross cloud test, visible, thermal IR tests, and low stratus test) designed to eliminate different cloud signatures at several stages. The techniques implemented by the SeaSpace Corporation, San Diego, Calif. (SEA) remove clouds from AVHRR scenes based on a series of tests and thresholds similar to (but not the same as) those used by NOAA NESDIS. These latter procedures, which allow for further relaxation of the input parameter, produce less conservative (i.e. more cloud contaminated pixels are present in the satellite image) images than the NOAA NESDIS techniques. Simpson and Humphrey, infra methods rely on a combination of visible and infrared data, empirical albedo model as function of solar zenith angle, a radiance threshold, and sampling techniques across the image to insure that all Sun-pixel-satellite angles associated with that image are included in the calculations. See Simpson, J. J. and Humphrey, C., An automated cloud screening algorithm for daytime advanced very high resolution radiometer imagery, J. Geophys. Res., 95, 13, 459-13, 481 (1990), incorporated by reference herein in its entirety for all purposes. Although the methods may vary, standard cloud techniques rely on spectral rather than on texture signatures of the scene.

Ebert demonstrated that the introduction of pattern recognition techniques to AVHRR cloud removal methods produced cloud cover and cloud field structure information with great accuracy and less computation time than techniques that used spectral information exclusively. See Ebert, E., A pattern recognition technique for distinguishing surface and cloud types in the polar regions, J. Clin. Appl. Meteorol., 26, 1412-1427 (1987), incorporated by reference herein in its entirety for all purposes. Texture information alone extracted from a single visible channel was sufficient to detect, segment and classify different cloud types in polar regions. The gray level co-occurrence (GLC) matrix is one approach for extracting texture information of a scene based on the overall spatial relationship existing among its gray levels (i.e. intensity). See Haralick, R. M., Shanmugam, K., and Denstein, I., Textural features for image classification, IEEE Trans. Syst. Man Cybernt., 3, 610-621 (1973), incorporated by reference herein in its entirety for all purposes.

The Sea Surface Temperature Analysis and Composite (SSTAC) module is part of the Navy Tactical Environmental Support System (TESS) created to provide ships at sea with the capability to produce SST imagery from real-time 1-km-resolution AVHRR data. Currently, the module includes time efficient routines to identify and remove clouds from these data. Unfortunately, the daytime version of the cloud-masking technique is highly conservative. It overestimates cloud cover and removes a large number of cloud free pixels from each scene. The nighttime version of the algorithm, although less rigorous than the daytime version, produces noisy data sets and allows for a large number of cloud-contaminated pixels to enter the SST calculations. See Phegley, L. and Crosiar, C., The third phase of TESS, Bull. Am. Meteoroi. Soc., 72, 954-960 (1991), incorporated by reference herein in its entirety for all purposes. None of the techniques reviewed met the accuracy, time efficiency and geographical independence required for the production of accurate SST images.

Current methods of cloud identification and removal in satellite imagery of the oceans rely on spectral features of the scene exclusively. In some cases the methods are overly conservative and tend to remove large numbers of cloud free pixels. In other cases, the methods perform satisfactorily for very cold clouds with high albedo values but underestimate the low-level clouds whose temperature and albedos are similar to that of the ocean. Cirrus, low stratus and small cumulus are examples of such clouds. In addition, the current methods of cloud removal do not perform consistently over a wide range of latitudes and environmental conditions. In fact, they do not perform consistently within the same image. Accurate retrievals of oceanographic parameters such as the SST and ocean color depend largely on an accurate method of cloud removal. For a method of cloud detection to be useful it needs to be: (1) automated, (2) time efficient, (3) accurate, and (4) independent of geographic location and environmental conditions. Current methods of cloud detection do not meet these requirements.

It is therefore an object of the present invention to provide an apparatus and method for the retrieval of SST images from bodies of water such as the oceans and for the accurate removal of clouds from the retrieved SST images.

It is another object of the present invention to detect and remove cloud contaminated pixels from SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and remove cloud contaminated pixels from daytime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and remove cloud contaminated pixels from nighttime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and remove cirrus, low stratus and small cumulus clouds whose signatures can introduce large errors into the SST measurements from SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and remove cirrus, low stratus and small cumulus clouds whose signatures can introduce large errors into the SST measurements from daytime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and remove cirrus, low stratus and small cumulus clouds whose signatures can introduce large errors into the SST measurements from nighttime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to provide an apparatus and method for the formation of a cloud mask from SST images of bodies of water such as the oceans.

It is another object of the present invention to detect and form a cloud mask of cloud contaminated pixels from SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and form a cloud mask of cloud contaminated pixels from daytime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and form a cloud mask of cloud contaminated pixels from nighttime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and form a cloud mask of clouds such as cirrus, low stratus and small cumulus clouds whose signatures can introduce large errors into the SST measurements from SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and form a cloud mask of clouds such as cirrus, low stratus and small cumulus clouds whose signatures can introduce large errors into the SST measurements from daytime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

It is another object of the present invention to detect and form a cloud mask of clouds such as cirrus, low stratus and small cumulus clouds whose signatures can introduce large errors into the SST measurements from nighttime SST images automatically, time efficiently, accurately, and independent of geographic or environmental conditions.

These and other objects are accomplished by obtaining spectral information of the oceans using visible and infrared channels. Texture information is extracted from the visible and infrared (i.e. thermal channels) channels using the cluster shade texture measure. The cluster shade texture measure is utilized to determine edge pixels. The edge pixel determination is repeated until all edge identified pixels are joined to other edge identified pixels forming a boundary around contiguous nonedge pixel polygons having an outer shape defined by the edge pixel boundary. Having determined the shape of the contiguous nonedge pixel polygons and having identified the edge pixel boundary, for each polygon, the determination of whether the polygon is cloud free or cloud contaminated is made on the basis of the assumption that clouds are colder and brighter than the surrounding ocean. Polygons determined to be clouds are removed from the initial SST image to produce a nearly cloud free SST image. The step of removing polygons identified as clouds may be carried out at a later time. One may practice the presently disclosed invention by merely forming the cloud mask, according to the presently disclosed invention, for subsequent use by another party.

The file of this patent contains at least one drawing executed in color. Copies of this patent with color drawing(s) will be provided by the Patent and Trademark Office upon request and payment of the necessary fee.

FIG. 1 (plate 1) is a daytime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 15, 1987. Clouds, shown as black pixels, are removed from the image using constant radiance threshold (CRT) methods computed from visible and infrared data. For comparison this photograph is representative of the daytime photograph having used CRT methods.

FIG. 2 (plate 2) is a daytime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 15, 1987. Clouds, shown as black pixels, are removed from the image using the SeaSpace routines. For comparison this photograph is representative of the daytime photograph having used the SeaSpace routines.

FIG. 3 (plate 3) is a daytime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 15, 1987. Clouds, shown as black pixels, are removed from the image using the Simpson and Humphrey routines. For comparison this photograph is representative of the daytime photograph having used the Simpson and Humphrey routines.

FIG. 4 (plate 4) is a daytime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 15, 1987. Clouds, shown as black pixels, are removed from the image using the NOAA NESDIS techniques. For comparison this photograph is representative of the daytime photograph having used the NOAA NESDIS routines.

FIG. 5 (plate 5) is a daytime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 15, 1987. Clouds, shown as black pixels, are removed from the image using the apparatus and method of the currently disclosed invention (the NRL-SSC cloud masking technique). For comparison this photograph is representative of the daytime photograph having used the NRL-SSC apparatus and routines of the presently disclosed invention.

FIG. 6 (plate 6) is a nighttime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 10, 1989. Clouds, shown as black pixels, are removed from the image using constant radiance threshold (CRT) methods computed from visible and infrared data. For comparison this photograph is representative of the nighttime photograph having used CRT methods.

FIG. 7 is a nighttime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 10, 1989. Clouds, shown as black pixels, are removed from the image using the NOAA NESDIS routines. For comparison this photograph is representative of the nighttime photograph having used the NOAA NESDIS nighttime routines.

FIG. 8 is a nighttime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 10, 1989. Clouds, shown as black pixels, are removed from the image using the apparatus and method of the present invention (the NRL-SSC routines). For comparison this photograph is representative of the nighttime photograph having used the NRL-SSC apparatus and methods of the presently disclosed invention.

FIG. 9 is a daytime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 15, 1987. Clouds, shown as black pixels, are removed from the image using a mask produced through local thresholds and general human interpretation. For comparison this photograph is representative of the daytime photograph using general human interpretation.

FIG. 10 is a nighttime AVHRR SST scene of the waters of the east coast of the United States acquired on Apr. 10, 1989. Clouds, shown as black pixels, are removed from the image using a mask produced through local thresholds and general human interpretation. The apparatus and method of the currently disclosed invention (the NRL-SSC cloud masking technique). For comparison this photograph is representative of the nighttime photograph using general human interpretation.

Color photographs of FIGS. 1-10 are incorporated by reference herein in their entirety for all purposes. To view the color photographs see the following article submitted with the attached information disclosure statement accompanying this application. See Gallegos, S. C., Hawkins, J. D., and Cheng, C. F., A New Automated Method of Cloud Masking for Advanced Very High Resolution Radiometer Full-Resolution Data Over the Ocean, Journal of Geophysical Research, Vol. 98, NO. C5 (1993), incorporated by reference herein in its entirety for all purposes.

FIG. 11 depicts a flow chart of the various method steps for forming the cloud mask according to the presently disclosed invention.

FIG. 12 is a schematic representation of an AVHRR scene with the intensity of the pixel indicated by a number from 0 to (L-1) at the spatial location of the pixel. The pixels are separated into local windows of size N.sub.0 .times.M.sub.0, N.sub.0 .times.M.sub.1, N.sub.0 .times.M.sub.2, N.sub.0 .times.M.sub.3, N.sub.1 .times.M.sub.0, N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, . . . as indicated in FIG. 12 wherein the local size of each window is 6 pixels by 6 pixels, as an example. Further, adjacent local windows N.times.M are separated by a distance of Dx=6, or -6 and Dy=0, 6, or -6, as indicated in exemplary FIG. 12.

In FIG. 13, the cluster shade values calculated for each of the local windows or neighborhood of pixels identified as N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 .times.M.sub.1, N.sub.3 .times.M.sub.2, and N.sub.3 .times.M.sub.3 in FIG. 12, are as indicated, respectively.

The GLC matrices for local windows N.sub.0 .times.M.sub.0, N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 33 M.sub.1, N.sub.3 .times.M.sub.2, and N.sub.3 .times.M.sub.3 are calculated and shown in FIGS. 14, 17, 19, 15, 21, 23, 25, 27, 29 and 31, respectively.

For the exemplary local windows N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 .times.M.sub.1, N.sub.3 .times.M.sub.2, and N.sub.3 .times.M.sub.3 of FIG. 12, the calculations of the cluster shade values, S (.DELTA.x, .DELTA.y), are given in FIGS. 18, 20, 16, 22, 24, 26, 28, 30 and 32, respectively.

The necessary calculations of u.sub.i and u.sub.j for the exemplary local windows N.sub.1 .times.M.sub.1, N.sub.1.times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 .times.M.sub.1, N.sub.3 .times.M.sub.2, and N.sub.3 .times.M.sub.3 of FIG. 12, are given in FIGS. 17, 19, 15, 21, 23, 25, 27, 29 and 31, respectively.

An exemplary cluster shade texture measure image derived from FIG. 12 is depicted as previously indicated in FIG. 13. FIG. 13 has been incorporated into a larger exemplary cluster shade texture measure image in FIG. 33 wherein the respective cluster shade texture values are indicated.

Exemplary FIG. 34 depicts a binary file image corresponding to FIG. 33 wherein the 6.times.6 edge pixel windows are indicated with shading and with an intensity value of 0 and the nonedge pixel windows are indicated without shading and an intensity value of 1.

The following detailed description of the invention is provided to aid those skilled in the art in practicing the present invention. However, the following detailed description of the invention should not be construed to unduly limit the present invention. Variations and modifications in the embodiments discussed may be made by those of ordinary skill in the art without departing from the scope of the present inventive discovery.

Several terms and symbols are used in this application which have meaning as described below:

L=the number of intensity levels. In some of the examples that follow L=3. In obtaining FIGS. 5 and 8, L was set to equal 1024. L may vary from 0 to infinity. Typically L is between 0-100,000, inclusive. More typically, L is between 0-10,000, inclusive. Most typically L is between 0-5,000, inclusive. Preferably, L is between 0-2,500, inclusive. More preferably, L is from 800-1,500, inclusive. Most preferably, L is from 0-1024, inclusive.

N=the number of pixels (i.e. columns of pixels) in the local neighborhood. N is an integer between 2 to infinity. In some of the examples that follow N=6. Typically, N is an integer between 2-2,000, inclusive. More typically, N is an integer between 2-1,000, inclusive. Most typically, N is an integer between 2-500, inclusive. Preferably, N is an integer between 2-100, inclusive. More preferably, N is an integer between 2-6, inclusive. Most preferably, N is an integer between 3-6, inclusive.

M=the number of lines (i.e. rows of pixels) in the local neighborhood. M is an integer between 2 to infinity. In some of the examples that follow M=6. Typically, M is an integer between 2-2, 000, inclusive. More typically, M is an integer between 2-1,000, inclusive. Most typically, M is an integer between 2-500, inclusive. Preferably, M is an integer between 2-100, inclusive. More preferably, M is an integer between 2-6, inclusive. Most preferably, M is an integer between 3-6, inclusive. Note that both M and N cannot equal 2 at the same time.

Dx=the vertical distance between local neighborhoods of size N.times.M. Typically, overlapping local neighborhoods of size 1000.times.1000 wherein Dx=300 pixels and Dy=300 pixels may be used. More typically, overlapping local neighborhoods of size 500.times.500 wherein Dx=150 pixels and Dy=150 pixels may be used. Most typically, overlapping local neighborhoods of size 200.times.200 wherein Dx=60 pixels and Dy=60 pixels may be used. Preferably, overlapping local neighborhoods of size 100.times.100 wherein Dx=30 pixels and Dy=30 pixels may be used. More preferably, overlapping local neighborhoods of size 64.times.64 wherein Dx=1 pixel and Dy=1 pixels may be used. Most preferably, overlapping local neighborhoods of size 3.times.3 wherein Dx=1 pixel and Dy=1 pixel may be used.

.DELTA.x=the horizontal component of the displacement vector within each local neighborhood of size N.times.M used to compute cluster shade. Typically, .DELTA.x is between 0-100, inclusive. More typically, .DELTA.x is between 0-50, inclusive. Most typically, .DELTA.x is between 0-25, inclusive. Preferably, .DELTA.x is between 1-10, inclusive. More preferably, .DELTA.x is between 1-5, inclusive. Most preferably, .DELTA.x is 1.

.DELTA.y=the vertical component of the displacement vector within each local neighborhood of size N.times.M used to compute cluster shade. Typically, .DELTA.y is between 0-100, inclusive. More typically, .DELTA.y is between 0-50, inclusive. Most typically, .DELTA.y is between 0-25, inclusive. Preferably, .DELTA.y is between 0-10, inclusive. More preferably, .DELTA.y is between 0-5, inclusive. Most preferably, .DELTA.y is 1.

One kilometer AVHRR scenes are obtained using a conventional satellite data capturing techniques well known in the industry. See Kidwell, K. N., Polar Orbiter data user's guide, NOAA-NCDC, SDIS, U.S. Dept. of Commerce, Washington, D.C. (1983), incorporated by reference herein in its entirety for all purposes. See Holyer, R. J. and Peckinpaugh, S. H., Edge detection applied to satellite imagery of the oceans, IEEE Trans. Geosci. Remote Sens., 27, 46-56 (1989), incorporated by reference herein in its entirety for all purposes. See Parikh, J., A comparative study of cloud classification techniques, Remote Sens. Environ., 6, 67-81 (1977), incorporated by reference herein in its entirety for all purposes. See Welch, R. M., Sengupta, S. K., and Chen, D. W. Cloud field classification based upon high spatial resolution textural features. 1. Gray level co-occurrence matrix approach, J. Geophys. Res., 93, 12, 663-12, 681 (1988), incorporated by reference herein in its entirety for all purposes. See Welch, R. M., Sengupta, S. K., and Navar, M. S., The effect of spatial resolution upon texture-based cloud field classifications, J. Geophys. Res., 94, 14, 767-14, 781 (1989), incorporated by reference herein in its entirety for all purposes. See Welch, R. M., Sengupta, S. K., Goroch, K., Rabindra, P., Rangaraj, N. and Navar, M. S., Polar cloud surface classification using AVHRR imagery: An intercomparison of methods, J. Appl. Meteorol., 31, 405-419 (1992), incorporated by reference herein in its entirety for all purposes. A flow chart depicting the various method steps are presented in FIG. 11. Spectral information is derived from one visible and three infrared (IR) channels. These channels are channel 1 (585-685 nm), channel 3 (3,575-3,983 nm), channel 4 (10,362-11,299 nm) and channel 5 (11,450-12,424 nm). The first channel is utilized for resolving clouds during the day, and the last three channels are used for detecting clouds at night. The visible channels of the AVHRR instrument are customarily used to identify clouds in AVHRR scenes at various remote sensing facilities. Channel 2 provides high contrast between clouds and water. However, channel 1 was chosen because it appears to have better texture information than channel 2. At night, in the absence of reflective data, differences in thermal channel responses ((channel 3-channel 4) and (channel 4-channel 5)) are utilized. Spectral differences at these wavelengths ((channel 3-channel 4) and (channel 4-channel 5)) measure the opacity of the atmosphere.

Texture information is extracted from the visible (channel 1) and thermal channels ((channel 3-channel 4) and (channel 4-channel 5)) using the cluster shade texture measure, a local segmentation scheme based on the GLC matrix.

Initial steps in the resolution of cloud edges in a scene requires the computation of the GLC matrix. The (i, j)th element of the matrix is the relative frequency of occurrence of gray level i and gray level j when the ith gray level is separated by a distance or displacement vector (.DELTA.x, .DELTA.y) from the jth gray level within a window or local neighborhood of size M.times.N. Given a range of intensity levels or gray levels from 0 to (L-1), the GLC matrix is calculated from: ##EQU1## if f (m,n)=i and f(m+.DELTA.x, n+.DELTA.y)=j and

otherwise, f(m,n) is the intensity level of the pixel at sample (m, n), the ith gray level, and f(m+.DELTA.x, n+.DELTA.y) is the intensity of the pixel at sample (m+.DELTA.x, n+.DELTA.y), the jth gray level. FIG. 12 is a schematic representation of an AVHRR scene with the intensity of the pixel indicated by a number from 0 to (L-1) at the spatial location of the pixel. The pixels arc separated into local windows of size N.sub.0 .times.M.sub.0, N.sub.0 .times.M.sub.1, N.sub.0 .times.M.sub.2, N.sub.0 .times.M.sub.3, N.sub.1 .times.M.sub.0, N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, . . . as indicated in FIG. 12 wherein the local size of each window is 6 pixels by 6 pixels, as an example. Further, adjacent local windows N.times.M are separated by a distance of Dx=6, or -6 and Dy=0, 6, or -6, as indicated in exemplary FIG. 12. As an example, the GLC matrices for local windows N.sub.0 .times.M.sub.0, N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 .times.M.sub.1, N.sub.3 .times.M.sub.2, and N.sub.3 .times.M.sub.3 are calculated and shown in FIGS. 14, 17, 19, 15, 21, 23, 25, 27, 29 and 31, respectively.

FIGS. 5 and 8 were obtained wherein the local windows were of size N.times.M=3.times.3 and wherein local windows were separated by a distance Dx=1, or -1 and Dy=1 or -1. Thus, while the exemplary local windows N.sub.0 .times.M.sub.0, N.sub.0 .times.M.sub.1, N.sub.0 .times.M.sub.2, N.sub.0 .times.M.sub.3, N.sub.1 .times.M.sub.0, N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, . . . in FIG. 12 do not overlap, local windows used to obtain FIGS. 5 and 8 did overlap. The small window size and displacement vectors used in obtaining FIGS. 5 and 8 were geared towards detecting changes in intensity occurring over small distances in the AVHRR scene. The small window size and small displacement vectors are preferred for locating small clouds and small cloud fringes. As is apparent from FIGS. 14, 17, 19, 15, 21, 23, 25, 27, 29 and 31, the GLC matrix computations result in a second-order probability matrix of size L.times.L from which cluster shade, a measure of the local variation in intensity between pixels, is calculated for each local neighborhood N.times.M according to the formula: ##EQU2## where u.sub.i and u.sub.j are estimates of mean intensity calculated from the weighted summations of rows and columns of the GLC matrix, according to the formulas: ##EQU3## For the exemplary local windows N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 .times.M.sub.1, N.sub.3 .times.M.sub.2, and N.sub.3 .times.M.sub.3 of FIG. 12, the calculations of the cluster shade values, S (.DELTA.x, .DELTA.y), are given in FIGS. 18, 20, 16, 22, 24, 26, 28, 30 and 32, respectively. The necessary calculations of u.sub.i and u.sub.i for the exemplary local windows N.sub.1.times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 .times.M.sub.1, N.sub.3 .times.M.sub.2, and N.sub.3 .times.M.sub.3 of FIG. 12, are given in FIGS. 17, 19, 15, 21, 23, 25, 27, 29 and 31, respectively.

From the computation of the cluster shade values, a new image is produced in which the center pixel of each local neighborhood (i.e. note that a local neighborhood N.times.M wherein N=M=an odd integer, a center pixel is readily identified whereas if N=M=an even integer, then one must choose one of four pixels to consistently to be the center pixel) is replaced with a cluster shade value calculated for that neighborhood. Thus, in exemplary FIG. 13, the cluster shade values calculated for each of the local windows or neighborhood of pixels identified as N.sub.1 .times.M.sub.1, N.sub.1 .times.M.sub.2, N.sub.1 .times.M.sub.3, N.sub.2 .times.M.sub.1, N.sub.2 .times.M.sub.2, N.sub.2 .times.M.sub.3, N.sub.3 .times.M.sub.1, N.sub.3 33 M.sub.2, and N.sub.3 .times.M.sub.3 in FIG. 12, are as indicated, respectively. The calculated cluster shade values are utilized to generate edges using the zero-crossing test.

With an exemplary local neighborhood size of 3.times.3 and Dx=1 and Dy=1, each pixel is the center of a given local overlapping neighborhood. After calculation of the cluster shade values, as indicated by EQUATION 2 above, for the entire AVHRR scene, each pixel now has associated with it a calculated cluster shade value, S(.DELTA.x, .DELTA.y). Due to the overlapping of local neighborhoods of size 3.times.3, and Dx=1 and Dy=1 between local neighborhoods, each pixel of the AVHRR scene has a computed cluster shade value associated with a given pixel. The entire set of cluster shade values forms a cluster shade texture image.

Edges are generated by a zero crossing test. Edges are generated at the site of the neighboring pixels whose cluster shade values are opposite in sign (zero crossing) and where both have cluster shade intensities greater than or equal to a chosen threshold. To avoid producing edges for noncloud targets, thresholds are introduced. The daytime version uses an exemplary cluster shade threshold of 3, and the nighttime version uses a threshold of 100.

The zero crossing test is implemented as follows. The user selects an initial threshold to be used for this test. Each exemplary 3.times.3 overlapping neighborhood of the computed cluster shade texture image is tested as follows. If the center pixel of the 3.times.3 overlapping neighborhood has an absolute cluster shade value.gtoreq.the initial threshold value, then the 8 pixel neighbors of the center pixel (i.e. in a 3.times.3 neighborhood) are tested as follows. If any of the 8 neighbors have an absolute cluster shade value.gtoreq.the initial threshold and have an opposite sign from the center pixel cluster shade value, then that center pixel is considered an edge pixel. Although most of the edges (i.e. pixels identified as edges using the zero crossing test described above) are in place after conducting the zero crossing test for each pixel one time, some of the pixels are determined to be edges forming isolated strings (i.e. the pixels defined as edges that do not form closed polygons). Since clouds occupy a finite space, individual clouds are defined by an outer edge boundary taking the shape of a polygon. Since some of the polygons defined by the edge pixels represent clouds, it is necessary to enable the isolated strings (i.e. a series of adjacent edge pixels, the ends of which are adjacent to nonedge pixels) to join with other strings to form polygons, by extending the original edge strings using a procedure similar to the one utilized in the first zero crossing test.

The extended zero crossing test is optional. This extension computation is used to extend edges detected in the first zero crossing test. The parameters involved are a second or minimum threshold and maximum allowed passes through the cluster shade texture measure image.

Thus, an extended zero erossing test (extension computation) is performed. New thresholds are set. An exemplary daytime threshold is 1 and an exemplary nighttime threshold is 40. For the edge extension computation to be successful, the exemplary 3.times.3 local neighborhoods of nonedge pixels are examined. For an exemplary local 3.times.3 window or neighborhood to enter the edge extension computations, the cluster shade value of the center pixel must be greater than the new threshold and one of the exemplary 9 pixel members of the exemplary 3.times.3 local window must have been defined as an edge in the previous zero crossing test. New edge pixels are generated at the site of those pixels whose cluster shade value is greater than and opposite in sign to the new threshold. While the extended zero crossing test may be repeated 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, or more times, the extended zero crossing test is exemplarily repeated 17 times, after which most of the connections have been established. Further runs of the extended zero crossing test (greater than 17) only increases processing time and yields very few additional connections.

After the zero crossing test and the extended zero crossing test are performed on the cluster shade texture measure image, edge pixels are identified as described above. Edge pixels are assigned a zero value, and nonedge pixels are set to 1. Areas of polygons in the context of overlapping 3.times.3 local windows are clusters of contiguous nonedge pixels (intensity=1) which are surrounded by one-pixel-wide boundaries of edge pixels (intensity=0). Thus, a binary file is generated wherein edges are represented by 0s and the polygons within the edges are represented by 1s.

An exemplary cluster shade texture measure image derived from FIG. 12 is depicted as previously indicated in FIG. 13. FIG. 13 has been incorporated into a larger exemplary cluster shade texture measure image in FIG. 33 wherein the respective cluster shade texture values are indicated. In exemplary FIG. 33, 6.times.6 nonoverlapping neighborhoods are subjected to the zero crossing test once. An exemplary threshold of 0.100 was used. After one pass of the zero crossing test, those 6.times.6 windows identified as edges are indicated by an X in the local window and a darker shade. Note that because 6.times.6 nonoverlapping windows were tested, the nonedge pixels are identified by 6 or more pixel wide boundaries of edge pixels. Exemplary FIG. 34 depicts a binary file image corresponding to FIG. 33 wherein the 6.times.6 edge pixel windows are indicated with shading and an intensity value of 0 and the nonedge pixel windows are indicated without shading and an intensity value of 1.

Areas in the AVHRR scene that fail to produce edges occur in homogeneous regions associated with cloud free oceanic areas or with the tops of very cold clouds. Heterogeneous regions, which produce large numbers of edges are found at the interphase of the clouds and the ocean, and also between cloud layers. Because the edge information produced in the binary file image was not sufficient to permit discrimination of cloud masses (i.e. tops of very cold clouds, the interphase of clouds and the ocean, and between cloud layers), an additional routine designed to separate cloudy from clear areas in the AVHRR scene was implemented. This additional routine automatically searches, identifies and labels areas or polygons generated in the binary files as clouds or as ocean. Areas identified as clouds become part of a cloud mask.

Each line of data from the binary files is scanned to search for strings of contiguous nonedge pixels (i.e. strings of 1s). Note that nonedge strings are being scanned. Upon encountering the beginning of a nonedge string, the algorithm first assigns the nonedge string a unique exemplary 16-bit identification number and then sets out to find links between that string and neighboring strings. The linkage information is stored in an exemplary 16-bit pointer array. When two nonedge strings connect, the resulting nonedge string takes on the lower of the two identification numbers. The process of searching for nonedge strings and linking them continues until all possible connections have been established and the nonedge strings are clustered into labeled polygons.

A given image can contain from a few polygons to thousands of polygons. To reduce the number of polygons, and thereby reduce processing time for the scene, the smaller polygons are eliminated from further calculations by assimilating them into the cloud mask (i.e. similar to a template that displays where all the clouds are located within the AVHRR scene). The criterion used for this elimination is that if the total number of pixels within the polygon is smaller than the total number of edge pixels enclosing it, the polygon is considered to be cloud noise. Thus, for example, 6.times.6 pixel window N.sub.2 .times.M.sub.4 of FIG. 34 would be assimilated into the cloud mask.

The larger remaining polygons are separated on the basis of the assumption that clouds are colder and brighter than the surrounding ocean. Mean calibrated albedos (channel 1 daytime pixel intensity data is used) or brightness temperatures (channel 4 nighttime pixel intensity--temperature--data is used) are calculated for each polygon and its boundary. The albedos are calculated from channel 1 data and the brightness temperatures are calculated from channel 4 data.

A cloud free polygon should have a lower mean albedo than that of the edge boundary or a higher mean brightness temperature than that of the edge boundary. This assumption holds because clouds are more white (higher albedo) and colder (lower temperature) than the oceans below. Thus, the opposite is true for a cloud filled polygon. A polygon that has been declared to be cloud filled by these procedures is assimilated into the cloud mask and all of its pixels are set to 0. Pixels in a cloud free polygon are set to 1.

By applying the aforementioned steps, a procedure is designed to produce a refined binary image identifying regions of clouds and regions of cloud free areas. The regions of clouds identified form a cloud mask. The cloud mask is used to remove SST data corresponding to the spatial location within an AVHRR image that is occupied by a cloud as determined by the cloud mask.

Note that according to the flowchart of FIG. 11, channel 1 albedo data is used to form a GLC matrix (not shown). The GLC matrix is used to form a cluster shade texture measure image (clustershade). The cluster shade texture measure image is subjected to zero crossing tests to generate edge pixels (edge generation). Nonedge pixels are located and linked to form polygons (find areas and link areas). The polygons are subjected to mean albedo comparison with the mean albedo of the edge boundary. If a polygon has a higher mean albedo than the mean albedo of the edge boundary, the polygon is determined to be a cloud (statistics). If the polygon is determined to be a cloud, the polygon is filled with Os (fill areas). After all the polygons have been subjected to statistics, all cloud polygons are combined to form a cloud mask. When SST AVHRR data spatially corresponds to the same location as a cloud (i.e. occupying the same space as indicated by the cloud mask), the SST AVHRR data is discarded and is depicted as the black cloud region as indicated in FIG. 5 (apply mask).

Note that according to the flowchart of FIG. 11, opacity data obtained at nighttime by the subtraction of (channel 3-channel 4 data) and (channel 4-channel 5 data) is utilized individually to form a GLC matrix (not shown) for each of (channel 3-channel 4 data) and (channel 4-channel 5 data). A cluster shade texture measure image is formed for each of (channel 3-channel 4 data) and (channel 4-channel 5 data) (i.e. clustershade). The cluster shade texture measure image for each of (channel 3-channel 4 data) and (channel 4-channel 5 data) is subjected to zero crossing tests to generate edge pixels (edge generation). The two edge generated matrices from each of (channel 3-channel 4 data) and (channel 4-channel 5 data) are combined to form a single matrix defining generated edge pixels (merge masks). Thereafter, the nonedge pixels are located and linked to form polygons (find areas and link areas). The polygons of nonedge pixels are subjected to mean temperature comparison with the mean temperature of the edge boundary. If a polygon has a lower mean temperature than the mean temperature of the edge boundary, the polygon is determined to be a cloud (statistics). If the polygon is determined to be a cloud the polygon is filled with 0s (fill areas). After all the polygons have been subjected to statistics, all cloud polygons are combined to form a cloud mask. When SST AVHRR data spatially corresponds to the same location as a cloud (i.e. occupying the same space as indicated by the cloud mask), the SST AVHRR data is discarded and is depicted as the black cloud region as indicated in FIG. 8 (apply mask).

Having described the invention, the following examples are given to illustrate specific applications of the invention, including the best mode now known to perform the invention. These specific examples are not intended to limit the scope of the invention described in this application.

The NRL-SSC algorithm was designed to be an operational procedure for shipboard use, and as such, its development was governed by three requirements: accuracy, geographical independence, and time efficiency. The extent to which this algorithm met the first two requirements is demonstrated in the examples below. The time efficiency factor was assessed via CPU times for each run. Two systems were used for this purpose: A VAX/8800 and a CONCURRENT/5450. The first system was chosen because it is widely used and can give the reader a frame of reference, and the second system is an early version of the TESS computer installed on board Navy ships. In both systems, the CPU time was less than 55 seconds for the daytime algorithm and less than 2 minutes for the nighttime algorithm.

The AVHRR image collections were performed following recommendations by the National Environmental Satellite Data and Information Service of the National Oceanic and Atmospheric Administration (NOAA/NESDIS). See Kidwell, K. N., (1983), supra.

The actual code of the NRL-SSC algorithm in FORTRAN used is given on the following pages.

Evaluation of the capabilities of this technique to detect clouds from AVHRR data was accomplished through comparisons of the NRL-SSC algorithm with cloud mask results of (1) other experimental and operational cloud removal techniques and (2) human interpretation of images of the cloud cover at the time of the satellite pass. These evaluations were performed with 1-km full-resolution AVHRR data.

Both daytime and nighttime evaluations were accomplished in the same area of the central and north Atlantic, from 24.18.degree. N to 45.00.degree. N and from 49.42.degree. W to 75.00.degree. W. This area was chosen because of availability of data and for its wide range of atmospheric and oceanographic gradients. Both the daytime and nighttime sets include cold water regimes in the higher latitudes and warm water in the lower latitudes. There is also a variety of cloud types, ranging from cold cloud masses to the thin warmer clouds of the lower latitudes.

The daytime evaluation uses data from 15 Apr. 1987 and five cloud masks: CRT, SEA, S&H, NOAA and NRL-SSC. Scenes containing cloud masks from these five techniques appear in FIGS. 1-5, respectively.

These scenes as well as the nighttime scenes (FIGS. 6-10) were processed with a color table in which colder waters are depicted by beige, yellow and red hues while warmer waters are shown as greens and blues. This color assignment is opposite to normal convention and was chosen for no other reason than to show contrast between cloud edges and the surrounding ocean. Clouds identified by the individual cloud removal techniques are shown as black. Pixels that are cloud contaminated but that failed to be included in a given cloud mask appear as scattered yellow and red regions around the black cloud masses. Cloud contaminated pixels, as defined here, are pixels that exhibit anomalous radiances in reference to the surrounding pixels. These occur in the fringes of the large cloud masses and in the form of isolated patches over the more homogeneous ocean.

The CRT mask (FIG. 1) is a cloud mask produced by constant radiance thresholds (CRT). In this procedure clouds are removed from the scene by choosing arbitrary threshold values in the visible and infrared data. One drawback of this method is that it changes from image to image and depending on the area chosen, the thresholds may or may not apply for the entire image.

The SEA mask (FIG. 2) is a mask generated with the SeaSpace operational techniques. These routines use most but not all of the NOAA cloud removal routines and allow for further relaxation of the input parameters to the cloud and uniformity tests. Because the procedures are less restrictive, they yield a larger number of "cloud-free" pixels in the scenes. This mask was processed by W. Kastak of the Center for Space Research of the University of Texas at Austin using software commercially available from the SeaSpace Corporation, San Diego, Calif.

The S&H cloud mask (FIG. 3) is a cloud mask processed with the Simpson Humphrey daytime cloud masking algorithm. This mask was provided by J. Simpson of the Scripps Institution of Oceanography, La Jolla, Calif.

The NOAA cloud mask (FIG. 4) includes all the techniques used by NOAA NESDIS to remove clouds operationally from AVHRR data sets. This mask was provided by D. May of NRL-SSC.

The NRL-SSC cloud mask (FIG. 5) is the subject of the presently disclosed invention.

For the daytime evaluations, data from channels 1, 4 and 5 were processed. Channel 1 was utilized exclusively for cloud-masking purposes. Channels 4 and 5 were used to calculate multichannel sea surface temperature (MCSST). Channel 4 and MCSST were chosen to evaluate the performance of the various cloud masking procedures. Clouds were removed from the data with the various techniques under consideration, ensuring that for each comparison, the same algorithm was used to simultaneously remove clouds from channel 4 and the MCSST data set. The remaining cloud-free pixels were used in a regression analysis of channel 4 versus MCSST. Five regressions were carried out, one for each cloud mask: CRT, NOAA, SEA, S&H and NRL-SSC.

The best cloud mask will produce the least variability between channel 4 and MCSST data sets while yielding the largest number of cloud-free pixels in the AVHRR scene. The results of the daytime evaluations are presented in Table Ia below:

The nighttime evaluations used the same procedures as the daytime evaluations. These comparisons include nighttime data from Apr. 10, 1989 and three cloud masks: CRT, NOAA, and NRL-SSC. Scenes containing cloud masks from these three techniques appear in FIGS. 6-8, respectively.

For the nighttime evaluations, AVHRR image data from three thermal channels of the AVHRR instrument were processed. Channels 3, 4 and 5 were used for the MCSST computations and for the cloud masking procedures. The results of the nighttime evaluations are presented in Table Ib below:

Each of the cloud masks used in Table la were compared with cloud masks produced by an experienced operator using a combination of visual interpretation of the cloud cover, manual removal of cloud contaminated pixels, and local thresholds of reflective (daytime) and infrared (nighttime) data. Each cloud mask was stored in a binary data file where cloud covered pixels were assigned values of 0 and cloud free pixels were assigned values of 1. Areas of agreement between a cloud mask and the manual cloud mask used as a "control mask" were those where both cloud masks had values of 0 and where both had values of 1. Conversely, areas of disagreement occurred where one mask had values of 1 while the other exhibited values of 0. Percentage agreement and disagreement were calculated for all masks. The results are presented in Table 2 below. The control masks for the daytime and nighttime data sets appear in FIGS. 9 and 10, respectively.

The NRL-SSC cloud mask was tested over a wide range of geographical, atmospheric and seasonal conditions. A total of 138 1-km AVHRR image data sets were used to evaluate global performance of the NRL-SSC mask. Four AVHRR image data sets from four different areas of the world's oceans were randomly selected. The areas were the Gulf of Mexico, the Persian Gulf, the North Atlantic, and the central Pacific. Two of the sets were acquired at night and two were acquired during the day. For each set, a "control mask" was produced by experienced technical personnel. As before, the cloud masks generated with NRL-SSC technique were compared with four control masks. The results of these evaluations are presented in Table 3 below: