CN104574346A - Optical remote sensing image decomposition algorithm - Google Patents

Optical remote sensing image decomposition algorithm Download PDF

Info

Publication number
CN104574346A
CN104574346A CN201310502774.XA CN201310502774A CN104574346A CN 104574346 A CN104574346 A CN 104574346A CN 201310502774 A CN201310502774 A CN 201310502774A CN 104574346 A CN104574346 A CN 104574346A
Authority
CN
China
Prior art keywords
image
wavelet
decomposition
remote sensing
low
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.)
Pending
Application number
CN201310502774.XA
Other languages
Chinese (zh)
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.)
Beijing Research Institute of Uranium Geology
Original Assignee
Beijing Research Institute of Uranium Geology
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 Beijing Research Institute of Uranium Geology filed Critical Beijing Research Institute of Uranium Geology
Priority to CN201310502774.XA priority Critical patent/CN104574346A/en
Publication of CN104574346A publication Critical patent/CN104574346A/en
Pending legal-status Critical Current

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/40Analysis of texture
    • G06T7/41Analysis of texture based on statistical description of texture
    • G06T7/48Analysis of texture based on statistical description of texture using fractals
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/40Analysis of texture
    • G06T7/41Analysis of texture based on statistical description of texture
    • G06T7/42Analysis of texture based on statistical description of texture using transform domain methods
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20048Transform domain processing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20048Transform domain processing
    • G06T2207/20064Wavelet transform [DWT]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30181Earth observation

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Probability & Statistics with Applications (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Image Analysis (AREA)

Abstract

The invention belongs to the field of remote sensing image processing and particularly relates to an optical remote sensing image decomposition algorithm. The optical remote sensing image decomposition algorithm comprises the following steps: step 1, reading the stored original image data; step 2, performing logarithmic transformation; taking logarithm of F(x, y), and F (x, y)=1g(F(x, y)); step 3, performing wavelet decomposition; carrying out discrete wavelet transformation, wherein the wavelet basis is haar wavelet; step 4, performing inverse transformation; using zero to substitute {Hp, t} for wavelet inverse transformation with low-frequency subgraphs Lt to obtain low-frequency image Li in the image f spatial domain; step 5, conducting fractal computation; step 6, performing cyclic judgment; if Di is less than or equal to Dd, actuating step 7, otherwise, adding 1 to i, and then actuating the step 3 (wavelet decomposition); step 7, performing optimal dimension judgment; step 8, realizing image generation. The optical remote sensing image decomposition algorithm has the effect of solving the problem that the traditional image decomposition algorithm cannot accurately judge the optimal decomposition dimension.

Description

A kind of remote sensing image decomposition algorithm
Technical field
The invention belongs to field of remote sensing image processing, be specifically related to a kind of remote sensing image decomposition algorithm.
Background technology
At present, at remote sensing fields in order to obtain the feature of terrestrial materials component space distribution more accurately, the method for topographic correction is often adopted to eliminate the impact of landform.For this reason, wavelet multiresolution analysis method is utilized to become a kind of generally acknowledged thought to study different scale clarification of objective on remote sensing images, the method realized is the subgraph utilizing wavelet function original image to be decomposed into different scale, then utilizes subgraph to carry out target classification and the identification of different scale feature.But, how to determine that the rank of decomposing is a difficult problem.Existing picture breakdown algorithm utilizes different stage decomposition result usually, passes through the test compared with actual characteristic and determines.This method is mainly experimental, not only time-consuming, and lacks theoretical foundation, is difficult to promote.
Summary of the invention
The object of the invention is for prior art defect, a kind of remote sensing image decomposition algorithm is provided.
The present invention is achieved in that a kind of remote sensing image decomposition algorithm, comprises the steps:
Step one: digital independent
Read the raw image data stored, the image read out is gray-scale value, and the data that namely this step obtains are three-dimensional array, wherein bidimensional is the transverse and longitudinal coordinate of image, and the third dimension is the gray-scale value that coordinate is corresponding, the result (x that this step obtains, y, F(x, y)) represent, wherein x, y is respectively the transverse and longitudinal coordinate of image, F(x, y) be point (x, y) gray-scale value
Step 2: log-transformation
To F(x, y) take the logarithm, i.e. f (x, y)=lg (F (x, y)),
The result that this step obtains represents with (x, y, f (x, y)), also represents the set of the image that f (x, y) is formed in follow-up expression with f,
Step 3: wavelet decomposition
Carry out wavelet transform to the gray-scale value Mallat algorithm that step 2 obtains, wavelet basis is Haar wavelet transform,
This step is decomposed f, obtains wavelet conversion coefficient { L t, H p,t, L trepresent the low frequency subgraph picture of yardstick t hypograph f, H p,trepresent the high frequency subimage in p direction under yardstick t, p=1 here, 2,3, p=1 represents horizontal direction, and p=2 represents vertical direction, and p=3 represents angular direction,
Step 4: inverse transformation
{ H is replaced with zero p,t, carry out wavelet inverse transformation with low frequency subgraph Lt, obtain the low-frequency image L of image f spatial domain i, i=t, L when first time calculates irepresent the low-frequency image of the spatial domain under i yardstick, in successive iterations calculates, the value of i is determined according to subsequent step, by L in subsequent step iimage judges image L as waiting 0,
Step 5: fractal calculation
Low-frequency image L is read in by the Fraclab module of Matlab 0and DEM, and calculate with box meter dimension instrument, under ensureing that cross-correlation coefficient equals the prerequisite of 1, getting continuous 5 box meter dimensions that matching maximum error and match point span are obtained than regression equation calculation time minimum is calculated box meter dimension value, obtains low-frequency image L respectively 0with the box meter dimension D of DEM iand D d,
Step 6: cycle criterion
If D i≤ D d, perform step 7, otherwise make i value add 1, then perform step 3 wavelet decomposition,
Step 7: best scale judges
Judge D i-1-D tabsolute value and D i-D dthe size of absolute value, if Abs (D i-Dd)≤Abs (D i-1-D d), best decomposition scale is i, otherwise best scale is i-1, and best scale assignments is to I best,
Described Abs represents absolute value,
Step 8: Computer image genration
{ H is replaced with zero p,t, t=1,2 ..., I best, p=1,2,3, with yardstick I bestunder low frequency subgraph carry out wavelet inverse transformation, then carry out antilogarithm change and obtain low-frequency image---the landform subgraph of spatial domain,
With high frequency subgraph set { H p,t, t=1,2 ..., I best, p=1,2,3 and zero carry out wavelet inverse transformation, then carry out high frequency imaging---the lithology subgraph that antilogarithm change obtains spatial domain,
Export and represent the high frequency subgraph of lithology and represent the low frequency subgraph of landform.
A kind of remote sensing image decomposition algorithm as above, wherein, described DEM is the earth's surface elevation map picture being defined in the same space territory with image f, three-dimensional array equally, (x, y, D(x can be used, y)) represent, wherein x, y are respectively the transverse and longitudinal coordinate of image, D(x, y) be the height value of point (x, y).
Effect of the present invention is used to be: first the present invention utilizes the computing of taking the logarithm, and the information that remote sensing images comprise is converted to additive operation result by multiplicative operation result; Then wavelet function is utilized to decompose image.In order to determine to decompose rank, by calculating the method for the fractal dimension whether equal (within the scope of certain error) of DEM and low frequency subgraph, judge optimum wavelet transform dimension, be finally the lithology subgraph of the low frequency landform subgraph of spatial domain and the high frequency of spatial domain picture breakdown, overcome in traditional picture breakdown algorithm without the problem accurately differentiating optimal Decomposition yardstick.Therefore, by the method that fractal dimension calculates and wavelet decomposition combines, overcome and need by experience to judge the shortcoming of wavelet decomposition scales and topographic correction DeGrain, greatly can improve the effect of picture breakdown, have important meaning and practical value for remote sensing Objects recognition and Remote Sensing Data Fusion Algorithm.The thought of this decomposition algorithm can also be generalized among the decomposition application of other optical images.
Accompanying drawing explanation
Fig. 1 is the basic flow sheet of the inventive method;
Fig. 2 is DS mining area TM image wavelet decomposition result.
Embodiment
A kind of remote sensing image decomposition algorithm, comprises the steps:
Step one: digital independent
Read the raw image data stored, the image read out is gray-scale value, and the data that namely this step obtains are three-dimensional array, and wherein bidimensional is the transverse and longitudinal coordinate of image, and the third dimension is the gray-scale value that coordinate is corresponding.Result that this step obtains is with (x, y, F(x, y)) represent, wherein x, y are respectively the transverse and longitudinal coordinate of image, F(x, y) be the gray-scale value of point (x, y).
Step 2: log-transformation
To F(x, y) take the logarithm, i.e. f (x, y)=lg (F (x, y)).
This step will affect two factors---atural object and the landform of remote sensing images, be converted to additive operation by multiplicative computing.
The result that this step obtains represents with (x, y, f (x, y)), also represents the set of the image that f (x, y) is formed in follow-up expression with f.
Step 3: wavelet decomposition
Carry out wavelet transform to the gray-scale value Mallat algorithm that step 2 obtains, wavelet basis is Haar wavelet transform.
This step is decomposed f, obtains wavelet conversion coefficient { L t, H p,t, L trepresent the low frequency subgraph picture of yardstick t hypograph f, H p,trepresent the high frequency subimage in p direction under yardstick t.Here p=1,2,3, p=1 represents horizontal direction, and p=2 represents vertical direction, and p=3 represents angular direction.Result of calculation as shown in Figure 2.
Step 4: inverse transformation
{ H is replaced with zero p,t, carry out wavelet inverse transformation with low frequency subgraph Lt, obtain the low-frequency image L of image f spatial domain i.I=t, L when first time calculates irepresent the low-frequency image of the spatial domain under i yardstick, in successive iterations calculates, the value of i is determined according to subsequent step.The application is by L iimage judges image L as waiting 0.
Step 5: fractal calculation
Low-frequency image L is read in by the Fraclab module of Matlab 0and DEM, and calculate with box meter dimension instrument, under ensureing that cross-correlation coefficient equals the prerequisite of 1, get at continuous 5 and make matching maximum error and match point span be calculated box meter dimension value than the box meter dimension that regression equation calculation time minimum obtains.Obtain low-frequency image L respectively 0with the box meter dimension D of DEM iand D d.
Described DEM is the earth's surface elevation map picture being defined in the same space territory with image f, is three-dimensional array equally, can use (x, y, D(x, y)) represent, wherein x, y are respectively the transverse and longitudinal coordinate of image, D(x, y) be the height value of point (x, y).
Step 6: cycle criterion
If D i≤ D d, perform step 7, otherwise make i value add 1, then perform step 3 wavelet decomposition.
In certain calculates, when decomposing the 5th grade, L 0fractal dimension be 2.16599, be less than the fractal dimension 2.27822 of DEM, so stop at the 5th grade.
Step 7: best scale judges
Judge D i-1-D tabsolute value and D i-D dthe size of absolute value, if Abs (D i-Dd)≤Abs (D i-1-D d), best decomposition scale is i, otherwise best scale is i-1.Best scale assignments is to I best.
Described Abs represents absolute value.
Step 8: Computer image genration
{ H is replaced with zero p,t, t=1,2 ..., I best, p=1,2,3.With yardstick I bestunder low frequency subgraph carry out wavelet inverse transformation, then carry out antilogarithm change and obtain low-frequency image---the landform subgraph of spatial domain.
With high frequency subgraph set { H p,t, t=1,2 ..., I best, p=1,2,3 and zero carry out wavelet inverse transformation, then carry out high frequency imaging---the lithology subgraph that antilogarithm change obtains spatial domain.
Export and represent the high frequency subgraph of lithology and represent the low frequency subgraph of landform.

Claims (2)

1. a remote sensing image decomposition algorithm, is characterized in that, comprises the steps:
Step one: digital independent
Read the raw image data stored, the image read out is gray-scale value, and the data that namely this step obtains are three-dimensional array, wherein bidimensional is the transverse and longitudinal coordinate of image, and the third dimension is the gray-scale value that coordinate is corresponding, the result (x that this step obtains, y, F(x, y)) represent, wherein x, y is respectively the transverse and longitudinal coordinate of image, F(x, y) be point (x, y) gray-scale value
Step 2: log-transformation
To F(x, y) take the logarithm, i.e. f (x, y)=lg (F (x, y)),
The result that this step obtains represents with (x, y, f (x, y)), also represents the set of the image that f (x, y) is formed in follow-up expression with f,
Step 3: wavelet decomposition
Carry out wavelet transform to the gray-scale value Mallat algorithm that step 2 obtains, wavelet basis is Haar wavelet transform,
This step is decomposed f, obtains wavelet conversion coefficient { L t, H p,t, L trepresent the low frequency subgraph picture of yardstick t hypograph f, H p,trepresent the high frequency subimage in p direction under yardstick t, p=1 here, 2,3, p=1 represents horizontal direction, and p=2 represents vertical direction, and p=3 represents angular direction,
Step 4: inverse transformation
{ H is replaced with zero p,t, carry out wavelet inverse transformation with low frequency subgraph Lt, obtain the low-frequency image L of image f spatial domain i, i=t, L when first time calculates irepresent the low-frequency image of the spatial domain under i yardstick, in successive iterations calculates, the value of i is determined according to subsequent step, by L in subsequent step iimage judges image L as waiting 0,
Step 5: fractal calculation
Low-frequency image L is read in by the Fraclab module of Matlab 0and DEM, and calculate with box meter dimension instrument, under ensureing that cross-correlation coefficient equals the prerequisite of 1, getting continuous 5 box meter dimensions that matching maximum error and match point span are obtained than regression equation calculation time minimum is calculated box meter dimension value, obtains low-frequency image L respectively 0with the box meter dimension D of DEM iand D d,
Step 6: cycle criterion
If D i≤ D d, perform step 7, otherwise make i value add 1, then perform step 3 wavelet decomposition,
Step 7: best scale judges
Judge D i-1-D tabsolute value and D i-D dthe size of absolute value, if Abs (D i-Dd)≤Abs (D i-1-D d), best decomposition scale is i, otherwise best scale is i-1, and best scale assignments is to I best,
Described Abs represents absolute value,
Step 8: Computer image genration
{ H is replaced with zero p,t, t=1,2 ..., I best, p=1,2,3, with yardstick I bestunder low frequency subgraph carry out wavelet inverse transformation, then carry out antilogarithm change and obtain low-frequency image---the landform subgraph of spatial domain,
With high frequency subgraph set { H p,t, t=1,2 ..., I best, p=1,2,3 and zero carry out wavelet inverse transformation, then carry out high frequency imaging---the lithology subgraph that antilogarithm change obtains spatial domain,
Export and represent the high frequency subgraph of lithology and represent the low frequency subgraph of landform.
2. a kind of remote sensing image decomposition algorithm as claimed in claim 1, it is characterized in that: described DEM is the earth's surface elevation map picture being defined in the same space territory with image f, is three-dimensional array equally, can (x be used, y, D(x, y)) represent, wherein x, y is respectively the transverse and longitudinal coordinate of image, D(x, y) be the height value of point (x, y).
CN201310502774.XA 2013-10-23 2013-10-23 Optical remote sensing image decomposition algorithm Pending CN104574346A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310502774.XA CN104574346A (en) 2013-10-23 2013-10-23 Optical remote sensing image decomposition algorithm

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310502774.XA CN104574346A (en) 2013-10-23 2013-10-23 Optical remote sensing image decomposition algorithm

Publications (1)

Publication Number Publication Date
CN104574346A true CN104574346A (en) 2015-04-29

Family

ID=53090323

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310502774.XA Pending CN104574346A (en) 2013-10-23 2013-10-23 Optical remote sensing image decomposition algorithm

Country Status (1)

Country Link
CN (1) CN104574346A (en)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2019205725A1 (en) * 2018-04-28 2019-10-31 中国石油大学(华东) Lithology recognition method for volcanic rock based on fractal dimensions of electrical imaging logging
CN111060985A (en) * 2019-12-20 2020-04-24 核工业北京地质研究院 Construction method of sandstone-type uranium ore multi-element data sample set
CN116452439A (en) * 2023-03-29 2023-07-18 中国工程物理研究院计算机应用研究所 Noise reduction method and device for laser radar point cloud intensity image

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1892698A (en) * 2006-05-19 2007-01-10 中国科学院研究生院 Remote-sensing picture interpolation method based on small wave fractal
CN101599173A (en) * 2009-07-24 2009-12-09 武汉大学 The radiation correction method of a kind of satellite or aviation remote sensing image
US20110069905A1 (en) * 2008-05-30 2011-03-24 Ge Healthcare Bio-Sciences Corp System and method for detecting and eliminating one or more defocused or low contrast-to-noise ratio images
CN102521831A (en) * 2011-12-02 2012-06-27 南京信息工程大学 Robot vision image segmentation method based on multi-scale fractal dimension and neural network

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1892698A (en) * 2006-05-19 2007-01-10 中国科学院研究生院 Remote-sensing picture interpolation method based on small wave fractal
US20110069905A1 (en) * 2008-05-30 2011-03-24 Ge Healthcare Bio-Sciences Corp System and method for detecting and eliminating one or more defocused or low contrast-to-noise ratio images
CN101599173A (en) * 2009-07-24 2009-12-09 武汉大学 The radiation correction method of a kind of satellite or aviation remote sensing image
CN102521831A (en) * 2011-12-02 2012-06-27 南京信息工程大学 Robot vision image segmentation method based on multi-scale fractal dimension and neural network

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
潘蔚等: "基于遥感图像地形结构-岩性组分分解的岩类多重分形特征研究", 《地学前源(中国地质大学(北京);北京大学)》 *
韩子媛: "基于小波域高斯混合模型的绝缘子红外图像去噪方法研究", 《中国优秀硕士学位论文全文数据库 信息科技辑》 *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2019205725A1 (en) * 2018-04-28 2019-10-31 中国石油大学(华东) Lithology recognition method for volcanic rock based on fractal dimensions of electrical imaging logging
CN111060985A (en) * 2019-12-20 2020-04-24 核工业北京地质研究院 Construction method of sandstone-type uranium ore multi-element data sample set
CN111060985B (en) * 2019-12-20 2022-03-22 核工业北京地质研究院 Construction method of sandstone-type uranium ore multi-element data sample set
CN116452439A (en) * 2023-03-29 2023-07-18 中国工程物理研究院计算机应用研究所 Noise reduction method and device for laser radar point cloud intensity image

Similar Documents

Publication Publication Date Title
Tahmasebi et al. MS-CCSIM: accelerating pattern-based geostatistical simulation of categorical variables using a multi-scale search in Fourier space
CN103654789B (en) Fast magnetic resonance parametric formation method and system
CN111161199B (en) Space spectrum fusion hyperspectral image mixed pixel low-rank sparse decomposition method
CN110763187B (en) Robust ground subsidence monitoring method based on radar distributed targets
CN105006003A (en) Random projection fern based real-time target tracking algorithm
CN105354841A (en) Fast matching method and system for remote sensing images
Song et al. Extraction and reconstruction of curved surface buildings by contour clustering using airborne LiDAR data
CN103226825B (en) Based on the method for detecting change of remote sensing image of low-rank sparse model
CN103454636A (en) Differential interferometric phase estimation method based on multi-pixel covariance matrixes
CN103514599B (en) A kind of segmentation of the image optimum based on neighborhood total variation scale selection method
CN104574346A (en) Optical remote sensing image decomposition algorithm
CN114265062B (en) InSAR phase unwrapping method based on phase gradient estimation network
CN102831621B (en) Video significance processing method based on spectral analysis
CN103325104B (en) Based on the face image super-resolution reconstruction method of iteration sparse expression
CN106056551A (en) Local similarity sample learning-based sparse de-noising method
CN101908211B (en) High spectral image fusion method based on variational method
Kang et al. Image registration based on harris corner and mutual information
CN103236041B (en) A kind of image super-resolution rebuilding method based on contourlet transformation
CN104200472A (en) Non-local wavelet information based remote sensing image change detection method
CN101000651B (en) Method for recognizing multiple texture image
CN105389798A (en) Synthetic aperture radar (SAR) image segmentation method based on deconvolution network and mapping inference network
CN103018728B (en) Laser radar real-time imaging and building characteristic extracting method
CN104392209A (en) Evaluation model for image complexity of target and background
CN103606189A (en) Track base selection method facing non-rigid body three-dimensional reconstruction
CN103310456A (en) Multi-temporal/multi-mode remote sensing image registration method based on Gaussian-Hermite moments

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
WD01 Invention patent application deemed withdrawn after publication

Application publication date: 20150429

WD01 Invention patent application deemed withdrawn after publication