CN113203399A - Underground space resource quantity analysis method - Google Patents

Underground space resource quantity analysis method Download PDF

Info

Publication number
CN113203399A
CN113203399A CN202110411708.6A CN202110411708A CN113203399A CN 113203399 A CN113203399 A CN 113203399A CN 202110411708 A CN202110411708 A CN 202110411708A CN 113203399 A CN113203399 A CN 113203399A
Authority
CN
China
Prior art keywords
building
shadow
extraction result
height
satellite
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110411708.6A
Other languages
Chinese (zh)
Other versions
CN113203399B (en
Inventor
董杰
徐红兵
夏伟强
王建收
刘海松
郝胜优
于晓曦
付佳妮
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Qingdao Geological Engineering Survey Institute
Original Assignee
Qingdao Geological Engineering Survey Institute
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Qingdao Geological Engineering Survey Institute filed Critical Qingdao Geological Engineering Survey Institute
Priority to CN202110411708.6A priority Critical patent/CN113203399B/en
Publication of CN113203399A publication Critical patent/CN113203399A/en
Application granted granted Critical
Publication of CN113203399B publication Critical patent/CN113203399B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C11/00—Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
    • G01C11/04—Interpretation of pictures
    • G01C11/06—Interpretation of pictures by comparison of two or more pictures of the same area
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C11/00—Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
    • G01C11/04—Interpretation of pictures
    • G01C11/30—Interpretation of pictures by triangulation
    • G01C11/34—Aerial triangulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00—Scenes; Scene-specific elements
    • G06V20/10—Terrestrial scenes
    • G06V20/176—Urban or other man-made structures

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Multimedia (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Theoretical Computer Science (AREA)
  • Image Processing (AREA)

Abstract

The application relates to a method for analyzing underground space resource quantity. The method comprises the following steps: obtaining a satellite image of an area to be analyzed; preprocessing the satellite image to obtain a processed satellite image; extracting buildings from the processed satellite images by adopting morphological building indexes and deep learning; extracting shadows from the processed satellite images according to the characteristic indexes; analyzing according to the sun azimuth line and the shadow extraction result to obtain the shadow length of each shadow in the shadow extraction result; building height inversion work is carried out by utilizing the building height inversion model, and the building height of each building in the building extraction result is determined; the building height of each building is converted into the building influence depth of each building, the utilization condition of the underground space resource quantity of each layer is analyzed by utilizing the building influence depth and the area of each building, the analysis result of the underground space resource quantity is obtained, and the problems that the utilized resource quantity is counted on the spot and the investigation cost is high are solved.

Description

Underground space resource quantity analysis method
Technical Field
The application relates to the technical field of remote sensing, in particular to an underground space resource quantity analysis method.
Background
The planning and construction of underground space resources in a city need to be realized, and reasonable layout and development are achieved.
At present, most scholars take resource quality evaluation as a research target in the aspect of research on the amount of resources which utilize underground space; some scholars develop specific evaluation means and technical method researches; some researchers develop evaluation system research; however, the amount of resources used is rarely investigated, and at present, the amount of resources used is mainly counted in the field in a manual visit investigation mode, which causes a problem of high investigation cost.
Disclosure of Invention
In view of the above, it is necessary to provide an underground space resource amount analysis method that can solve the problem of high investigation cost.
A method of underground spatial resource volume analysis, the method comprising:
acquiring a satellite image of an area to be analyzed;
preprocessing the satellite image to obtain a processed satellite image;
extracting buildings from the processed satellite images by adopting morphological building indexes and deep learning to obtain building extraction results;
shadow extraction is carried out on the processed satellite images according to the characteristic indexes to obtain shadow extraction results;
analyzing according to the sun azimuth line and the shadow extraction result to obtain the shadow length of each shadow in the shadow extraction result;
according to the building extraction result and the shadow length of each shadow in the shadow extraction result, utilizing a building height inversion model to perform building height inversion work, and determining the building height of each building in the building extraction result;
and converting the building height of each building into the building influence depth of each building, and carrying out hierarchical underground space resource utilization analysis by using the building influence depth and area of each building to obtain an underground space resource analysis result.
In one embodiment, the step of extracting a building from the processed satellite image by using a morphological building index and deep learning to obtain a building extraction result includes:
extracting buildings from the processed satellite images by adopting a morphological building index to obtain a first extraction result;
extracting buildings from the processed satellite images by adopting deep learning to obtain a second extraction result;
fusing the first extraction result and the second extraction result to obtain a fused extraction result;
and optimizing the roof shape of the fused extraction result to obtain a building extraction result.
In one embodiment, the step of extracting shadows from the processed satellite images according to the feature indexes to obtain a shadow extraction result includes:
comparing the intensity of the shadow expression force of each characteristic index, and determining the characteristic index with the highest intensity of the shadow expression force as the characteristic index of shadow extraction;
adopting the characteristic index of the shadow extraction to carry out shadow extraction work, and obtaining a preliminary shadow extraction result;
and optimizing the preliminary shadow extraction result by adopting clustering processing to obtain a shadow extraction result.
In one embodiment, the step of determining the building height of each building in the building extraction result by performing building height inversion work using a building height inversion model according to the building extraction result and the shadow length of each shadow in the shadow extraction result includes:
determining a building to be subjected to building height inversion according to the building extraction result;
and based on the spatial position relation among the sun, the satellite and the building and the relation between the building height and the shadow of the building, performing the inversion work of the building height by using a building height inversion model, and determining the building height of each building in the building extraction result.
In one embodiment, when the azimuth of the sun and the satellite are the same, the altitude angle of the sun is smaller than that of the satellite, and the projection of the sun ray on the ground is perpendicular to the building trend, the height inversion formula of the building height inversion model is as follows:
H＝S·(tanα·tanβ)/(tanα-tanβ)
where α is the altitude of the satellite, β is the altitude of the sun, S is the shadow length of the building in the shadow extraction result, and H is the building height of the building.
In one embodiment, when the building trend is not perpendicular to the projection of the solar ray on the ground, and the projection of the solar ray on the ground and the projection of the satellite sight line on the ground are located on two sides of the building, the height inversion formula of the building height inversion model is as follows:
H＝S·tanβ/sinγ
wherein gamma is the included angle between the building trend and the solar ray projection direction.
In one embodiment, when the building trend is not perpendicular to the projection of the solar ray on the ground, and the projection of the solar ray and the projection of the satellite sight line are located on the same side of the building, the height inversion formula of the building height inversion model is as follows:
H＝S·sinα·sinβ/(sinα·cosβ·sinγ-cosα·sinβ·sinδ)
wherein, delta is an included angle between the trend of the building and the projection direction of the satellite sight.
In one implementation, the step of converting the building height of each building into the building influence depth of each building, and performing hierarchical underground space resource utilization analysis by using the building influence depth and area of each building to obtain the underground space resource analysis result includes:
converting the building height of each building into the building influence depth of each building according to the conversion relation table of the building height and the building influence depth;
and analyzing the utilization condition of the layered underground space resource quantity according to the building influence depth and area of each building to obtain an underground space resource quantity analysis result.
In one implementation, the step of performing hierarchical underground space resource utilization analysis according to the building influence depth and area of each building to obtain an underground space resource analysis result includes:
according to the building influence depth of each building, respectively taking the building influence depth as 0m-10m, 10m-30m, 30m-50m and 50m-100m as calculation units;
and performing superposition analysis on the utilization conditions of the underground space resources of each computing unit to obtain an underground space resource analysis result.
According to the method for analyzing the underground space resource quantity, the satellite image of the area to be analyzed is obtained; preprocessing the satellite image to obtain a processed satellite image; extracting buildings from the processed satellite images by adopting morphological building indexes and deep learning to obtain building extraction results; shadow extraction is carried out on the processed satellite images according to the characteristic indexes to obtain shadow extraction results; analyzing according to the sun azimuth line and the shadow extraction result to obtain the shadow length of each shadow in the shadow extraction result; according to the building extraction result and the shadow length of each shadow in the shadow extraction result, utilizing a building height inversion model to perform building height inversion work, and determining the building height of each building in the building extraction result; the building height of each building is converted into the building influence depth of each building, the utilization condition of the underground space resource quantity of each layer is analyzed by utilizing the building influence depth and the area of each building, the analysis result of the underground space resource quantity is obtained, and the problems that the utilized resource quantity is counted on the spot and the investigation cost is high are solved.
Drawings
FIG. 1 is a schematic flow chart of a method for analyzing the amount of resources in a subterranean space according to one embodiment;
FIG. 2 is a schematic diagram illustrating a building extraction result of a partially processed satellite image according to an embodiment;
FIG. 3 is a schematic diagram of a partially processed satellite image according to one embodiment;
FIG. 4 is a diagram illustrating shadow extraction results of partially processed satellite images according to an embodiment;
FIG. 5 is a schematic diagram illustrating a geometrical relationship between the same side of the star and the same side of the sun under ideal conditions in the method for analyzing the amount of underground space resources in one embodiment;
FIG. 6 is a schematic diagram illustrating an ideal case of geometric relationships between the different sides of the satellite and the sun in the method for analyzing the amount of underground space resources in one embodiment;
FIG. 7 is a schematic diagram illustrating a geometrical relationship between the same sides of the star and the sun under non-ideal conditions in the method for analyzing the amount of underground space resources in one embodiment;
fig. 8 is a schematic diagram of geometric relationships of the satellite-sun different-side time under the non-ideal condition in the underground space resource amount analysis method in one embodiment.
Detailed Description
In order to make the objects, technical solutions and advantages of the present application more apparent, the present application is described in further detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the present application and are not intended to limit the present application.
In one embodiment, as shown in fig. 1, there is provided a method for analyzing the amount of underground space resources, comprising the steps of:
in step S220, a satellite image of the area to be analyzed is obtained.
The region to be analyzed refers to a geographical region in which underground space resource amount analysis is required. The satellite image is a ground image obtained by a satellite, has the characteristics of high resolution, wide coverage area, short update period, rich spectral characteristics and the like, and is called as a senyei eye and a perspective eye in space.
Step S240, preprocessing the satellite image to obtain a processed satellite image.
The satellite image is preprocessed in at least one of radiometric calibration, atmospheric correction, orthometric correction, image fusion and image cropping.
And step S260, extracting the building from the processed satellite image by adopting the morphological building index and the deep learning to obtain a building extraction result.
Wherein, the morphological building index and the deep learning are adopted to respectively carry out the building extraction work from the processed satellite images, the two extraction results are fused, the fused result is optimized in the roof shape, and finally the obvious, clear and complete building single body extraction result with clear boundary is obtained, the building single body extraction result is the building extraction result, and is shown as the building extraction result schematic diagram of the partially processed satellite images shown in fig. 2.
In one embodiment, the step of extracting the building from the processed satellite image by using the morphological building index and the deep learning to obtain the building extraction result comprises: extracting buildings from the processed satellite images by adopting a morphological building index to obtain a first extraction result; extracting buildings from the processed satellite images by adopting deep learning to obtain a second extraction result; fusing the first extraction result and the second extraction result to obtain a fused extraction result; and optimizing the roof shape of the fused extraction result to obtain a building extraction result.
The basic idea of the Morphological Building Index (MBI) is to establish a relationship between implicit features of a building (such as brightness, size and contrast) and morphological operators (such as top hat obtained through reconstruction, granularity determination and directionality), and extract the building from the processed satellite image according to the relationship established between the implicit features of the building and the morphological operators of the building to obtain a first extraction result.
The morphological building index is based on the fact that the edge spectrum of a building varies greatly, while the interior spectrum of the building varies little. The morphological building index is constructed by considering the characteristics of the shape, the direction, the brightness, the contrast ratio and the like of a building, and the method comprises the following concrete implementation steps:
the first step, calculating the brightness b of the processed satellite image, wherein the brightness of the processed satellite image comprises the brightness value of each pixel x in the processed satellite image, and the brightness value of the pixel x is calculated as follows: the maximum spectral value of the multispectral band for pixel x is recorded as the luminance value of pixel x:
b(x)＝max1≤k≤K(bandk(x)) (1)
wherein b (x) is the brightness value of pixel x, and K is the number of visible light spectrum bands, bandk(x) The spectral value of the kth spectral band pixel x is represented, and since the visible light band in the multispectral bands has a large influence on the spectral information of the building, the maximum value of each pixel in the visible light band is selected as the brightness value of the pixel x.
Step two, MBI construction: DMP using top-hat represents the spectral structure characteristics (e.g., contrast, size, directionality) of the building through a series of linear SE reconstructions. The morphological operators used to construct the MBI are as follows:
1. morphological white cap reconstruction (W-TH)
Figure BDA0003024087590000071
Wherein W-TH (d, s) is morphology of the processed satellite image after reconstructionThe white cap is a white cap, and the white cap,
Figure BDA0003024087590000072
the luminance b of the processed satellite image is subjected to morphological opening operation, d is the direction of the linear structural element, s is the scale of the linear structural element, and b is the luminance of the processed satellite image. Because the spectral information of buildings and roads are relatively similar, the directions of the roads are less, the roads generally extend along one or two directions, the directions of the buildings are more, and linear structural elements with multiple directions and scales are selected for well distinguishing the roads from the buildings.
2. Calculating the morphological section MP:
Figure BDA0003024087590000073
wherein, MPW-TH(d, s) is the morphological section of the processed satellite image, MPW-THAnd (d,0) is a morphological section of the processed satellite image when the scale of the linear structural element is 0.
3. Calculating Differential Morphology Profile (DMP):
DMPW-TH(d,s)＝|MPW-TH(d,(s+Δs))-MPW-TH(d,s)| (4)
wherein, DMPW-TH(d, s) represents the differential morphological section of the processed satellite image, Δ s is the profile interval and the range is smin≤s≤smax。
4. Calculating a morphological building index MBI:
Figure BDA0003024087590000081
where MBI is the morphological building index of the processed satellite image, D represents the number of directions of the contour, and S represents the number of scales of the contour. Since the increase in the number of directions D of the contour does not improve the accuracy of building detection, the present study considered four directions (0 °, 45 °, 90 °, 135 °). Structural element smin、smaxAnd ΔThe size of s should be determined based on the spatial characteristics and spatial resolution of the building, and in one embodiment smax＝52、smin2 and 5. The dimension S of the profile is calculated as:
Figure BDA0003024087590000082
the construction of MBI is based on building structures having large values in most directions of the white top hat DMP histogram because they have high local contrast in these directions, and therefore, the building extraction from the processed satellite imagery is performed with often reference to large MBI values, and the first extraction result is obtained.
Extracting the building from the processed satellite image by adopting deep learning to obtain a second extraction result, specifically: the total area of 5488 buildings is 84.0095Km2The method comprises the steps of establishing a training set for satellite images in the area, enhancing the training set through operations such as rotation, mirroring, blurring and illumination adjustment, and finally obtaining a data set with the quantity of 10 ten thousand 256 × 256 pictures, wherein the data set only comprises two types including buildings and non-buildings, the satellite images cover building areas with different densities, the building areas comprise areas with dense buildings and areas with sparse buildings, training according to the training set by using a Deeplab-V3 network to obtain a building extraction model, inputting the processed satellite images into the building extraction model, and extracting corresponding buildings, namely second extraction results, according to the input building extraction model and the input processed satellite images.
Deep is an algorithm designed by google on the basis of FCN. FCN to obtain a more dense position sensitivity map, a 500x500 input image is directly padded with 100 on the first convolutional layer, and finally a 16x16 position sensitivity map is obtained on the seventh convolutional layer. Deeplab uses a more subtle approach here: the step size of the 14 and 15 pooling layers of the VGG network is changed from 2 to 1, plus 1 padding. It is a modification that the total step of the VGG network is changed from 32 to 8, so that the seventh convolutional layer can obtain a 67x67 position sensitive map when the input image is 514x514, and the map is denser than the FCN.
Deeplab-V3 is still mainly researched in cavity convolution, and is designed to capture a multi-scale background by adopting multi-proportion perforated convolution cascade or parallel, a building extraction model obtained by utilizing Deeplab-V3 network training is utilized to transfer images through a plurality of convolutions, the convolutions are all set to have the step length of 1, the same filling is used instead of using pooling and full connection layers, and through the processing, each convolution retains the input spatial dimension, and a segmentation model (namely the building extraction model) for extracting the building is finally obtained through a plurality of stacked convolutions.
Step S280, shadow extraction is carried out from the processed satellite images according to each characteristic index, and a shadow extraction result is obtained.
Wherein, each characteristic index is an index such as a shadow extraction index, a green channel vegetation index (GNDVI), a fourth wave band of the image and the like. The shadow expressive forces of the characteristic indexes are compared, and the characteristic index with the clearest shadow expressive force is selected for shadow extraction. However, the method for extracting the shadow is based on the pixel, so that a certain salt and pepper phenomenon exists in the extraction result, small pattern spots exist in the extracted primary shadow extraction result, meanwhile, a part of small shadows can be generated due to the fact that the roof is not flat and a part of small buildings such as staircases exist, the small pattern spots and the small shadow shadows are removed or reclassified, and clustering processing is adopted for classification post-processing of the shadow extraction result. As shown in the schematic diagram of the partially processed satellite image shown in fig. 3, shadows irradiated by the sun exist in the processed satellite image, and shadow extraction is performed from the processed satellite image according to each feature index, so as to obtain a schematic diagram of a partial shadow extraction result shown in fig. 4.
In one embodiment, the step of extracting the shadow from the processed satellite image according to each feature index to obtain the shadow extraction result includes:
comparing the intensity of the shadow expression force of each characteristic index, and determining the characteristic index with the highest intensity of the shadow expression force as the characteristic index of shadow extraction; adopting the characteristic index of shadow extraction to carry out shadow extraction work, and obtaining a preliminary shadow extraction result; and optimizing the preliminary shadow extraction result by adopting clustering treatment to obtain a shadow extraction result.
Among them, the characteristic index having the highest intensity of the shadow expression is the characteristic index in which the shadow expression is most clearly apparent. And optimizing the extracted shadow by adopting clustering (column). Clustering (clump) is the use of mathematical morphological operators (erosion and dilation) to cluster and merge adjacent similarly classified regions. Classified images often lack spatial continuity (the presence of blobs or holes in the classified region). Although low-pass filtering can be used for smoothing the images, the class information is often interfered by codes of adjacent classes, and the clustering process solves the problem.
And step S300, analyzing according to the sun azimuth line and the shadow extraction result to obtain the shadow length of each shadow in the shadow extraction result.
And drawing a sun azimuth line according to the satellite image data, and obtaining the shadow length of each shadow in the shadow extraction result according to the intersection of the sun azimuth line and each shadow in the shadow extraction result.
And S320, performing building height inversion work by using the building height inversion model according to the building extraction result and the shadow length of each shadow in the shadow extraction result, and determining the building height of each building in the building extraction result.
The method comprises the steps of utilizing a building height inversion model to conduct building height inversion work according to a building extraction result and shadow lengths of shadows in the shadow extraction result, establishing corresponding models for simulation according to four conditions of different position relations among the sun, a satellite and a building, utilizing data such as a sun azimuth angle, a sun altitude angle, a satellite azimuth angle, a building trend line and the like to conduct calculation according to the building height inversion model to complete height inversion, correspondingly endowing building height to each building of the building extraction result, and determining the building height of each building in the building extraction result.
In one embodiment, the step of determining the building height of each building in the building extraction result by performing building height inversion work using the building height inversion model according to the building extraction result and the shadow length of each shadow in the shadow extraction result includes:
determining a building to be subjected to building height inversion according to a building extraction result; and based on the spatial position relation among the sun, the satellite and the building and the relation between the building height and the shadow of the building, performing the inversion work of the building height by using the building height inversion model, and determining the building height of each building in the building extraction result.
The building height inversion method is based on the spatial position relationship among the sun, the satellite and the building, and the relationship between the building height and the shadow, and can be divided into four conditions that the directions of the sun and the satellite are the same or different, and the sun and the satellite are on the different side and the same side of the building. Ideally, when the sun and the satellite are in the same plane in the vertical plane, there is no shadow in the case when the azimuth and the altitude of the sun and the satellite are the same, or the azimuth of the sun and the satellite are the same and the altitude of the sun is greater than or equal to the altitude of the satellite.
Establishing different height inversion models according to different conditions to realize accurate inversion of buildings, wherein the height inversion models comprise the following height inversion formulas:
in one embodiment, when the azimuth of the sun and the satellite are the same, the altitude of the sun is smaller than that of the satellite, and the projection of the sun ray on the ground is perpendicular to the building direction, the height inversion formula of the building height inversion model is as follows:
H＝S·(tanα·tanβ)/(tanα-tanβ)
where α is the altitude of the satellite, β is the altitude of the sun, S is the shadow length of the building in the shadow extraction result, and H is the building height of the building.
As shown in fig. 5, when the sun and the satellite have the same azimuth and the altitude of the sun is smaller than that of the satellite, and the projection of the sun ray on the ground is perpendicular to the building trend, the relationship among the sun, the satellite, the building and the shadow is as follows:
S＝L2-L1＝H/tanβ-H/tanα (7)
where α is the altitude of the satellite, β is the altitude of the sun, L2Is the actual shadow length of the building, L1The length of the shadow occluded by the building, S the length of the shadow of the building in the shadow extraction result, and H the building height of the building. Once the shadow length S of the building in the shadow extraction result is obtained, the building height of the building can be calculated according to the following formula:
H＝S·(tanα·tanβ)/(tanα-tanβ) (8)
if the satellite is in a front view posture during shooting, the satellite can observe the whole shadow of the building, and shadow parts which cannot be observed by the satellite do not exist, so that L1＝0，L2When it is S
H＝S·tanβ (9)
As shown in fig. 6, when the sun and the satellite are located on different sides of the building, and the projection of the sun ray on the ground is perpendicular to the direction of the building, at this time, the whole shadow of the building can be seen, and the building has the shadow, then:
building shadow length:
L2＝S＝H/tanβ (10)
building height:
H＝S·tanβ (11)
the two situations only consider two special situations that the satellite and the sun are positioned on different sides or the same side of the building when the satellite, the building and the sun are positioned in the same vertical plane, and actually, the satellite and the sun have a certain included angle, and the trend of the building is not necessarily perpendicular to the projection of the sun ray on the ground. Taking into account the sun and satellite orientations and their relationship to the building strike, the building height calculation can take advantage of the following non-idealities:
specifically, the method comprises the following steps: when the trend of the building is not perpendicular to the projection of the solar ray on the ground, and the projection of the solar ray on the ground and the projection of the satellite sight line on the ground are positioned on two sides of the building, the height inversion formula of the building height inversion model is as follows:
H＝S·tanβ/sinγ
wherein gamma is the included angle between the building trend and the solar ray projection direction.
As shown in fig. 7, when the projection of the sun ray on the ground and the projection of the satellite sight line on the ground are located on two sides of the building, most of shadows (including the main shadow and the falling shadow) can be clearly displayed on the remote sensing image, and at this time, the shadows used for calculating the height of the building are all the falling shadows, and the geometric relationship is as follows:
sinγ＝L2s (two type)
H＝L2Tan beta (Shisansan type)
Wherein S is the shadow length of the building in the shadow extraction result, L2The actual shadow length of the building, H the building height and gamma the included angle between the building trend and the solar ray projection direction.
H ═ S · tan β/sin γ (shijia shima)
When the building trend is not perpendicular to the projection of the solar ray on the ground, and the projection of the solar ray and the projection of the satellite sight line are positioned on the same side of the building, the height inversion formula of the building height inversion model is as follows:
H＝S·sinα·sinβ/(sinα·cosβ·sinγ-cosα·sinβ·sinδ)
wherein, delta is an included angle between the trend of the building and the projection direction of the satellite sight.
As shown in fig. 8, when the sun ray projection and the satellite sight line projection are located on the same side of the building, the light-facing surface of the building is imaged, so that a part of the shadow is blocked, and the visible shadow on the satellite image is only a part of the home shadow of the building, and the geometric relationship is as follows:
H＝S1·tanβ＝S2·tanα (15)
L2＝S1sin gamma (land picking up type)
L1＝S2Sin delta (Ji Qi type)
Wherein, delta is the included angle between the building trend and the satellite sight line projection direction, S1Shadow length for sun casting direction, S2Casting a shadow length for the satellite, again because
S＝L2-L1(Shibai type)
Therefore, it is not only easy to use
H ═ S · sin α · sin β/(sin α · cos β · sin γ -cos α · sin β · sin δ) (jijiujiujiu)
In the formula: s is the shadow length of the building in the shadow extraction result (i.e. the visible shadow length perpendicular to the building line on the satellite image).
Figure BDA0003024087590000131
Figure BDA0003024087590000132
In the formula: a. the1Is the sun azimuth, A2Is the satellite azimuth, (x)1，y1)、(x2，y2) The coordinates of two points on the building trend line or the parallel line of the building trend line are shown.
Step S340, converting the building height of each building into the building influence depth of each building, and analyzing the utilization condition of the layered underground space resource quantity by using the building influence depth and area of each building to obtain an underground space resource quantity analysis result.
Taking a high-resolution second image near the beginning of 9 months at the bottom of 2016 and 8 months in Qingdao city of six Qingdao cities, which is based on central longitude and latitude of (120.3E, 36.2N), (120.5E, 36.1N), (120.2E, 36.0N), (120.5E, 35.9N), (120.3E, 36.4N), (120.6E, 36.3N), as a satellite image of an area to be analyzed, preprocessing the satellite image to obtain a processed satellite image, performing shadow extraction from the processed satellite image according to each characteristic index to obtain a shadow extraction result, analyzing according to a solar square line and the shadow extraction result to obtain shadow lengths of each shadow in the shadow extraction result, performing building height inversion work by using a building height inversion model according to the shadow lengths in the building extraction result and the shadow extraction result, and determining the building height of each building in the building extraction result, and comparing the building height with the actual building height of each building to obtain a building height inversion accuracy evaluation table 1.1.
Figure BDA0003024087590000141
Therefore, the building height inversion method can finish the building height inversion work more accurately, and obtain more accurate building height.
In one embodiment, the step of converting the building height of each building into the building influence depth of each building, and performing hierarchical underground space resource utilization analysis by using the building influence depth and area of each building to obtain the underground space resource analysis result includes: converting the building height of each building into the building influence depth of each building according to the conversion relation table of the building height and the building influence depth; and analyzing the utilization condition of the layered underground space resource quantity according to the influence depth and the influence area of the buildings of each building to obtain an underground space resource quantity analysis result.
In one embodiment, the step of analyzing the utilization condition of the layered underground space resource amount according to the building influence depth and the area of each building to obtain the analysis result of the underground space resource amount comprises the steps of respectively taking the building influence depths of 0-10m, 10-30 m, 30-50 m and 50-100m as calculation units according to the building influence depth of each building; and carrying out superposition analysis on the utilization conditions of the underground space resources on each computing unit to obtain an underground space resource analysis result.
The method comprises the steps of converting the building height of each building into the building influence depth according to a conversion relation table 1.2 of the building height and the building influence depth, then carrying out superposition analysis by using the building influence depth and area and taking 4 layers of 0m-10m, 10m-30m, 30m-50m and 50m-100m as calculation units, carrying out layered underground space resource utilization statistical calculation, and obtaining an underground space resource analysis result.
Table 1.2: conversion relation table of building height and building influence depth
Figure BDA0003024087590000151
Specifically, the resource amount available for the underground space is subjected to superposition calculation, and the resource amount of the urban underground space has different purposes under the conditions of different depths. This needs to be taken into account when the available resource amount is analyzed in superposition. 4 layers of 0m-10m, 10m-30m, 30m-50m and 50m-100m are taken as calculation units to carry out superposition analysis and calculate the available resource amount. Taking a superposition calculation method of available resource amount of 0m-10m underground space as an example: the calculation method is that the used resource amount of 10m and the used part of 0m-10m of underground space resource with the depth exceeding 10m are deducted from the total resource amount of 0m-10m, and the rest parts are similar. Through the steps, a spatial resource utilization distribution diagram of 4 layers of computing units of 0m-10m, 10m-30m, 30m-50m and 50m-100m is obtained, and finally, an underground spatial resource quantity analysis result is obtained.
The superposition calculation method of the available resource amount of the 0m-10m underground space is that the total resource amount of 0m-10m deducts the used resource amount with the influence depth of 10m and the used part of 0m-10m of the underground space resource with the influence depth exceeding 10m, and the calculation formula is as follows:
K0-10＝(S10+S30+S50+S100) 10 (two and two)
V(0-10m available resources)＝V(0-10m Total resource amount)-V(0-10m amount of used resources)-K0-10(Ershisanx type)
In the formula: k0-10Is a part of 0m to 10m, S, with the depth of more than 10m, where underground space resources are utilized10Is the area of the land with the influence of the depth of 10m, S30Is the area of the land with an influence of a depth of 30m, S50Is the area of the land, S, affecting the depth of 50m100Is the area of the land with an influence of the depth of 100m, V(0-10m available resources)Is 0m-10m of underground space available resource quantity, V(0-10m Total resource amount)Is the total resource amount of the underground space of 0m-10m, V(0-10m amount of used resources)Is the amount of the utilized resources of the underground space of 0m to 10 m.
The method for calculating the superposition of the available resource amount of the underground space of 10m-30m comprises the following steps: the 10m-30m total resource amount deducts the part of 10m-30m used by underground space resources with the influence depth exceeding 10m, and the calculation formula is as follows:
K10-30＝(S30+S50+S100) 20 (two and four)
V(10m-30m available resources)＝V(10m-30m Total resource amount)-K10-30(two Shi Wu)
In the formula: k10-30A 10m-30m part, V, of underground space resources having a depth of more than 10m(10m-30m available resources)Is 10m-30m of underground space available resource quantity, V(10m-30m Total resource amount)Is the total resource amount of the underground space of 10m-30 m.
The method for calculating the superposition of the available resource quantity of the underground space of 30m-50m comprises the following steps: and deducting the utilized part of 30m-50m of underground space resources with the influence depth exceeding 30m by the total resource amount of 30m-50m, and calculating according to the following formula:
K30-50＝(S50+S100) 20 (two shilelu)
V(30m-50m available resources)＝V(total resource amount of 30m-50 m)-K30-50(two ten seven)
In the formula: k30-50A 30m-50m part, V, of underground space resources with a depth of more than 30m(30m-50m available resources)Is 30m-50m of underground space available resource quantity, V(total resource amount of 30m-50 m)Is the total resource amount of the underground space of 30m-50 m.
The superposition calculation method of the available resource amount of the underground space with the size of 50m-100m comprises the following steps: the total resource amount of 50m-100m is deducted from the part of 50m-100m used by underground space resources with the influence depth exceeding 50m, and the calculation formula is as follows:
K50-100＝S10050 (eight type)
V(50-100m available resources)＝V(50-100m Total resource amount)-K50-100(Er Ji Jiu)
In the formula: k50-100A 50m-100m part, V, of underground space resources with a depth of more than 50m(50-100m available resources)Is 50m-100m of underground space available resource quantity, V(50-100m Total resource amount)Is the total resource amount of the underground space of 50m-100 m.
According to the method for analyzing the underground space resource quantity, the satellite image of the area to be analyzed is obtained; preprocessing the satellite image to obtain a processed satellite image; extracting buildings from the processed satellite images by adopting morphological building indexes and deep learning to obtain building extraction results; shadow extraction is carried out on the processed satellite images according to the characteristic indexes to obtain shadow extraction results; analyzing according to the sun azimuth line and the shadow extraction result to obtain the shadow length of each shadow in the shadow extraction result; according to the building extraction result and the shadow length of each shadow in the shadow extraction result, utilizing a building height inversion model to perform building height inversion work, and determining the building height of each building in the building extraction result; the building height of each building is converted into the building influence depth of each building, the utilization condition of the underground space resource quantity of each layer is analyzed by utilizing the building influence depth and the area of each building, the analysis result of the underground space resource quantity is obtained, and the problems that the utilized resource quantity is counted on the spot and the investigation cost is high are solved.
Building extraction is carried out on the area to be analyzed by using the high-score No. 2 remote sensing image and adopting a method combining building indexes and deep learning; by utilizing the high-resolution image fourth band threshold segmentation method, most of building shadows can be simply and effectively extracted, and the method is very suitable for actual work production; the height inversion models under different conditions are carefully analyzed in the aspect of building height inversion, and the height inversion work can be specifically and accurately completed; the resource amount of the urban underground space has different purposes under different depth conditions, and layered analysis is carried out when the resource amount available in the underground space is superposed and analyzed in consideration of the resource amount. And finally, combining the building extraction result with the building image depth result to complete the calculation of the underground space resource amount of the area to be analyzed.
It should be understood that, although the steps in the flowchart of fig. 1 are shown in order as indicated by the arrows, the steps are not necessarily performed in order as indicated by the arrows. The steps are not performed in the exact order shown and described, and may be performed in other orders, unless explicitly stated otherwise. Moreover, at least a portion of the steps in fig. 1 may include multiple sub-steps or multiple stages that are not necessarily performed at the same time, but may be performed at different times, and the order of performance of the sub-steps or stages is not necessarily sequential, but may be performed in turn or alternately with other steps or at least a portion of the sub-steps or stages of other steps.
The technical features of the above embodiments can be arbitrarily combined, and for the sake of brevity, all possible combinations of the technical features in the above embodiments are not described, but should be considered as the scope of the present specification as long as there is no contradiction between the combinations of the technical features.
The above-mentioned embodiments only express several embodiments of the present application, and the description thereof is more specific and detailed, but not construed as limiting the scope of the invention. It should be noted that, for a person skilled in the art, several variations and modifications can be made without departing from the concept of the present application, which falls within the scope of protection of the present application. Therefore, the protection scope of the present patent shall be subject to the appended claims.

Claims (9)

1. A method for analyzing the amount of underground spatial resources, the method comprising:
acquiring a satellite image of an area to be analyzed;
preprocessing the satellite image to obtain a processed satellite image;
extracting buildings from the processed satellite images by adopting morphological building indexes and deep learning to obtain building extraction results;
shadow extraction is carried out on the processed satellite images according to the characteristic indexes to obtain shadow extraction results;
analyzing according to the sun azimuth line and the shadow extraction result to obtain the shadow length of each shadow in the shadow extraction result;
according to the building extraction result and the shadow length of each shadow in the shadow extraction result, utilizing a building height inversion model to perform building height inversion work, and determining the building height of each building in the building extraction result;
and converting the building height of each building into the building influence depth of each building, and carrying out hierarchical underground space resource utilization analysis by using the building influence depth and area of each building to obtain an underground space resource analysis result.
2. The method of claim 1, wherein the step of extracting the building from the processed satellite images by using the morphological building index and the deep learning to obtain the building extraction result comprises:
extracting buildings from the processed satellite images by adopting a morphological building index to obtain a first extraction result;
extracting buildings from the processed satellite images by adopting deep learning to obtain a second extraction result;
fusing the first extraction result and the second extraction result to obtain a fused extraction result;
and optimizing the roof shape of the fused extraction result to obtain a building extraction result.
3. The method according to claim 1, wherein the step of extracting the shadow from the processed satellite image according to each feature index to obtain the shadow extraction result comprises:
comparing the intensity of the shadow expression force of each characteristic index, and determining the characteristic index with the highest intensity of the shadow expression force as the characteristic index of shadow extraction;
adopting the characteristic index of the shadow extraction to carry out shadow extraction work, and obtaining a preliminary shadow extraction result;
and optimizing the preliminary shadow extraction result by adopting clustering processing to obtain a shadow extraction result.
4. The method according to claim 1, wherein the step of determining the building height of each building in the building extraction result by performing a building height inversion operation using a building height inversion model according to the building extraction result and the shadow length of each shadow in the shadow extraction result comprises:
determining a building to be subjected to building height inversion according to the building extraction result;
and based on the spatial position relation among the sun, the satellite and the building and the relation between the building height and the shadow of the building, performing the inversion work of the building height by using a building height inversion model, and determining the building height of each building in the building extraction result.
5. The method of claim 4, wherein when the orientation of the sun and the satellite is the same, the altitude of the sun is smaller than the altitude of the satellite, and the projection of the sun ray on the ground is perpendicular to the building direction, the height inversion formula of the building height inversion model is as follows:
H＝S·(tanα·tanβ)/(tanα-tanβ)
where α is the altitude of the satellite, β is the altitude of the sun, S is the shadow length of the building in the shadow extraction result, and H is the building height of the building.
6. The method of claim 4, wherein when the building is not vertical to the projection of the sun rays on the ground, and the projection of the sun rays on the ground and the projection of the satellite sight line on the ground are located on two sides of the building, the height inversion formula of the building height inversion model is as follows:
H＝S·tanβ/sinγ
wherein gamma is the included angle between the building trend and the solar ray projection direction.
7. The method of claim 4, wherein when the building is not vertical to the projection of the sun rays on the ground and the projection of the sun rays and the projection of the satellite sight line are on the same side of the building, the height inversion formula of the building height inversion model is as follows:
H＝S·sinα·sinβ/(sinα·cosβ·sinγ-cosα·sinβ·sinδ)
wherein, delta is an included angle between the trend of the building and the projection direction of the satellite sight.
8. The method of claim 1, wherein the step of converting the building height of each building into a building impact depth of each building, and using the building impact depth and area of each building to perform a hierarchical geospatial resource usage analysis to obtain a geospatial resource usage analysis result comprises:
converting the building height of each building into the building influence depth of each building according to the conversion relation table of the building height and the building influence depth;
and analyzing the utilization condition of the layered underground space resource quantity according to the building influence depth and area of each building to obtain an underground space resource quantity analysis result.
9. The method according to claim 8, wherein the step of performing a hierarchical underground space resource utilization analysis according to the building influence depth and area of each building to obtain an underground space resource analysis result comprises:
according to the building influence depth of each building, respectively taking the building influence depth as 0m-10m, 10m-30m, 30m-50m and 50m-100m as calculation units;
and performing superposition analysis on the utilization conditions of the underground space resources of each computing unit to obtain an underground space resource analysis result.
CN202110411708.6A 2021-04-16 2021-04-16 Underground space resource quantity analysis method Active CN113203399B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110411708.6A CN113203399B (en) 2021-04-16 2021-04-16 Underground space resource quantity analysis method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110411708.6A CN113203399B (en) 2021-04-16 2021-04-16 Underground space resource quantity analysis method

Publications (2)

Publication Number Publication Date
CN113203399A true CN113203399A (en) 2021-08-03
CN113203399B CN113203399B (en) 2022-02-15

Family

ID=77027212

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110411708.6A Active CN113203399B (en) 2021-04-16 2021-04-16 Underground space resource quantity analysis method

Country Status (1)

Country Link
CN (1) CN113203399B (en)

Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2007033157A (en) * 2005-07-25 2007-02-08 Ntt Data Corp Image analysis system, image analysis method, and program
CN102855494A (en) * 2012-08-20 2013-01-02 中国测绘科学研究院 Method and device for extracting water body of satellite remote sensing image
CN106650689A (en) * 2016-12-30 2017-05-10 厦门理工学院 Coastal city time sequence land utilization information extracting method
CN106886584A (en) * 2017-02-14 2017-06-23 北京市地质调查研究院 Based on the underground space development present situation evaluation method on the basis of the various geodatas in city
CN107516100A (en) * 2017-08-31 2017-12-26 北京航天绘景科技有限公司 A kind of image building extracting method based on elevation morphology building index
CN107679441A (en) * 2017-02-14 2018-02-09 郑州大学 Method based on multi-temporal remote sensing image shadow extraction City Building height
CN110287275A (en) * 2019-05-19 2019-09-27 中国地质调查局西安地质调查中心 A kind of Underground Space Resource information processing system and method based on negative inventory method
CN111476896A (en) * 2020-03-24 2020-07-31 清华大学 Panoramic three-dimensional modeling method for urban overground and underground buildings and personnel
CN111721267A (en) * 2020-07-22 2020-09-29 河南大学 Method for predicting building height by using satellite image
CN112149594A (en) * 2020-09-29 2020-12-29 同济大学 Urban construction assessment method based on deep learning and high-resolution satellite images

Patent Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2007033157A (en) * 2005-07-25 2007-02-08 Ntt Data Corp Image analysis system, image analysis method, and program
CN102855494A (en) * 2012-08-20 2013-01-02 中国测绘科学研究院 Method and device for extracting water body of satellite remote sensing image
CN106650689A (en) * 2016-12-30 2017-05-10 厦门理工学院 Coastal city time sequence land utilization information extracting method
CN106886584A (en) * 2017-02-14 2017-06-23 北京市地质调查研究院 Based on the underground space development present situation evaluation method on the basis of the various geodatas in city
CN107679441A (en) * 2017-02-14 2018-02-09 郑州大学 Method based on multi-temporal remote sensing image shadow extraction City Building height
CN107516100A (en) * 2017-08-31 2017-12-26 北京航天绘景科技有限公司 A kind of image building extracting method based on elevation morphology building index
CN110287275A (en) * 2019-05-19 2019-09-27 中国地质调查局西安地质调查中心 A kind of Underground Space Resource information processing system and method based on negative inventory method
CN111476896A (en) * 2020-03-24 2020-07-31 清华大学 Panoramic three-dimensional modeling method for urban overground and underground buildings and personnel
CN111721267A (en) * 2020-07-22 2020-09-29 河南大学 Method for predicting building height by using satellite image
CN112149594A (en) * 2020-09-29 2020-12-29 同济大学 Urban construction assessment method based on deep learning and high-resolution satellite images

Also Published As

Publication number Publication date
CN113203399B (en) 2022-02-15

Similar Documents

Publication Publication Date Title
CN107067415B (en) A kind of object localization method based on images match
CN107909607B (en) A kind of year regional vegetation coverage computational methods
CN103971115B (en) Automatic extraction method for newly-increased construction land image spots based on NDVI and PanTex index
Wang et al. Segmentation of high spatial resolution remote sensing imagery based on hard-boundary constraint and two-stage merging
CN112381013B (en) Urban vegetation inversion method and system based on high-resolution remote sensing image
CN108053408B (en) High-automation land use updating method based on remote sensing satellite image
CN108195736B (en) Method for extracting vegetation canopy clearance rate through three-dimensional laser point cloud
Zhu et al. Robust registration of aerial images and LiDAR data using spatial constraints and Gabor structural features
CN110263716B (en) Remote sensing image super-resolution land cover mapping method based on street view image
CN110992366A (en) Image semantic segmentation method and device and storage medium
CN116030352B (en) Long-time-sequence land utilization classification method integrating multi-scale segmentation and super-pixel segmentation
Berveglieri et al. Identification of successional stages and cover changes of tropical forest based on digital surface model analysis
CN115512247A (en) Regional building damage grade assessment method based on image multi-parameter extraction
Anders et al. Rule set transferability for object-based feature extraction: An example for cirque mapping
CN113203399B (en) Underground space resource quantity analysis method
CN114842356B (en) High-resolution earth surface type sample automatic generation method, system and equipment
Ding et al. Study on building extraction from high-resolution images using Mbi
CN110929739A (en) Automatic impervious surface range remote sensing iterative extraction method
Trekin et al. Deep neural networks for determining the parameters of buildings from single-shot satellite imagery
CN107220615B (en) Urban impervious surface information extraction method fusing interest point big data
Gan et al. A random forest based method for urban object classification using lidar data and aerial imagery
Panigrahi et al. Image pan-sharpening and sub-pixel classification enabled building detection in strategically challenged forest neighborhood environment
Yu et al. An integrative object-based image analysis workflow for UAV images
Li Urban Remote Sensing Using Ground‐Based Street View Images
Niederheiser et al. Mapping alpine vegetation location properties by dense matching

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant