Geospatial Data Processing Characteristics for Environmental Monitoring Tasks

As standard, processing is reduced to localizing various natural and man-made objects on the underlying surface and evaluating their properties according to the satellite imagery parameters. The classic way to obtain information about environmental objects from satellite imagery when creating digital terrain models to form a coherent picture of the studied environmental situation on the ground is performed at image decoding stage. As such, the time taken to decode using traditional technologies ranges from 20 to 40 percent of the total. This is due to the fact that the most complex and time-consuming stages of decryption are performed using partially automated or even non-automated visual-based ways. Considerable attention is given to the development of methods and means of automated decryption and analysis of data obtained using remote sensing of the Earth [1]. In particular, attention is paid to the creation of software that implements image processing methods for classifying and encoding digital information about the terrain, the formation of structures for terrain data storage and delivery to consumers as well as pattern recognition based on remote aerospace research data. The results of similar studies are described in sufficient detail in the works of V.I. Lyalko, M.A. Popov [2], V.V. Kozoderov [3], Yu.I. Gurevich [4], R.E. Pashchenko [5], D. Lu, O. Weng [6] and others. However, an analysis of the considered works showed that there isn’t enough research pertaining to the tasks related to the creation of new methods for joint processing of multi-temporal image data obtained in different spectral ranges, as well as tasks related to improving the accuracy of satellite imagery materials geospatial referencing while automating the processing of large information amounts to increase the efficiency of decision-making on the current condition of monitored objects.

Currently, due to the constant increase in the speed and volume of transmitted remote sensing data, the tasks related to improving the efficiency of automated information processing due to new approaches to the processing and storage of significant amounts of geospatial data are relevant. It is known that at the moment, the speed of reception for satellite information reaches 800 Mbit/s (WorldView – 1, 2, 3, 4) [7] – 2000 Mbit/s (Ofeq-5) [8], and the information capacity recorded in one file communication session reaches 10 GB [7, 8]. In order for the results of the analysis of the obtained data to be effective and timely, complex processing of such a volume of information should be performed for at least 30 minutes after receiving it. In this case, the decisions made about the current state of the monitored objects and forecasts of the dynamics of its possible changes will be relevant.

An analysis of known methods showed that computer analysis of satellite images, as a rule, is based on the use of formalized natural features of the image or artificial features, that were designed on the basis of formalized natural features. In addition, it is shown that the most difficult to formalize structural feature is texture. In visual decryption, as standard, the texture is described by one or two adjectives, for example: filiform, spongy, radial, etc. [9]. The preliminary stage of visual decryption involves the use of the whole complex of decryption features involving various auxiliary information, namely, the results of previously performed decryptions, various a priori and reference information. At the same time, brightness and structural ones are usually considered as the main decoding features. Moreover, the role of luminance features when identifying objects increases significantly in the case of using data with low spatial resolution. Structural features are most often used to identify objects. They play a major role in working with high spatial resolution data. Often objects are identified based on a combination of luminance and structural features. This approach underlies various GIS technologies. The use of geographic information systems allows thematic processing of image data of high and ultra-high resolution, etc. Even the most modern satellite images will not allow to obtain the necessary information without the fast and reliable decryption methods found in most modern specialized software products. Techniques for automated decoding of images developed by a number of international research institutes and research institutes of Ukraine, including classification, contouring and recognition of images of topographic objects [9, 10, 11, 12], allow to maximize the efficiency of decryption, in comparison with traditional non-automated methods. However, the developers of these techniques emphasize that the presence of a number of interactive processes, such as contouring, searching for reference information, classification of images by texture and brightness characteristics, significantly reduces labor productivity during decryption. Poor adaptation of existing methods of automated decryption to changing shooting and imagery processing conditions, as well as the duration of the processes of preparing and performing decryption, which is sufficiently long, is a constraining factor for effective creation and updating of topographic, situational, digital and electronic maps. Thus, the use of well-known GIS packages, such as, for example, ArcGIS, ERDAS IMAGINE, ENVI, etc. to solve a number of problems associated with the processing of multi-temporal space monitoring data is not always justified. This is primarily due to the large number of functions performed by these packages, a large number of variables, which leads to the accumulation of errors to a loss of accuracy during geometric and sometimes topological operations. Their use gives good results when making decisions on visual assessments, when constructing maps of various scales according to the data of remote sensing of the Earth, etc. Moreover, the systems existing at the moment, in fact, like all multidisciplinary programs, do not allow predicting the further distribution of various types of changes in the state of the monitored object, leading to possible technological and environmental disasters. For this, it is necessary to modernize the process of forming the knowledge base and methods for determining not only the attribute space, but also the vector of parameters, which determines the trend of further changes.

Thus, the result of the analysis of the considered works and the analysis of the existing methods and models used showed that the issues of creating a functionally complete mathematical software that comprehensively implements all stages of processing, from analyzing the source information to obtaining products of output information, have not yet received their comprehensive conclusion. This is especially true when taking into account the features of the analyzed geospatial data and the specifics of their reception during the joint processing of significant amounts of heterogeneous a priori information about the monitoring object. Thus, there is a need to develop additional software and technology tools as well as introduce new methods in the field of automated decryption aimed at the joint processing of space monitoring data, ground-based research and statistical data with elimination of redundant data.

The analysis of the problem that is directly related to the geoinformation processing [13, 14] made it possible to identify the prerequisites for the integrated approach need to create a complex of geoinformation models that would allow to study and analyze the cause-effect relationships of the emergence and development of dangerous environmental situations:

The integration of remote sensing data processing technologies and GIS technologies is one of the methods for the operational analysis and forecasting of the development of anthropogenic, natural and social factors interaction processes [15].

Thus, to solve the aforementioned problems, a concept is proposed for predicting the dynamics of changes according to environmental monitoring, built on the basis of generalizing the theory of forecasting natural, social and economic consequences caused by both natural and anthropogenic changes in ecosystems, using interdisciplinary approaches (geographical, biological, geodetic, chemical and photogrammetric) in the joint analysis of field observation data and the decryption of aerospace sensing. This gives additional opportunities when studying the properties of observation objects to establish cause-effect relationships of the occurrence of anomalies of various origins when building a forecast for their further distribution only on a priori data (in the absence of subsequent images).

For this, when defragmenting raster imagery, a complex approach to the organization and storage of remote sensing data is necessary, taking into account the specifics of their acquisition. This approach allows saving information resources and obtain quick ways to access spatial data, as well as automatically detect new changes in images with their automatic photogrammetric normalization and georeferencing. When carrying out photogrammetric processing, despite the fact that all the received images stored in the database have undergone the procedure of geometric and radiometric correction, it is necessary to additionally carry out the data georeferencing procedure to convert the images to a view independent of the data source [16]. The advantage of geo-referenced images is the ability to conduct various geometric and overlay operations with “criterion trees” [13, 15, 17] built on their basis, and combining satellite imagery data with other data determined on the basis of geographic reference of the studied area. In particular, it is possible to combine data obtained at different times from the same satellite, or data obtained by different remote sensing systems.

Alignment of the images to a single scale, a single coordinate system and establishing quantitative characteristics of distortions caused by both external, random factors, and the image geometry is based on mathematical dependencies – that express the basic geometric properties of a single image, namely, using the formulas (1) for the relationship of the point coordinates on an arbitrary image with the coordinates of the corresponding points on the horizontal image and with the coordinates of the same points on the ground [18].

where x, y –aerial photograph point coordinates;

f – focal length;

a₁, a₂, a₃, b₁, b₂, b₃, c₁, c₂, c₃ – directional cosines;

X, Y – ground point coordinates;

X_s, Y_s, Z_s – photographing pointcoordinates;

It is taken into account that in the general case, the scale changes not only when moving from one point of the image to another, but also when changing directions in the same point. Thus, the scale of a photograph at its different points can be determined by the formula (2) [18]:

where 1/m – vertical aerial photograph average scale on an arbitrary point;

H – photographing height;

α – inclination angle of the aerial photograph;

φ – image direction by which the aerial photograph average scale is determined.

The task of linking the characteristic (nodal) and tying points of the image to the control pointscoordinateson the ground is carried out on the basis of approximate orientation element values indicated in the comments on the image, or in the absence of them, as a result of elementary geometric calculations based on available parameters about the motion of the spacecraft.

However, when implementing the proposed method, a number of difficulties arise associated with the information sources, namely, when the requested images are presented in a synthesized form. First of all, the difficulty is due to the absence of a number of a priori data indicated in the comments to the images, such as: the altitude of the spacecraft, elements of internal and external orientation, allowing to restore a bunch of projection beams that existed at the time of photographing [19].

Thus, when transforming an image for its conversion into a screen coordinate system in a given map projection, it is proposed to perform standard operations – spatial interpolation and brightness pixel values interpolation [20].

Spatial interpolation is performed according to the adjusted values of reference points in a rectangular coordinate system as a result of the photogrammetric processing described above with fitting by the least squares method [19].

In the general case, such a transformation is described by an affine transformation with six independent parameters that correspond to six elementary transformations of the image, i.e., displacement along the x axis, displacement along the y axis, zooming along each of these axes, parallel to the shift of the image boundaries and rotation [13, 19].

This conversion is applied to all pixels in the image.

The final step in image processing is the use of the “nearest neighbor” method for interpolating the brightness values, i.e. recalculation of the brightness of the pixels of the original image into the brightness of transformed pixels. At the same time, 3 main resampling methods were considered: the “nearest neighbor” method, the bilinear interpolation method, the cubic convolution method. The use of the second and third methods has allowed to get very clear pictures that do not have a block structure, which manifests itself when using the method of “nearest neighbor”. However, taking into account the further construction of “criteria trees”, one of the properties of which is the level of homogeneity in each of its fragments, the most appropriate is the use of the “nearest neighbor” method [13, 15, 21].

To achieve the level of homogeneity in automatic segmentation, the use of the brightness segmentation method is proposed. At the same time, other methods for clustering data were considered, determining the differences in the values of segmentation parameters (in accordance with the established criteria when constructing the “criterion trees”), the corresponding boundaries between the selected segments, and contour segmentation (combining adjacent image elements with similar values of the parameter vector (increasing the areas) – clarification of the analyzed area boundaries). It is determined that the use of image segmentation algorithms based on minimizing or maximizing the selected criterion for clustering quality, for example, algorithms that use k-means and dynamic condensation methods, are not always justified. The use of this algorithm (k-intragroup means) in its pure form has a significant drawback – the segmentation results depend on the viewing image elements sequence, as well as the vector of initial approximations of the center of the segments, and the number of distinguished segments. The most appropriate is the use of the ISOMAD algorithm, which is a modification of the k-means algorithm, supplemented by a number of heuristic techniques that allow to adjust the number of allocated segments depending on the intermediate segmentation results. This method ensures the absence of empty classes that do not contain a single vector image and thereby increases the speed [13, 22].

It is known that thematic processing of satellite images is a complex of operations with images, which makes it possible to extract information from them that is of interest in terms of solving various thematic problems. In particular, the effective identification of various kinds of changes on the ground according to the remote sensing of the Earthdata, the identification of these changes, first of all, attributing them to anomalies, and taking the topological characteristics of the object based on the results of decryption. Their automatic detection by assigning a code for each homogeneous fragment during segmentation, according to the hierarchical structure of the “criterion tree” in conformity with a certain set of rules, allows you to identify the appurtenance of the analyzed image fragment to a particular class of objects. In other words, the “criterion tree” is a hypergraphic structure constructed as a result of applying the considered algorithms.

The implementation of this approach allows to generate a rule set that defines segmentation on a relatively small fragment, and then use the resulting set of rules for the entire image or set of images. This approach can significantly reduce processing time and save memory [19].

Additionally, when segmenting, it is necessary to evaluate the intermediate segmentation results. Evaluation can be made based on the analysis of the distances between the centers of segments, intra-segment dispersions and other quantitative criteria, as well as using the connectedness of the selected subsets of image elements and the correspondence of the object localizations allocated during segmentation with a priori information.

This approach to the satellite imagery analysis is the most appropriate, since it involves the analysis of not the entire image, but only the part of it, which, according to a priori data, is of the greatest interest. However, to compare the two images obtained from different sources, it was decided to break the image into cells by geographical coordinates. Then, touse “criterion trees” for mutual analysis.

Using the above algorithm, the images were segmented in accordance with the criteria for determining the potential flooding danger of objects that are potential pollution sources,with the Western Bug River basin used as the example.

Sentinel-2 images with a spatial resolution of 10 m/pixel for the period from 2016 to 2019 were used as initial data for flood zones localization (Fig. 1). In a joint analysis of Earth remote sensing data, as a rule, “criterion trees” are used, built on the basis of time-differentiated imagery and often from several spectral ranges. For the given territory, combinations of the following channels were used: NIR (8 band), Red (4 band), Green (3 band), Blue (2 band) [23].

Figure 1.

In order to determine the maximum possible flooding area during catastrophic floods, data on the hydrological regime of the Western Bug basin and the TIN relief model (Fig. 1) were used, obtained on the basis of radar survey data (SRTM data) [24].

According to indirect decryption features, 14 potential pollution sources are localized within the boundaries of the presented site. Also, according to indirect decryption features, 8 slagheaps, 3 coal processing plant’ssedimentation tanks and 3 illegal landfills were identified (Fig. 2).

Figure 2.

The indicated objects (Fig. 2) are potential pollution sources, polluted by heavy metals, sulfur oxides, phenols, oil products, etc. One of the possible options for pollution to enter the Western Bug river basin is their flooding or partial destruction during high or low floods.

Based on the results of a comprehensive analysis of contact and remote sensing data, forecasted flooding zones during a catastrophic flood were constructed (Fig. 2), a flood which occurs once every 100 years. 2 of the coal factory’s sedimentation tanks, 4 slagheaps and 2 unauthorized landfills fall into the flood zone.

A result of applying this approach to working with spatial data is a significant time gain, obtained when constructing timing series that allow to optimize RAM usage and work not with the whole image, but with individual fragments (homogeneous sheets), as well as use and acquire faster algorithms to access spatial data.

In a joint analysis of Earth remote sensing data, “criterion trees” are used, which are built on the basis of multi-time images and often from several spectral ranges.

A comparative analysis of the classical and proposed processing methods as applied to the satellite image fragment considered above, corresponding to the criterion element of 20×20 pixels. A comparative analysis of the results allows us to draw the following conclusions. When processing the image in accordance with standard approaches, information corresponding to 1224 bytes is used and contains: the size of the processed fragment in three RGB color channels – 20×20×3=1200 pixels, the coordinates of the selected fragment – 16 bytes, the number of rows and columns of the fragment – 8 bytes. The total capacity is 1224 bytes. The sheet of the “criterion tree” corresponding to the same fragment contains the address of the final sheet — approximately 50 bytes, the mathematical expectation of a homogeneous fragment in each RGB color channel is 3×8 bytes, the standard deviation in the color channels is 3×8 bytes, the brightness change gradient inside the fragment is 3×8 bytes, as a texture attribute, the number of brightness drops within the fragment is considered to be 1×8 bytes (analyzed in gradations of grey). The total size is 150 bytes. Thus, the application of the “criterion tree” method uses 8 times less RAM when processing a minimal element, which allows to save RAM and get faster algorithms for accessing spatial data due to the hierarchical organization of data and working only with homogeneous sheets. The above example illustrates this point. It was determined that the area of the selected fragment in the image occupies 245×245=60515 pixels, and the leaf of the “tree” corresponding to this fragment occupies 3595 pixels. Figure 3 presents a comparative histogram of the amount of occupied memory when working with the image according to the classical methods versus the “criterion trees” method.

Figure 3.

In the Figure 3, stage 1 corresponds to the initial selection of homogeneous segments, when processing is performed throughout the entire field of the image, stage 2 corresponds to work with segments selected in the first stage. As a result, taking into account the fact that the first time when isolating homogeneous fragments the processing is carried out with the whole image, the time gain when constructing time series and performing overlay operations with time-differentiated “criterion trees” reaches 53% compared to standard processing methods: (100+3595/60515)t/200=0.53t (t is the processing time using standard methods). The new trees formed as a result of the transformation, built only in the area of interest revealed that, as a result of overlay operations and changes with geo-referencing, which allow to highlight previously existing objects with better clarity.

Thus, the method of conducting overlay operations with data obtained by means of different survey equipment and at different times, indexing all combined fragments of superimposed images that was considered in this work allows the creation of a general “tree” that tracks the dynamics of changes in the characteristics of the observed object and contains the attributes of all the analyzed fragments.

The analysis of determining the natural and anthropogenic processes development dynamics according to remote sensing data should be carried out using the overlay operation with “criteria trees”, built for each shot. The differences obtained as a result of pixel-by-pixel subtraction of localized geo-referenced fragments of images of the observed territory on the analyzed image and images from the existing archive. At the same time, changes in the decoding attributes are associated with changes in attribute information at a fixed point in time. By superimposing similar fragments over the entire fixed time series, a “criterion tree” is formed, which displays the combined areas in its structure in different-time images. Finding these fragments is carried out by traversing all branches of each “tree” up to the minimum fragments of the image with the greatest number of attributes. In the absence of branching in some “criterion tree”, the quantitative value of the attribute is transferred to all subsequent levels. This technique allows you to get a single “tree”, which displays the attributes of all the studied fragments.

One of the tasks of space monitoring is to establish the causes of the dynamics of objects of different genesis in space and time, as well as to determine the types and degrees of interconnections between neighboring areas. The above method using “criteria trees” allows you to solve this problem. So, with thematic processing of remote sensing data, vector layers that occupy a significantly smaller amount of memory than raster ones are created. In these layers, the boundaries of various objects with attribute information are displayed. To identify the relationships between the objects, the spatial indices of the selected fragments obtained during the overlay are analyzed. This allows you to speed up the search for features with a complex shape and it is especially effective when processing polygon features. This approach allows you to quickly identify various relationships between objects and predict changes in their dynamics.

When determining the area of an object which has the shape of a polygon, defined by a sequence of vertices, it is advisable to use an algorithm based on its partition into trapezoidal lines bounded by the polygon segment line, perpendiculars dropped from the vertices of the segment to the corresponding coordinate axes. In this case, the areas of all trapezoids that make up the object for all segments of the polygon are calculated and then summed. It should be noted that in real conditions, due to the use of rectangular rectilinear coordinate systems with large numerical values, when the summed areas of the trapeziums are repeatedly added, a computational error will accumulate. The relative error is especially high for small polygons. In this case, the linear transformation of the polygon is applied to the new coordinate system in which there are no large values, and the calculation of the area within it. Next, the area in the source system is calculated. When determining the border of two adjacent polygons, the area is calculated once. Then, for the right polygons, this area is summed with the plus sign, and for the left – with the minus sign. This is due to the fact that the proposed method of splitting image fragments into a “criterion tree” is a graph structure and polygons are formed from arcs. Coding is relative to the direction of digitization of the arc.

Since all database images are subject to the same photogrammetric processing and stored as “criterion trees” before the thematic analysis, the analysis of the reference image trees and the tree built on the new image data makes it possible to detect various discrepancies, calculate the area of these changes, accurately determine the coordinates, etc. New trees formed as a result of this transformation, built only in the area of interest resulting from overlay operations, with geo-referencing, allow to identify the preexisting objects with much greater clarity. To assess the state changes in the monitored object, a pixel-by-pixel subtraction algorithm can be used, which consists in subtracting from the brightness value of each pixel in a fragment of one image the brightness value of the corresponding pixel in a fragment of another image, which is then combined with the first. Objects, reflectance value of which has changed slightly between the two shots will be painted in light gray tones in the differential image, and dark and bright areas of the image will correspond to significant changes in reflectance.

The algorithms for working with “criterion trees” discussed above are the basis for more complex algorithms that allow reconstructing spatial relationships between fragments of “trees” superimposed on each other. An analysis of the problems that have arisen touches on the issues of combining and creating new analysis criteria for the already constructed differential “criterion tree”, which contains a large number of variables and factors that interact with each other and respond to changes in each other variables, etc., which makes it possible to track the change dynamics in the monitoring object.

With this approach, the main idea of carrying out overlay operations on “criterion trees”, constructed from time-differentiated satellite images and the corresponding attribute information, is to combine similar tree fragments in the analyzed (input) image and on the tree from the images archive and obtain a new “criterion tree”, which indexes all aligned fragments of overlaid image snapshots. To do this, all the trees are evaluated simultaneously, following the branches existing in all the trees, and corresponding to the minimal fragments of images with the greatest number of attributes.

All images from their respective database, including the current one before the thematic analysis, undergo the same photogrammetric processing and are stored in the form of “criterion trees”. Analysis of the reference image trees and the tree built on the new image and various operations with them (overlaying, subtracting, etc.) makes it possible to detect various changes (discrepancies), calculate the area of these changes, accurately determine the coordinates, etc.

The new trees formed as a result of the transformation, built only in the area of interest revealed by changes in overlay operations, with a geographical reference, allow for a much clearer distinction of previously existing objects [10].

To assess the change in the state of the monitored object, a new differential “criterion tree” is formed as a result of using the pixel-by-pixel subtraction algorithm, which consists of subtracting from each pixel’s brightness value in one image fragment the brightness value of the corresponding pixel in a fragment of another image, which is then combined with the first. Objects whose reflectance have changed slightly between the two shots will be painted in light gray tones in the differential image, and dark and bright areas of the image will correspond to significant changes in reflectance.

When working with a differential “criterion tree” constructed at the previous stage, which is considered the reference at the current stage working with the image and the newly formed one, the addition operation is expedient, which is averaging and allows to further reduce the noise level of the selected fragment containing the analyzed object. This operation is especially often used in spatial filtering to select objects of a certain scale [15].

When forming a vector of parameters characterizing the monitored object state, it is recommended to use the division operation. This allows to calculate many spectral features characterizing the slightest changes in the spectral reflectivity of different ground covers and changing a number of decoding features. The ratio of values from different spectral ranges gives important information about objects, since various objects have high reflectivity and absorb radiation well in different regions of the spectrum. The result of the division operation application is characterized by a vegetation index, which is also an integral part of the parameter vector which characterizes the observation object. An additional advantage of using this algorithm is that the resulting data set contains only relative, not absolute, brightness values. Due to this, data analysis to a much lesser extent depends on such nuances as the illumination of the scene, which arise due to the features of the terrain.

The algorithms that work with “criterion trees” discussed above are the basis for more complex algorithms that allow to restore spatial relationships between stacked tree fragments, in particular, perform overlaying operations of polygons to clarify the boundaries of the identified anomaly and determine the zone of changes in the period of time between acquired imagery.

To clarify the boundaries of the revealed anomaly, it is recommended to use polygon operations inside the tree. Overlaying two polygons results in a graphical interpretation of the union or intersection of sets.

To estimate the relationship between the sections of the criterion tree, characterizing the identified changes which cannot be classified without additional analysis, with indicators of adjacent fragments, it is necessary to determine the correlation coefficient. This stage is especially important for the case of violation of topological characteristics after automated processing or in conditions of uncertainty caused by the use of smoothing filters in the primary processing of satellite images. For example, if the threshold value is incorrectly selected during the smoothing procedure, there is a chance of removing small decryption features. Further, when making decisions about the state of the monitoring object and analyzing the causes that led to various changes in its geometric or topological characteristics, the loss of such signs may lead to incorrect conclusions. To do this, according to well-known formulas, the correlation coefficient is calculated. Correlation is determined from the inequality -1 ≤ r_xy ≤1 by the correlation coefficient value: if r_xy < 0, then there is a negative correlation between the variables, and if r_xy > 0, it is positive. The values |r_xy|= 1 correspond to the existence of a functional linear relationship of variables, direct at r_xy = 1 and inverse at r_xy = -1. The closer |r_xy| is to one, the stronger the connection, and the closer to zero the weaker. The value r_xy characterizes the intensity of only the linear relationship between the variables [15, 19].

In the case of the existence of significant nonlinear correlation relationships, it makes sense to study the correlation only for those pairs of variables whose linear correlation coefficient absolute value exceeds 0.6.

Similar operations to determine the degree of correlation are carried out not only for satellite imagery data, but also for revealing stochastic relationships between the zonal brightness levels of satellite imagery and environmental parameters, probabilistic characteristics of the factors affecting the anomaly, determining the nature of these relationships and assessing their significance.

After revealing significant correlations between any analyzed indicators, dependencies are constructed for them by regression analysis methods.

At the same time, on the first stage, regression dependencies are constructed for the case when the parameters of satellite images are dependent variables, and explanatory variables – indicators of the environmental state and regression dependences of the satellite images parameters of various fragments at the junction of the old anomalyboundary and the new one for the constructed differential “criterion tree”.