CN111681207A - Remote sensing image fusion quality evaluation method - Google Patents
Remote sensing image fusion quality evaluation method Download PDFInfo
- Publication number
- CN111681207A CN111681207A CN202010385384.9A CN202010385384A CN111681207A CN 111681207 A CN111681207 A CN 111681207A CN 202010385384 A CN202010385384 A CN 202010385384A CN 111681207 A CN111681207 A CN 111681207A
- Authority
- CN
- China
- Prior art keywords
- band
- image
- column vector
- sub
- equal
- 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
Links
- 230000004927 fusion Effects 0.000 title claims abstract description 71
- 238000000034 method Methods 0.000 title claims abstract description 52
- 238000013441 quality evaluation Methods 0.000 title claims abstract description 31
- 239000013598 vector Substances 0.000 claims abstract description 279
- 238000012360 testing method Methods 0.000 claims abstract description 81
- 238000012549 training Methods 0.000 claims abstract description 48
- 230000003595 spectral effect Effects 0.000 claims abstract description 23
- 238000005070 sampling Methods 0.000 claims description 30
- 239000011159 matrix material Substances 0.000 claims description 18
- 238000005520 cutting process Methods 0.000 claims description 15
- 238000003708 edge detection Methods 0.000 claims description 4
- 238000011156 evaluation Methods 0.000 abstract description 6
- 230000008447 perception Effects 0.000 abstract description 4
- 239000003550 marker Substances 0.000 description 3
- 238000010586 diagram Methods 0.000 description 2
- 230000006870 function Effects 0.000 description 2
- 238000007500 overflow downdraw method Methods 0.000 description 2
- GCKMFJBGXUYNAG-HLXURNFRSA-N Methyltestosterone Chemical compound C1CC2=CC(=O)CC[C@]2(C)[C@@H]2[C@@H]1[C@@H]1CC[C@](C)(O)[C@@]1(C)CC2 GCKMFJBGXUYNAG-HLXURNFRSA-N 0.000 description 1
- AUQMGYLQQPSCNH-UHFFFAOYSA-L NIR-2 dye Chemical compound [K+].[K+].C1=CC2=C(S([O-])(=O)=O)C=C(S([O-])(=O)=O)C=C2C(C2(C)C)=C1[N+](CC)=C2C=CC=CC=C1C(C)(C)C2=CC(C(O)=O)=CC=C2N1CCCCS([O-])(=O)=O AUQMGYLQQPSCNH-UHFFFAOYSA-L 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 229940085503 testred Drugs 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/13—Edge detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20081—Training; Learning
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30168—Image quality inspection
Landscapes
- Engineering & Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Quality & Reliability (AREA)
- Image Analysis (AREA)
Abstract
The invention discloses a remote sensing image fusion quality evaluation method, which comprises the steps of extracting LBP characteristic statistical histograms, edge histograms and spectral characteristics of subblocks in images of different wave bands as characteristic vectors, and constructing an original multivariate Gaussian model in a training stage; in the testing stage, constructing an original multivariate Gaussian model in the testing stage, and calculating a full-resolution quality evaluation predicted value according to two original multivariate Gaussian models constructed in training and testing; calculating the spectral similarity and the spatial similarity of all band images of the multispectral image referenced by the test remote sensing fusion image and all band images of the test remote sensing fusion image after downsampling operation to obtain a reduced resolution quality evaluation predicted value; and then an objective quality evaluation predicted value is obtained, and the obtained characteristic vector information can better reflect the quality change condition of the remote sensing fusion image, so that the correlation between the objective evaluation result and the subjective perception of human eyes is effectively improved.
Description
Technical Field
The invention relates to an image quality evaluation method, in particular to a remote sensing image fusion quality evaluation method.
Background
Panchromatic multispectral fusion is a research hotspot in the field of remote sensing image fusion at present. Multispectral images are generally of lower spatial resolution, while panchromatic images generally have higher spatial resolution but lack spectral information. The purpose of panchromatic multispectral fusion is to fuse the spatial information of a panchromatic image into a multispectral image to obtain the multispectral image with high spatial resolution. However, different fusion methods have different fusion effects on different panchromatic images and multispectral images, and how to evaluate the quality of fused images is a technical problem to be solved at present.
In the existing panchromatic multispectral fusion image quality evaluation field, the following two problems mainly exist: (1) the spatial resolution of the panchromatic multispectral fusion image is inconsistent with that of the multispectral image, so that the panchromatic multispectral fusion image cannot be evaluated by directly adopting the existing full-reference quality evaluation method; (2) because the distortion problem of the spectral domain is difficult to measure by a non-reference quality evaluation method, and the problem of the distortion of the spectral domain is too serious, which can cause difficulty in identifying ground features, how to obtain evaluation reference information from a panchromatic multispectral fusion image with full resolution and reduced resolution and how to obtain the characteristic capable of measuring the distortion of the spectral domain is a problem which needs to be solved for evaluating the quality of the panchromatic multispectral fusion image.
Disclosure of Invention
The invention aims to provide a method for evaluating the fusion quality of remote sensing images, which can effectively improve the correlation between objective evaluation results and subjective perception of human eyes.
The technical scheme adopted by the invention for solving the technical problems is as follows: a method for evaluating fusion quality of remote sensing images is characterized by comprising a training stage and a testing stage;
the specific steps of the training phase process are as follows:
step ① _1, selecting N original multispectral images with width W and height H, each original multispectral image containing B different waveband images, and recording the gamma waveband image of the u-th original multispectral image as
1 st wave band image of u original multispectral image
The 2 nd wave band image of the u nd original multispectral image is the wave band image of the red wave bandImage of a person
Is a wave band image of a green wave band, and a 3 rd wave band image of a u original multispectral image
The wave band image of the blue wave band, the 4 th wave band image of the u original multispectral image
The band image is a near infrared band; wherein N is a positive integer, N is more than 1, B is a positive integer, u is more than or equal to 1 and less than or equal to N, gamma is a positive integer, gamma is more than or equal to 1 and less than or equal to B, and gamma is simultaneously used as a waveband mark of the waveband image;
① _2, cutting a sub-block with size of 96 × 96 from the same coordinate position in all band images of the original multispectral image, marking the sub-blocks cut from the same coordinate position in all band images of the original multispectral image with the same block index number, marking the sub-blocks cut from the same coordinate position in all band images of the original multispectral image with band marks to indicate the source of the sub-blocks, repeating the cutting process in a random cutting mode, cutting M times, then cutting M sub-blocks from all band images of the original multispectral image, arranging the pixel values of all pixel points in each cut sub-block to form a column vector corresponding to the sub-block, forming a column vector set with M 'sub-blocks of the same cut band mark band, and marking the M' sub-block corresponding to gamma vector set as a column set
The intercepted wave band is marked as a column vector set formed by column vectors corresponding to M' sub-blocks with the gamma being 1
Defining as training red wave band column vector set, marking the intercepted wave band as column vector set formed from column vectors correspondent to M' sub-blocks whose gamma is 2
Defining as training green band column vector set, marking the intercepted band as column vector set formed from column vectors correspondent to M' sub-blocks whose gamma is 3
Defining a column vector set formed by column vectors corresponding to M' sub-blocks with the intercepted wave band marked as gamma being 4 as a training blue wave band column vector set
Defining as training near infrared wave band column vector set, where M is positive integer, M is greater than or equal to 1, M 'is N × M, M is positive integer, M is greater than or equal to 1 and less than or equal to M', M is used as block index number of sub-block in training phase,
to represent
The m-th column vector in (e), i.e. the column vector representing the sub-block with the band labeled gamma and the block index m,
are (96 × 96) × 1,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose denoted band is γ ═ 1 and whose block index is m,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose denoted band is γ ═ 2 and block index m,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose representation band is denoted by γ ═ 3 and whose block index is m,
to represent
The mth column vector in (1), that is, the column vector corresponding to the sub-block whose band is denoted by γ ═ 4 and whose block index is m;
① _3, according to the training red band column vector set
Training a green band column vector set
And training a blue band column vector set
The training sample set of luminance components is obtained by calculation and is marked as { Ym|1≤m≤M'}，
wherein ,YmHas the dimension of (96 × 96) × 1 and symbols
Is a rounding down operation symbol;
step ① _4, calculating { YmThe statistical histogram of LBP features of all pixel points in the sub-block corresponding to each column vector in the sub-block with the value of 1 ≦ M ≦ M' } is calculated according to the statistical histogram of the LBP features of all the pixel points in the sub-block corresponding to each column vector in the sub-block with the value of YmThe statistical histogram of LBP features of all pixel points in the sub-block corresponding to the mth column vector in the M 'is | < M ≦ 1 ≦ M' } is marked as Fm； wherein ,FmHas a dimension of 8 × 1;
step ① _5, extracting { Y ] by adopting Canny edge detection operatormThe M is more than or equal to 1 and less than or equal to M' } is the edge map of the sub-block corresponding to each column vector; then calculate { YmThe edge histograms of all pixel points in the edge graph of the subblock corresponding to each column vector in the |1 ≦ M ≦ M' } are calculated according to the edge histograms of all the pixel points in the subblock corresponding to each column vector in the column vectormThe edge histograms of all pixel points in the edge graph of the sub-block corresponding to the mth column vector in the |1 ≦ M ≦ M' } are marked as Gm； wherein ,GmHas a dimension of 2 × 1;
① _6, according to the training red band column vector set
Training a green band column vector set
Training blue band column vector set
And training a near infrared band column vector set
Calculating the spectral characteristics of all subblocks with the same block index number, and recording the spectral characteristics of all subblocks with the block index number m as Hm； wherein ,HmHas a dimension of 6 × 1;
step ① _7, according to { Y }mI1 is less than or equal to M is less than or equal to M '} and the LBP feature statistical histogram of all pixel points in the sub-block corresponding to each column vector in the { Y ≦ M' } is obtainedmThe edge histograms of all pixel points in the edge graph of the subblock corresponding to each column vector in the |1 ≦ M ≦ M' } and the spectral characteristics of the subblocks with the same block index number are obtained to obtain all the same block index numbersThe feature vector of the subblock of which the index number of all the blocks is m is marked as Vm，Vm＝[Fm T,Gm T,Hm T]T； wherein ,VmHas a dimension of 16 × 1 [ F ]m T,Gm T,Hm T]Is shown asm T、Gm T and Hm TAre connected to form a vector, Fm TIs FmTranspose of (G)m TIs GmTranspose of (H)m TIs HmIs transposed, [ F ]m T,Gm T,Hm T]TIs [ F ]m T,Gm T,Hm T]Transposing;
step ① _8, constructing original multivariate Gaussian models of all band images of all original multispectral images according to the total M' eigenvectors obtained in the step ① _7, and marking as the original multivariate Gaussian models
wherein ,
to represent
The mean value vector of (a) is,
to represent
The covariance matrix of (a);
the specific steps of the test phase process are as follows:
in step ② _1, any test remote sensing fused image with width W 'and height H' includes B different waveband images, and the gamma waveband image of the test remote sensing fused image is marked as { S(γ)(x ', y') }; testing remote sensing fused image referenceThe multispectral image with the width W '/2 and the height H'/2 also comprises B different waveband images, and the gamma waveband image of the multispectral image referenced by the test remote sensing fusion image is marked as { R(γ)(x ", y") }; wherein x 'is more than or equal to 1 and less than or equal to W', y 'is more than or equal to 1 and less than or equal to H', x is more than or equal to 1 and less than or equal to W '/2, y is more than or equal to 1 and less than or equal to H'/2, gamma is more than or equal to 1 and less than or equal to B, gamma is simultaneously used as a waveband mark of a waveband image, S is(γ)(x ', y') denotes { S }(γ)(x ', y') } pixel value of pixel point with coordinate position (x ', y'), R(γ)(x ", y") denotes { R }(γ)(x ", y") } the pixel value of the pixel point whose coordinate position is (x ", y");
step ② _2, calculating and testing the full resolution quality evaluation predicted value of the remote sensing fusion image, and recording the predicted value as QfullThe specific process is as follows:
② _2a, using a moving window with a window size of 96 × 96 and a window moving step length of 1 pixel, moving in all band images of the tested remote sensing fused image to intercept subblocks with the same coordinate position and size of 96 × 96, marking the subblocks intercepted in the same coordinate position in all band images of the tested remote sensing fused image with the same block index number, marking each subblock intercepted in all band images of the tested remote sensing fused image with a band mark to indicate the source of the subblock, intercepting M ' subblocks from all band images of the tested remote sensing fused image after moving for M ' times, arranging pixel values of all pixel points in each intercepted subblock in sequence to form a column vector corresponding to the subblock, forming a column vector set corresponding to the M ' subblocks of each band mark intercepted to form a column vector set, and recording the column set formed by the column vectors corresponding to the M ' subblocks marked as gamma, intercepted M ' subblocks marked as gamma
Marking the intercepted wave band as a column vector set formed by column vectors corresponding to M' sub-blocks with gamma being 1
Defining as a test red wave band column vector set, and marking the intercepted wave bands asA column vector set composed of column vectors corresponding to M' sub-blocks with gamma being 2
Defining a column vector set formed by column vectors corresponding to M' sub-blocks for testing a green band column vector set and marking the intercepted band as gamma being 3
Defining a column vector set formed by column vectors corresponding to M' sub-blocks with the intercepted wave band marked as gamma being 4 as a test blue wave band column vector set
Is defined as testing near-infrared band column vector set, wherein M 'is positive integer, 1 is more than M and less than or equal to (W' -96) × (H '-96), M' is positive integer, 1 is less than or equal to M and less than or equal to M ', M' is simultaneously used as block index number of subblocks in the testing stage process,
has a dimension of (96 × 96) × 1,
to represent
The m "column vector in (1), i.e. the column vector representing the sub-block with the band marked as gamma and the block index number m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted γ and the block index m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted γ ═ 2 and the block index m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted as γ ═ 3 and the block index m",
to represent
The m "column vector in (1), that is, the column vector corresponding to the sub-block whose band is denoted by γ ═ 4 and whose block index is m";
step ② _2b, obtaining the feature vectors of all sub-blocks with the same block index number in the same way according to the process from step ① _3 to step ① _7, and recording the feature vectors of all sub-blocks with block index number m' as Vtest,m”； wherein ,Vtest,m”Has a dimension of 16 × 1;
step ② _2c, constructing original multivariate Gaussian models of all wave band images of the remote sensing fused image to be tested in the same way according to the M' eigenvectors obtained in the step ② _2b and the process of the step ① _8, and marking as the original multivariate Gaussian models
wherein ,
to represent
The mean value vector of (a) is,
to represent
The covariance matrix of (a);
step ② _2d according to
And
calculating Qfull，
wherein ,
is composed of
The transpose of (a) is performed,
is composed of
The inverse of (1);
step ② _3, calculating and testing the reduced resolution quality evaluation prediction value of the remote sensing fusion image, and recording the value as QreduThe specific process is as follows:
step ② _3a, performing down-sampling operation on each waveband image of the tested remote sensing fusion image to obtain a corresponding down-sampling waveband image with width W '/2 and height H'/2, and performing down-sampling on { S }(γ)(x ', y') } the down-sampling waveband image with the width W '/2 and the height H'/2 obtained after the down-sampling operation is recorded as the down-sampling waveband image
Wherein x is more than or equal to 1 and less than or equal to W '/2, y is more than or equal to 1 and less than or equal to H'/2,
to represent
The pixel value of the pixel point with the middle coordinate position of (x', y ");
step ② _3b, calculating the spectral similarity of the down-sampling waveband image obtained by the down-sampling operation of all waveband images of the multi-spectral image referenced by the test remote sensing fusion image and all waveband images of the test remote sensing fusion image, and recording as Qspe，
wherein ,Rx”,y”A vector S formed by sequentially connecting pixel values of pixel points with coordinate positions (x, y') in all band images of the multispectral image for testing remote sensing fusion image referencex”,y”The pixel values of pixel points with coordinate positions (x, y') in the downsampled waveband images obtained by downsampling all waveband images representing the tested remote sensing fusion image are connected in sequence to form a vector Rx”,y” and Sx”,y”All dimensions of (A) are B × 1, symbol "<>"is an inner product operation symbol, acrcos () is an inverse cosine function, the symbol" | | | | luminance2"is the 2-norm sign of the matrix;
step ② _3c, calculating the spatial similarity of the down-sampling waveband images obtained by the down-sampling operation of all waveband images of the multispectral image referenced by the test remote sensing fusion image and all waveband images of the test remote sensing fusion image, and recording as Qspa，
wherein ,
represents Rx”,y”The average of the values of all the elements in (a),
denotes Sx”,y”The average of the values of all the elements in (a),
represents Rx”,y”The standard deviation of the values of all elements in (a),
denotes Sx”,y”The standard deviation of the values of all elements in (a),
denotes Sx”,y”The values of all elements in (1) and Rx”,y”Of all elements in (c), ω1 and ω2Is a control parameter;
step ② _3d, according to Qspe and QspaCalculating Qredu，Qredu＝(Qspe)η×(Qspa)1-ηWherein η is a weight coefficient;
step ② _4, according to Qfull and QreduCalculating and testing the objective quality evaluation predicted value of the remote sensing fusion image and recording the predicted value as Qtest，
Wherein λ is a weight coefficient.
In the step ① _6, HmThe acquisition process comprises the following steps:
step ① _6a, calculation
And
spectral feature of (1), denoted as d12，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
The standard deviation of the values of all elements in (a),
representation calculation
Standard deviation of the values of all elements in (a);
step ① _6b, calculate
And
spectral feature of (1), denoted as d13，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
Standard deviation of the values of all elements in (a);
step ① _6c, calculate
And
spectral feature of (1), denoted as d14，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
Standard deviation of the values of all elements in (a);
step ① _6d, calculate
And
spectral feature of (1), denoted as d23，
wherein ,
representation calculation
All values of (1) and
value of all elements inThe covariance of (a);
step ① _6e, calculate
And
spectral feature of (1), denoted as d24，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a);
step ① _6f, calculate
And
spectral feature of (1), denoted as d34，
wherein ,
representation calculation
All ofValue of element and
the covariance of the values of all elements in (a);
step ① _6g, according to d12、d13、d14、d23、d24 and d34Obtaining Hm，Hm＝[d12,d13,d14,d23,d24,d34]T； wherein ,HmHas a dimension of 6 × 1 [ d ]12,d13,d14,d23,d24,d34]Is shown as12、d13、d14、d23、d24 and d34Are connected to form a vector, [ d ]12,d13,d14,d23,d24,d34]TIs [ d ]12,d13,d14,d23,d24,d34]The transposing of (1).
In the step ① _8, the data is sent,
the acquisition process comprises the following steps:
step ① _8a, obtaining the feature matrix of all sub-blocks with the same block index number according to the total M' feature vectors obtained in step ① _7, and marking as V,
wherein the dimension of V is M' × 16, (V)1)TIs a V1Transpose of (V)1Feature vector representing subblocks with all block indices of 1, (V)2)TIs a V2Transpose of (V)2Feature vector representing subblocks with all chunk indices of 2, (V)m)TIs a VmTranspose of (V)M')TIs a VM'Transpose of (V)M'A feature vector representing a subblock having all block indexes M';
step ① _8b, calculating the mean vector of the column vectors in V, and recording as
wherein ,
is 1 × 16, V (m,1) represents the value of the element in row m, column 1 in V, V (m,2) represents the value of the element in row m, column 2 in V, V (m,16) represents the value of the element in row m, column 16 in V;
and calculate the covariance matrix of the column vectors in V, noted
wherein ,
dimension 16 × 16, cov (a, b) represents calculating the covariance of the values of all elements in vector a and all elements in vector b;
step ① _8c, will
As
Mean vector of (2) will
As
Is constructed to obtain a covariance matrix
Compared with the prior art, the invention has the advantages that:
the method of the invention considers the influence of the spectral domain distortion and the spatial domain distortion on the quality of the remote sensing fusion image, and constructs the original multivariate Gaussian model of all the wave band images of all the original multispectral images by extracting LBP feature statistical histograms, edge histograms and spectral features of all sub-blocks of different wave band images as feature vectors; in the testing stage, original multivariate Gaussian models of all band images of the tested remote sensing fused image are constructed, and the full-resolution quality evaluation predicted value of the tested remote sensing fused image is calculated according to two original multivariate Gaussian models constructed in training and testing; the method comprises the steps of calculating the spectral similarity and the spatial similarity of all band images of a multispectral image referenced by a test remote sensing fusion image and all band images of the test remote sensing fusion image after downsampling operation by downsampling each band image of the test remote sensing fusion image to obtain a reduced resolution quality evaluation predicted value of the test remote sensing fusion image; the objective quality evaluation predicted value is obtained by combining the full-resolution quality evaluation predicted value and the deresolved quality evaluation predicted value of the tested remote sensing fusion image, and the correlation between the objective evaluation result and the subjective perception of human eyes is effectively improved because the obtained characteristic vector information can better reflect the quality change condition of the remote sensing fusion image.
Drawings
Fig. 1 is a block diagram of the overall implementation of the method of the present invention.
Detailed Description
The invention is described in further detail below with reference to the accompanying examples.
The general implementation block diagram of the method for evaluating the fusion quality of the remote sensing images provided by the invention is shown in fig. 1, and the method comprises a training stage and a testing stage.
The specific steps of the training phase process are as follows:
step ① _1, selecting N original multispectral images with width W and height H, each original multispectral image containing B different waveband images, and recording the gamma waveband image of the u-th original multispectral image as
1 st wave band image of u original multispectral image
The wave band image of the red wave band and the 2 nd wave band image of the u original multispectral image
Is a wave band image of a green wave band, and a 3 rd wave band image of a u original multispectral image
The wave band image of the blue wave band, the 4 th wave band image of the u original multispectral image
The band image is a near infrared band; where N is a positive integer, N > 1, where N is 100, B is a positive integer, the number of bands of the multispectral image is B, u is a positive integer, u is greater than or equal to 1 and less than or equal to N, γ is a positive integer, γ is greater than or equal to 1 and less than or equal to B, and γ is simultaneously used as a band marker of the band image, that is, the band markers of the band images of the same band of all the original multispectral images are the same, for example: the band label of the 1 st band image of all the original multispectral images is 1.
Here, the multispectral images of IKONOS, QuickBird, Gaofen-1 and WorldView-4 satellites used contain 4 band images of red, green, blue and near infrared, i.e. B is 4; the multispectral images of the adopted WorldView-2 and WorldView-3 satellites contain 8 band images of red, green, blue, near infrared, coast band, yellow band, red edge band and near infrared 2 band, namely B is 8.
① _2, cutting a sub-block with size of 96 × 96 from the same coordinate position in all band images of each original multispectral image, marking the sub-blocks cut from the same coordinate position in all band images of the same original multispectral image with the same block index number, marking the sub-blocks cut from the same coordinate position in all band images of the same original multispectral image with the same band mark to indicate the source of the sub-blocks, repeating the cutting process in a random cutting mode, cutting a total of M multispectral sub-blocks from all band images of each original multispectral image after cutting M times, namely, N × M sub-blocks are shared by the band images of the same band mark of N different multispectral images, and N × B × M sub-blocks are marked by the sub-blocks of the same band image of N different multispectral images, then arranging the pixels of each cut sub-block in sequence to form a corresponding vector column set of the vector column, and marking the M sub-blocks as a vector column set of the same vector column
That is, column vectors corresponding to M 'sub-blocks cut from band images of the same band marker of all multispectral images form a column vector set, and a column vector set formed by column vectors corresponding to M' sub-blocks cut from band images of the same band marker γ is recorded as
The intercepted wave band is marked as a column vector set formed by column vectors corresponding to M' sub-blocks with the gamma being 1
Defining as training red wave band column vector set, marking the intercepted wave band as column vector set formed from column vectors correspondent to M' sub-blocks whose gamma is 2
Defining as training green band column vector set, marking the intercepted band as column vector set formed from column vectors correspondent to M' sub-blocks whose gamma is 3
Defining a column vector set formed by column vectors corresponding to M' sub-blocks with the intercepted wave band marked as gamma being 4 as a training blue wave band column vector set
Defining as a training near-infrared band column vector set, where M is a positive integer, M is greater than or equal to 1, where M is 100 in this embodiment, arranging pixel values of all pixel points in each intercepted subblock in sequence to form a column vector corresponding to the subblock, where the sequence may be in a conventional row-by-row ordering manner, or in a conventional column-by-column ordering manner, or in a conventional zigzag scanning manner, or in a self-defined ordering manner, where M 'is N × M, M is a positive integer, M is greater than or equal to 1 and less than or equal to M', and M is simultaneously used as a block index number of the subblock in the training stage process,
to represent
The m-th column vector in (e), i.e. the column vector representing the sub-block with the band labeled gamma and the block index m,
are (96 × 96) × 1,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose denoted band is γ ═ 1 and whose block index is m,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose denoted band is γ ═ 2 and block index m,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose representation band is denoted by γ ═ 3 and whose block index is m,
to represent
I.e. the column vector representing the sub-block with the band denoted y 4 and the block index m.
① _3, according to the training red band column vector set
Training a green band column vector set
And training a blue band column vector set
The training sample set of luminance components is obtained by calculation and is marked as { Ym|1≤m≤M'}，
wherein ,YmHas the dimension of (96 × 96) × 1 and symbols
To round the operator down.
Step ① _4, calculate { Y Using Prior ArtmThe statistical histogram of LBP features of all pixel points in the sub-block corresponding to each column vector in the sub-block with the value of 1 ≦ M ≦ M' } is calculated according to the statistical histogram of the LBP features of all the pixel points in the sub-block corresponding to each column vector in the sub-block with the value of YmThe statistical histogram of LBP features of all pixel points in the sub-block corresponding to the mth column vector in the M 'is | < M ≦ 1 ≦ M' } is marked as Fm； wherein ,FmHas a dimension of 8 × 1.
Step ① _5, extracting { Y ] by adopting Canny edge detection operatormThe M is more than or equal to 1 and less than or equal to M' } is the edge map of the sub-block corresponding to each column vector; then using the prior art to calculate { Y }mThe edge histograms of all pixel points in the edge graph of the subblock corresponding to each column vector in the |1 ≦ M ≦ M' } are calculated according to the edge histograms of all the pixel points in the subblock corresponding to each column vector in the column vectormThe edge histograms of all pixel points in the edge graph of the sub-block corresponding to the mth column vector in the |1 ≦ M ≦ M' } are marked as Gm； wherein ,GmHas a dimension of 2 × 1.
① _6, according to the training red band column vector set
Training a green band column vector set
Training blue band column vector set
And training a near infrared band column vector set
Calculating the spectral characteristics of the subblocks with the same block index number in all the waveband images, and recording the spectral characteristics of the subblocks with the block index number m in all the waveband images as Hm； wherein ,HmHas a dimension of 6 × 1.
In the present embodiment, in step ① _6, HmThe acquisition process comprises the following steps:
step ① _6a, calculation
And
spectral feature of (1), denoted as d12，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
Of all elements inThe mean value of the values is determined,
representation calculation
The standard deviation of the values of all elements in (a),
representation calculation
Standard deviation of the values of all elements in (a).
Step ① _6b, calculate
And
spectral feature of (1), denoted as d13，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
Of the values of all elements inThe average value of the average value is calculated,
representation calculation
Standard deviation of the values of all elements in (a).
Step ① _6c, calculate
And
spectral feature of (1), denoted as d14，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
Of the values of all elements inStandard deviation.
Step ① _6d, calculate
And
spectral feature of (1), denoted as d23，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a).
Step ① _6e, calculate
And
spectral feature of (1), denoted as d24，
wherein ,
representation calculation
All elements inThe value of the element and
the covariance of the values of all elements in (a).
Step ① _6f, calculate
And
spectral feature of (1), denoted as d34，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a).
Step ① _6g, according to d12、d13、d14、d23、d24 and d34Obtaining Hm，Hm＝[d12,d13,d14,d23,d24,d34]T； wherein ,HmHas a dimension of 6 × 1 [ d ]12,d13,d14,d23,d24,d34]Is shown as12、d13、d14、d23、d24 and d34Are connected to form a vector, [ d ]12,d13,d14,d23,d24,d34]TIs [ d ]12,d13,d14,d23,d24,d34]The transposing of (1).
Step ① _7, according to { Y }mI1 is less than or equal to M is less than or equal to M '} and the LBP feature statistical histogram of all pixel points in the sub-block corresponding to each column vector in the { Y ≦ M' } is obtainedmThe edge histograms of all pixel points in the edge graph of the subblock corresponding to each column vector in the [ 1 ] M 'are larger than or equal to the [ M' ], the spectral characteristics of the subblocks with the same block index number in all the wave band images are obtained, the characteristic vectors of the subblocks with the same block index number in all the wave band images are recorded as V, and the characteristic vector of the subblock with the block index number M in all the wave band images is recorded as Vm，Vm＝[Fm T,Gm T,Hm T]T； wherein ,VmHas a dimension of 16 × 1 [ F ]m T,Gm T,Hm T]Is shown asm T、Gm T and Hm TAre connected to form a vector, Fm TIs FmTranspose of (G)m TIs GmTranspose of (H)m TIs HmIs transposed, [ F ]m T,Gm T,Hm T]TIs [ F ]m T,Gm T,Hm T]The transposing of (1).
Step ① _8, constructing original Multivariate Gaussian (MVG) models of all band images of all original multispectral images according to the total M' eigenvectors obtained in step ① _7, and marking as MVG models
wherein ,
to represent
The mean value vector of (a) is,
to represent
The covariance matrix of (2).
In this embodiment, in step ① — 8,
the acquisition process comprises the following steps:
step ① _8a, obtaining the feature matrix of the sub-blocks with the same block index number in all the waveband images according to the M' feature vectors obtained in step ① _7, and marking as V,
wherein the dimension of V is M' × 16, (V)1)TIs a V1Transpose of (V)1A feature vector representing a subblock with a block index of 1 in all band images, (V)2)TIs a V2Transpose of (V)2Feature vector (V) representing subblocks with block index 2 in all band imagesm)TIs a VmTranspose of (V)M')TIs a VM'Transpose of (V)M'The feature vector representing the subblock with the block index M' in all band images.
Step ① _8b, calculating the mean vector of the column vectors in V, and recording as
wherein ,
is 1 × 16, V (m,1) represents the value of the element in row m, column 1 in V, V (m,2) represents the value of the element in row m, column 2 in V, and V (m,16) represents the value of the element in row m, column 16 in V.
And calculates the covariance of the column vectors in VDifference matrix, is recorded as
wherein ,
dimension (b) 16 × 16, cov (a, b) represents calculating the covariance of the values of all elements in vector a and the values of all elements in vector b.
Step ① _8c, will
As
Mean vector of (2) will
As
Is constructed to obtain a covariance matrix
The specific steps of the test phase process are as follows:
step ② _1, selecting any test remote sensing fusion image packet with width W' and height HB different waveband images are contained, and the gamma waveband image of the tested remote sensing fused image is marked as { S(γ)(x ', y') }; the multispectral image with the width W '/2 and the height H'/2 referred by the test remote sensing fusion image also comprises B different waveband images, and the gamma waveband image of the multispectral image referred by the test remote sensing fusion image is marked as { R(γ)(x ", y") }; wherein, x ' is more than or equal to 1 and less than or equal to W ', y ' is more than or equal to 1 and less than or equal to H ', x is more than or equal to 1 and less than or equal to W '/2, y ' is more than or equal to 1 and less than or equal to H '/2, gamma is more than or equal to 1 and less than or equal to B, gamma is simultaneously used as the wave band mark of the wave band image, the wave band number of the tested remote sensing fusion image and the multispectral(γ)(x ', y') denotes { S }(γ)(x ', y') } pixel value of pixel point with coordinate position (x ', y'), R(γ)(x ", y") denotes { R }(γ)(x ", y") } the pixel value of the pixel point whose coordinate position is (x ", y").
Step ② _2, calculating and testing the full resolution quality evaluation predicted value of the remote sensing fusion image, and recording the predicted value as QfullThe specific process is as follows:
② _2a, using a moving window with a window size of 96 × 96 and a window moving step length of 1 pixel, moving in all band images of the tested remote sensing fused image to intercept subblocks with the same coordinate position and size of 96 × 96, marking the subblocks intercepted in the same coordinate position in all band images of the tested remote sensing fused image with the same block index number, marking each subblock intercepted in all band images of the tested remote sensing fused image with a band mark to indicate the source of the subblock, intercepting M ' subblocks from all band images of the tested remote sensing fused image after moving for M ' times, arranging pixel values of all pixel points in each intercepted subblock in sequence to form a column vector corresponding to the subblock, forming a column vector set corresponding to the M ' subblocks of each band mark intercepted to form a column vector set, and recording the column set formed by the column vectors corresponding to the M ' subblocks marked as gamma, intercepted M ' subblocks marked as gamma
Marking the intercepted wave band as a column formed by column vectors corresponding to M' sub-blocks with gamma being 1Vector collection
Defining a column vector set formed by column vectors corresponding to M' sub-blocks for testing red band column vector set and marking the intercepted band as gamma-2
Defining a column vector set formed by column vectors corresponding to M' sub-blocks for testing a green band column vector set and marking the intercepted band as gamma being 3
Defining a column vector set formed by column vectors corresponding to M' sub-blocks with the intercepted wave band marked as gamma being 4 as a test blue wave band column vector set
Is defined as testing near-infrared band column vector set, wherein M 'is positive integer, 1 is more than M and less than or equal to (W' -96) × (H '-96), M' is positive integer, 1 is less than or equal to M and less than or equal to M ', M' is simultaneously used as block index number of subblocks in the testing stage process,
has a dimension of (96 × 96) × 1,
to represent
The m "column vector in (1), i.e. the column vector representing the sub-block with the band marked as gamma and the block index number m",
to represent
The m "th column vector in (1), i.e. the column corresponding to the sub-block with the band denoted γ ═ 1 and the block index m ″The vector of the vector is then calculated,
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted γ ═ 2 and the block index m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted as γ ═ 3 and the block index m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block whose denoted band is γ ═ 4 and block index m".
Step ② _2b, obtaining the feature vectors of all sub-blocks with the same block index number in the same way according to the process from step ① _3 to step ① _7, and recording the feature vectors of all sub-blocks with block index number m' as Vtest,m”； wherein ,Vtest,m”Has a dimension of 16 × 1, i.e. according to
Calculating to obtain a test sample set of the brightness components; calculating LBP characteristic statistical histograms of all pixel points in the sub-blocks corresponding to each column vector in the test sample set; extracting an edge map of a subblock corresponding to each column vector in the test sample set by adopting a Canny edge detection operator; calculating edge histograms of all pixel points in the edge graph of the sub-block corresponding to each column vector in the test sample set; according to
Calculating the spectral characteristics of all subblocks with the same block index number; and obtaining the feature vectors of the subblocks with the same block index number according to the LBP feature statistical histogram of all pixel points in the subblock corresponding to each column vector in the test sample set, the edge histograms of all pixel points in the edge graph of the subblock corresponding to each column vector and the spectral features of the subblocks with the same block index number.
Step ② _2c, constructing original multivariate Gaussian (MVG) models of all band images of the remote sensing fused image in the same way according to the total M' eigenvectors obtained in the step ② _2b and the process of the step ① _8, and marking as the MVG models
wherein ,
to represent
The mean value vector of (a) is,
to represent
The covariance matrix of (2).
Step ② _2d according to
And
calculating Qfull，
wherein ,
is composed of
The transpose of (a) is performed,
is composed of
The inverse of (c).
Step ② _3, calculating and testing the reduced resolution quality evaluation prediction value of the remote sensing fusion image, and recording the value as QreduThe specific process is as follows:
step ② _3a, performing down-sampling operation on each waveband image of the tested remote sensing fusion image to obtain a corresponding down-sampling waveband image with width W '/2 and height H'/2, and performing down-sampling on { S }(γ)(x ', y') } the down-sampling waveband image with the width W '/2 and the height H'/2 obtained after the down-sampling operation is recorded as the down-sampling waveband image
Wherein x is more than or equal to 1 and less than or equal to W '/2, y is more than or equal to 1 and less than or equal to H'/2,
to represent
The pixel value of the pixel point with the middle coordinate position of (x ', y').
Step ② _3b, calculating the spectral similarity of the down-sampling waveband image obtained by the down-sampling operation of all waveband images of the multi-spectral image referenced by the test remote sensing fusion image and all waveband images of the test remote sensing fusion image, and recording as Qspe，
wherein ,Rx”,y”Coordinates in all band images of multispectral image representing reference of testing remote sensing fusion imageA vector, S, formed by connecting the pixel values of the pixels with the position (x', y ″)x”,y”The pixel values of pixel points with coordinate positions (x, y') in the downsampled waveband images obtained by downsampling all waveband images representing the tested remote sensing fusion image are connected in sequence to form a vector Rx”,y” and Sx”,y”All dimensions of (A) are B × 1, symbol "<>"is an inner product operation symbol, acrcos () is an inverse cosine function, the symbol" | | | | luminance2"is the 2-norm sign of the matrix.
Step ② _3c, calculating the spatial similarity of the down-sampling waveband images obtained by the down-sampling operation of all waveband images of the multispectral image referenced by the test remote sensing fusion image and all waveband images of the test remote sensing fusion image, and recording as Qspa，
wherein ,
represents Rx”,y”The average of the values of all the elements in (a),
denotes Sx”,y”The average of the values of all the elements in (a),
represents Rx”,y”The standard deviation of the values of all elements in (a),
denotes Sx”,y”The standard deviation of the values of all elements in (a),
denotes Sx”,y”The values of all elements in (1) and Rx”,y”Of all elements in (c), ω1 and ω2To control the parameters, in this example, take ω1＝0.0001，ω2＝0.0009。
Step ② _3d, according to Qspe and QspaCalculating Qredu，Qredu＝(Qspe)η×(Qspa)1-ηWherein η is a weight coefficient, and in this example η is 0.08.
Step ② _4, according to Qfull and QreduCalculating and testing the objective quality evaluation predicted value of the remote sensing fusion image and recording the predicted value as Qtest，
Where λ is a weight coefficient, and in this embodiment, λ is 0.1.
To further illustrate the feasibility and effectiveness of the method of the present invention, the method of the present invention was tested.
In this embodiment, the method of the present invention is used to test the remote sensing fusion image in the remote sensing fusion image database established by Ningbo university. 200 IKONOS satellites, 500 quick bird satellites, 500 Gaofen-1 satellites, 500 WorldView-4 satellites, 500 multispectral images and panchromatic images shot by the WorldView-2 satellites and 160 WorldView-3 satellites are collected in the remote sensing fusion image database, then the multispectral images and the panchromatic images are fused by 6 fusion methods to obtain corresponding remote sensing fusion images, and the average subjective score average value of each remote sensing fusion image is given. In this embodiment, 3 common objective parameters for evaluating the image quality are used as evaluation indexes, that is, a Pearson correlation coefficient (PLCC), a Spearman correlation coefficient (SROCC), a mean square error (RMSE), and the PLCC and the RMSE reflect the accuracy of the objective quality evaluation prediction value, and the SROCC reflects the monotonicity thereof under a nonlinear regression condition. The correlation between the objective quality evaluation prediction value obtained by the method and the average subjective score mean value is shown in the table 1.
TABLE 1 correlation between the objective quality evaluation prediction value obtained by the method of the present invention and the mean of the average subjective scores
| Satellite | PLCC | SROCC | RMSE |
| IKONOS | 0.6054 | 0.6085 | 11.38 |
| QuickBird | 0.5051 | 0.5043 | 11.90 |
| Gaofen-1 | 0.2966 | 0.2568 | 13.71 |
| WorldView-4 | 0.4697 | 0.4388 | 12.21 |
| WorldView-2 | 0.6738 | 0.6383 | 11.23 |
| WorldView-3 | 0.6383 | 0.6205 | 12.00 |
As can be seen from Table 1, the correlation between the objective quality evaluation prediction value obtained by the method of the present invention and the average subjective score mean is very high, which indicates that the objective evaluation result is more consistent with the result of human eye subjective perception, and is sufficient to illustrate the effectiveness of the method of the present invention.
Claims (3)
1. A method for evaluating fusion quality of remote sensing images is characterized by comprising a training stage and a testing stage;
the specific steps of the training phase process are as follows:
step ① _1, selecting N original multispectral images with width W and height H, each original multispectral image containing B different waveband images, and recording the gamma waveband image of the u-th original multispectral image as
1 st wave band image of u original multispectral image
The wave band image of the red wave band and the 2 nd wave band image of the u original multispectral image
Is a wave band image of a green wave band, and a 3 rd wave band image of a u original multispectral image
The wave band image of the blue wave band, the 4 th wave band image of the u original multispectral image
The band image is a near infrared band; wherein N is a positive integer, N is more than 1, B is a positive integer, u is more than or equal to 1 and less than or equal to N, gamma is a positive integer, gamma is more than or equal to 1 and less than or equal to B, and gamma is simultaneously used as a waveband mark of the waveband image;
① _2, cutting a sub-block with size of 96 × 96 from the same coordinate position in all band images of the original multispectral image, marking the sub-blocks cut from the same coordinate position in all band images of the original multispectral image with the same block index number, marking the sub-blocks cut from the same coordinate position in all band images of the original multispectral image with band marks to indicate the source of the sub-blocks, repeating the cutting process in a random cutting mode, cutting M times, then cutting M sub-blocks from all band images of the original multispectral image, arranging the pixel values of all pixel points in each cut sub-block to form a column vector corresponding to the sub-block, forming a column vector set with M 'sub-blocks of the same cut band mark band, and marking the M' sub-block corresponding to gamma vector set as a column set
The intercepted wave band is marked as a column vector set formed by column vectors corresponding to M' sub-blocks with the gamma being 1
Defining as training red wave band column vector set, marking the intercepted wave band as column vector set formed from column vectors correspondent to M' sub-blocks whose gamma is 2
Defining as training green band column vector set, marking intercepted band as M' sub-components of gamma-3Column vector set composed of column vectors corresponding to blocks
Defining a column vector set formed by column vectors corresponding to M' sub-blocks with the intercepted wave band marked as gamma being 4 as a training blue wave band column vector set
Defining as training near infrared wave band column vector set, where M is positive integer, M is greater than or equal to 1, M 'is N × M, M is positive integer, M is greater than or equal to 1 and less than or equal to M', M is used as block index number of sub-block in training phase,
to represent
The m-th column vector in (e), i.e. the column vector representing the sub-block with the band labeled gamma and the block index m,
are (96 × 96) × 1,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose denoted band is γ ═ 1 and whose block index is m,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose denoted band is γ ═ 2 and block index m,
to represent
The mth column vector in (1), i.e. the column vector corresponding to the sub-block whose representation band is denoted by γ ═ 3 and whose block index is m,
to represent
The mth column vector in (1), that is, the column vector corresponding to the sub-block whose band is denoted by γ ═ 4 and whose block index is m;
① _3, according to the training red band column vector set
Training a green band column vector set
And training a blue band column vector set
The training sample set of luminance components is obtained by calculation and is marked as { Ym|1≤m≤M'}，
wherein ,YmHas the dimension of (96 × 96) × 1 and symbols
Is a rounding down operation symbol;
step ① _4, calculating { YmThe statistical histogram of LBP features of all pixel points in the sub-block corresponding to each column vector in the sub-block with the value of 1 ≦ M ≦ M' } is calculated according to the statistical histogram of the LBP features of all the pixel points in the sub-block corresponding to each column vector in the sub-block with the value of YmThe sub-column vector corresponding to the mth column vector in the M is more than or equal to 1 and less than or equal to M' }The statistical histogram of LBP features of all the pixels in the block is marked as Fm； wherein ,FmHas a dimension of 8 × 1;
step ① _5, extracting { Y ] by adopting Canny edge detection operatormThe M is more than or equal to 1 and less than or equal to M' } is the edge map of the sub-block corresponding to each column vector; then calculate { YmThe edge histograms of all pixel points in the edge graph of the subblock corresponding to each column vector in the |1 ≦ M ≦ M' } are calculated according to the edge histograms of all the pixel points in the subblock corresponding to each column vector in the column vectormThe edge histograms of all pixel points in the edge graph of the sub-block corresponding to the mth column vector in the |1 ≦ M ≦ M' } are marked as Gm； wherein ,GmHas a dimension of 2 × 1;
① _6, according to the training red band column vector set
Training a green band column vector set
Training blue band column vector set
And training a near infrared band column vector set
Calculating the spectral characteristics of all subblocks with the same block index number, and recording the spectral characteristics of all subblocks with the block index number m as Hm； wherein ,HmHas a dimension of 6 × 1;
step ① _7, according to { Y }mI1 is less than or equal to M is less than or equal to M '} and the LBP feature statistical histogram of all pixel points in the sub-block corresponding to each column vector in the { Y ≦ M' } is obtainedmThe edge histograms of all pixel points in the edge graph of the subblock corresponding to each column vector in the I1 is more than or equal to M and less than or equal to M', the spectral characteristics of the subblocks with the same block index number are obtained, the characteristic vectors of the subblocks with the same block index number are obtained, and the characteristic vectors of the subblocks with the block index number M are marked as Vm，Vm＝[Fm T,Gm T,Hm T]T； wherein ,VmHas a dimension of 16 × 1 [ F ]m T,Gm T,Hm T]Is shown asm T、Gm T and Hm TAre connected to form a vector, Fm TIs FmTranspose of (G)m TIs GmTranspose of (H)m TIs HmIs transposed, [ F ]m T,Gm T,Hm T]TIs [ F ]m T,Gm T,Hm T]Transposing;
step ① _8, constructing original multivariate Gaussian models of all band images of all original multispectral images according to the total M' eigenvectors obtained in the step ① _7, and marking as the original multivariate Gaussian models
wherein ,
to represent
The mean value vector of (a) is,
to represent
The covariance matrix of (a);
the specific steps of the test phase process are as follows:
in step ② _1, any test remote sensing fused image with width W 'and height H' includes B different waveband images, and the gamma waveband image of the test remote sensing fused image is marked as { S(γ)(x ', y') }; the multispectral image with the width W '/2 and the height H'/2 referred by the test remote sensing fusion image also comprises B different waveband images, and the test remote sensing fusion image is fusedThe gamma band image of the multispectral image referenced by the image is marked as { R(γ)(x ", y") }; wherein x 'is more than or equal to 1 and less than or equal to W', y 'is more than or equal to 1 and less than or equal to H', x is more than or equal to 1 and less than or equal to W '/2, y is more than or equal to 1 and less than or equal to H'/2, gamma is more than or equal to 1 and less than or equal to B, gamma is simultaneously used as a waveband mark of a waveband image, S is(γ)(x ', y') denotes { S }(γ)(x ', y') } pixel value of pixel point with coordinate position (x ', y'), R(γ)(x ", y") denotes { R }(γ)(x ", y") } the pixel value of the pixel point whose coordinate position is (x ", y");
step ② _2, calculating and testing the full resolution quality evaluation predicted value of the remote sensing fusion image, and recording the predicted value as QfullThe specific process is as follows:
② _2a, using a moving window with a window size of 96 × 96 and a window moving step length of 1 pixel, moving in all band images of the tested remote sensing fused image to intercept subblocks with the same coordinate position and size of 96 × 96, marking the subblocks intercepted in the same coordinate position in all band images of the tested remote sensing fused image with the same block index number, marking each subblock intercepted in all band images of the tested remote sensing fused image with a band mark to indicate the source of the subblock, intercepting M ' subblocks from all band images of the tested remote sensing fused image after moving for M ' times, arranging pixel values of all pixel points in each intercepted subblock in sequence to form a column vector corresponding to the subblock, forming a column vector set corresponding to the M ' subblocks of each band mark intercepted to form a column vector set, and recording the column set formed by the column vectors corresponding to the M ' subblocks marked as gamma, intercepted M ' subblocks marked as gamma
Marking the intercepted wave band as a column vector set formed by column vectors corresponding to M' sub-blocks with gamma being 1
Defining a column vector set formed by column vectors corresponding to M' sub-blocks for testing red band column vector set and marking the intercepted band as gamma-2
Defining a column vector set formed by column vectors corresponding to M' sub-blocks for testing a green band column vector set and marking the intercepted band as gamma being 3
Defining a column vector set formed by column vectors corresponding to M' sub-blocks with the intercepted wave band marked as gamma being 4 as a test blue wave band column vector set
Is defined as testing near-infrared band column vector set, wherein M 'is positive integer, 1 is more than M and less than or equal to (W' -96) × (H '-96), M' is positive integer, 1 is less than or equal to M and less than or equal to M ', M' is simultaneously used as block index number of subblocks in the testing stage process,
has a dimension of (96 × 96) × 1,
to represent
The m "column vector in (1), i.e. the column vector representing the sub-block with the band marked as gamma and the block index number m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted γ and the block index m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted γ ═ 2 and the block index m",
to represent
The m "th column vector in (1), i.e. the column vector corresponding to the sub-block with the band denoted as γ ═ 3 and the block index m",
to represent
The m "column vector in (1), that is, the column vector corresponding to the sub-block whose band is denoted by γ ═ 4 and whose block index is m";
step ② _2b, obtaining the feature vectors of all sub-blocks with the same block index number in the same way according to the process from step ① _3 to step ① _7, and recording the feature vectors of all sub-blocks with block index number m' as Vtest,m”； wherein ,Vtest,m”Has a dimension of 16 × 1;
step ② _2c, constructing original multivariate Gaussian models of all wave band images of the remote sensing fused image to be tested in the same way according to the M' eigenvectors obtained in the step ② _2b and the process of the step ① _8, and marking as the original multivariate Gaussian models
wherein ,
to represent
The mean value vector of (a) is,
to represent
The covariance matrix of (a);
step ② _2d according to
And
calculating Qfull，
wherein ,
is composed of
The transpose of (a) is performed,
is composed of
The inverse of (1);
step ② _3, calculating and testing the reduced resolution quality evaluation prediction value of the remote sensing fusion image, and recording the value as QreduThe specific process is as follows:
step ② _3a, performing down-sampling operation on each waveband image of the tested remote sensing fusion image to obtain a corresponding down-sampling waveband image with width W '/2 and height H'/2, and performing down-sampling on { S }(γ)(x ', y') } the down-sampling waveband image with the width W '/2 and the height H'/2 obtained after the down-sampling operation is recorded as the down-sampling waveband image
Wherein x is more than or equal to 1 and less than or equal to W '/2, y is more than or equal to 1 and less than or equal to H'/2,
to represent
The pixel value of the pixel point with the middle coordinate position of (x', y ");
step ② _3b, calculating the spectral similarity of the down-sampling waveband image obtained by the down-sampling operation of all waveband images of the multi-spectral image referenced by the test remote sensing fusion image and all waveband images of the test remote sensing fusion image, and recording as Qspe，
wherein ,Rx”,y”A vector S formed by sequentially connecting pixel values of pixel points with coordinate positions (x, y') in all band images of the multispectral image for testing remote sensing fusion image referencex”,y”The pixel values of pixel points with coordinate positions (x, y') in the downsampled waveband images obtained by downsampling all waveband images representing the tested remote sensing fusion image are connected in sequence to form a vector Rx”,y” and Sx”,y”All dimensions of (A) are B × 1, symbol "<>"is an inner product operation symbol, acrcos () is an inverse cosine function, the symbol" | | | | luminance2"is the 2-norm sign of the matrix;
step ② _3c, calculating the spatial similarity of the down-sampling waveband images obtained by the down-sampling operation of all waveband images of the multispectral image referenced by the test remote sensing fusion image and all waveband images of the test remote sensing fusion image, and recording as Qspa，
wherein ,
represents Rx”,y”All elements inThe mean value of the values of the elements,
denotes Sx”,y”The average of the values of all the elements in (a),
represents Rx”,y”The standard deviation of the values of all elements in (a),
denotes Sx”,y”The standard deviation of the values of all elements in (a),
denotes Sx”,y”The values of all elements in (1) and Rx”,y”Of all elements in (c), ω1 and ω2Is a control parameter;
step ② _3d, according to Qspe and QspaCalculating Qredu，Qredu＝(Qspe)η×(Qspa)1-ηWherein η is a weight coefficient;
step ② _4, according to Qfull and QreduCalculating and testing the objective quality evaluation predicted value of the remote sensing fusion image and recording the predicted value as Qtest，
Wherein λ is a weight coefficient.
2. The remote sensing image fusion quality evaluation method according to claim 1, wherein in step ① _6, HmThe acquisition process comprises the following steps:
step ① _6a, calculation
And
spectral feature of (1), denoted as d12，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
The standard deviation of the values of all elements in (a),
representation calculation
Standard deviation of the values of all elements in (a);
step ① _6b, calculate
And
spectral feature of (1), denoted as d13，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
Standard deviation of the values of all elements in (a);
step ① _6c, calculate
And
spectral feature of (1), denoted as d14，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a),
representation calculation
The average of the values of all the elements in (a),
representation calculation
Standard deviation of the values of all elements in (a);
step ① _6d, calculate
And
spectral feature of (1), denoted as d23，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a);
step ① _6e, calculate
And
spectral feature of (1), denoted as d24，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a);
step ① _6f, calculate
And
spectral feature of (1), denoted as d34，
wherein ,
representation calculation
All values of (1) and
the covariance of the values of all elements in (a);
step ① _6g, according to d12、d13、d14、d23、d24 and d34Obtaining Hm，Hm＝[d12,d13,d14,d23,d24,d34]T； wherein ,HmHas a dimension of 6 × 1 [ d ]12,d13,d14,d23,d24,d34]Is shown as12、d13、d14、d23、d24 and d34Are connected to form a vector, [ d ]12,d13,d14,d23,d24,d34]TIs [ d ]12,d13,d14,d23,d24,d34]The transposing of (1).
3. The remote sensing image fusion quality evaluation method according to claim 1 or 2, characterized in that in step ① _8,
the acquisition process comprises the following steps:
step ① _8a, obtaining the feature matrix of all sub-blocks with the same block index number according to the total M' feature vectors obtained in step ① _7, and marking as V,
wherein the dimension of V is M' × 16, (V)1)TIs a V1Transpose of (V)1Feature vector representing subblocks with all block indices of 1, (V)2)TIs a V2Transpose of (V)2Feature vector representing subblocks with all chunk indices of 2, (V)m)TIs a VmTranspose of (V)M')TIs a VM'Transpose of (V)M'A feature vector representing a subblock having all block indexes M';
step ① _8b, calculating the mean vector of the column vectors in V, and recording as
wherein ,
has a dimension of 1 × 16, and V (m,1) represents the m-th row in VThe value of the element of column 1, V (m,2) denotes the value of the element of row m, column 2 in V, and V (m,16) denotes the value of the element of row m, column 16 in V;
and calculate the covariance matrix of the column vectors in V, noted
wherein ,
dimension 16 × 16, cov (a, b) represents calculating the covariance of the values of all elements in vector a and all elements in vector b;
step ① _8c, will
As
Mean vector of (2) will
As
Is constructed to obtain a covariance matrix
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010385384.9A CN111681207B (en) | 2020-05-09 | 2020-05-09 | Remote sensing image fusion quality evaluation method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010385384.9A CN111681207B (en) | 2020-05-09 | 2020-05-09 | Remote sensing image fusion quality evaluation method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111681207A true CN111681207A (en) | 2020-09-18 |
| CN111681207B CN111681207B (en) | 2023-10-27 |
Family
ID=72433468
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010385384.9A Active CN111681207B (en) | 2020-05-09 | 2020-05-09 | Remote sensing image fusion quality evaluation method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111681207B (en) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112258464A (en) * | 2020-10-14 | 2021-01-22 | 宁波大学 | Full-reference remote sensing image fusion quality evaluation method |
| CN113240758A (en) * | 2021-05-28 | 2021-08-10 | 珠江水利委员会珠江水利科学研究院 | Remote sensing image fusion method, system, equipment and medium based on fusion derivative index |
| CN115049576A (en) * | 2021-02-26 | 2022-09-13 | 北京小米移动软件有限公司 | Image quality evaluation method and device, equipment and storage medium |
| CN116129278A (en) * | 2023-04-10 | 2023-05-16 | 牧马人(山东)勘察测绘集团有限公司 | Land utilization classification and identification system based on remote sensing images |
Citations (14)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102542549A (en) * | 2012-01-04 | 2012-07-04 | 西安电子科技大学 | Multi-spectral and panchromatic image super-resolution fusion method based on compressive sensing |
| CN103383775A (en) * | 2013-07-02 | 2013-11-06 | 中国科学院东北地理与农业生态研究所 | Method for evaluating remote-sensing image fusion effect |
| CN104851091A (en) * | 2015-04-28 | 2015-08-19 | 中山大学 | Remote sensing image fusion method based on convolution enhancement and HCS transform |
| CN104881867A (en) * | 2015-05-13 | 2015-09-02 | 华中科技大学 | Method for evaluating quality of remote sensing image based on character distribution |
| CN106023111A (en) * | 2016-05-23 | 2016-10-12 | 中国科学院深圳先进技术研究院 | Image fusion quality evaluating method and system |
| CN106157317A (en) * | 2016-07-21 | 2016-11-23 | 武汉大学 | The high-resolution remote sensing image fusion rules method guided based on dispersion tensor |
| CN107040775A (en) * | 2017-03-20 | 2017-08-11 | 宁波大学 | A kind of tone mapping method for objectively evaluating image quality based on local feature |
| CN107862678A (en) * | 2017-10-19 | 2018-03-30 | 宁波大学 | A kind of eye fundus image reference-free quality evaluation method |
| CN108682005A (en) * | 2018-04-25 | 2018-10-19 | 西北工业大学 | Half based on covariance matrix feature refers to 3D composograph quality evaluating methods |
| CN109191450A (en) * | 2018-09-11 | 2019-01-11 | 宁波大学 | A kind of remote sensing image fusion quality evaluating method |
| CN109919899A (en) * | 2017-12-13 | 2019-06-21 | 香港纺织及成衣研发中心有限公司 | The method for evaluating quality of image based on multispectral imaging |
| CN109949270A (en) * | 2019-01-28 | 2019-06-28 | 西北工业大学 | Multispectral and full-colour image based on region convolutional network merges space quality evaluation method |
| CN110555843A (en) * | 2019-09-11 | 2019-12-10 | 浙江师范大学 | High-precision non-reference fusion remote sensing image quality analysis method and system |
| AU2020100179A4 (en) * | 2020-02-04 | 2020-03-19 | Huang, Shuying DR | Optimization Details-Based Injection Model for Remote Sensing Image Fusion |
-
2020
- 2020-05-09 CN CN202010385384.9A patent/CN111681207B/en active Active
Patent Citations (14)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102542549A (en) * | 2012-01-04 | 2012-07-04 | 西安电子科技大学 | Multi-spectral and panchromatic image super-resolution fusion method based on compressive sensing |
| CN103383775A (en) * | 2013-07-02 | 2013-11-06 | 中国科学院东北地理与农业生态研究所 | Method for evaluating remote-sensing image fusion effect |
| CN104851091A (en) * | 2015-04-28 | 2015-08-19 | 中山大学 | Remote sensing image fusion method based on convolution enhancement and HCS transform |
| CN104881867A (en) * | 2015-05-13 | 2015-09-02 | 华中科技大学 | Method for evaluating quality of remote sensing image based on character distribution |
| CN106023111A (en) * | 2016-05-23 | 2016-10-12 | 中国科学院深圳先进技术研究院 | Image fusion quality evaluating method and system |
| CN106157317A (en) * | 2016-07-21 | 2016-11-23 | 武汉大学 | The high-resolution remote sensing image fusion rules method guided based on dispersion tensor |
| CN107040775A (en) * | 2017-03-20 | 2017-08-11 | 宁波大学 | A kind of tone mapping method for objectively evaluating image quality based on local feature |
| CN107862678A (en) * | 2017-10-19 | 2018-03-30 | 宁波大学 | A kind of eye fundus image reference-free quality evaluation method |
| CN109919899A (en) * | 2017-12-13 | 2019-06-21 | 香港纺织及成衣研发中心有限公司 | The method for evaluating quality of image based on multispectral imaging |
| CN108682005A (en) * | 2018-04-25 | 2018-10-19 | 西北工业大学 | Half based on covariance matrix feature refers to 3D composograph quality evaluating methods |
| CN109191450A (en) * | 2018-09-11 | 2019-01-11 | 宁波大学 | A kind of remote sensing image fusion quality evaluating method |
| CN109949270A (en) * | 2019-01-28 | 2019-06-28 | 西北工业大学 | Multispectral and full-colour image based on region convolutional network merges space quality evaluation method |
| CN110555843A (en) * | 2019-09-11 | 2019-12-10 | 浙江师范大学 | High-precision non-reference fusion remote sensing image quality analysis method and system |
| AU2020100179A4 (en) * | 2020-02-04 | 2020-03-19 | Huang, Shuying DR | Optimization Details-Based Injection Model for Remote Sensing Image Fusion |
Non-Patent Citations (4)
| Title |
|---|
| 刘军;邵振峰;: "基于特征结构相似度的遥感影像融合质量评价指标" * |
| 孟祥超等;: "基于多分辨率分析的GF-5和GF-1遥感影像空―谱融合" * |
| 柳文等;: "ALOS全色波段与多光谱影像融合方法的比较研究" * |
| 陈鹰;郭睿;: "非负矩阵分解在遥感图像融合中的应用" * |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112258464A (en) * | 2020-10-14 | 2021-01-22 | 宁波大学 | Full-reference remote sensing image fusion quality evaluation method |
| CN112258464B (en) * | 2020-10-14 | 2023-04-28 | 宁波大学 | Full-reference remote sensing image fusion quality evaluation method |
| CN115049576A (en) * | 2021-02-26 | 2022-09-13 | 北京小米移动软件有限公司 | Image quality evaluation method and device, equipment and storage medium |
| CN113240758A (en) * | 2021-05-28 | 2021-08-10 | 珠江水利委员会珠江水利科学研究院 | Remote sensing image fusion method, system, equipment and medium based on fusion derivative index |
| CN116129278A (en) * | 2023-04-10 | 2023-05-16 | 牧马人(山东)勘察测绘集团有限公司 | Land utilization classification and identification system based on remote sensing images |
Also Published As
| Publication number | Publication date |
|---|---|
| CN111681207B (en) | 2023-10-27 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111681207A (en) | Remote sensing image fusion quality evaluation method | |
| CN111080629B (en) | Method for detecting image splicing tampering | |
| CN109191450B (en) | Remote sensing image fusion quality evaluation method | |
| CN111127374B (en) | Pan-sharing method based on multi-scale dense network | |
| Kotwal et al. | A novel approach to quantitative evaluation of hyperspectral image fusion techniques | |
| CN109978854B (en) | Screen content image quality evaluation method based on edge and structural features | |
| CN110363134B (en) | Human face shielding area positioning method based on semantic segmentation | |
| CN111932532B (en) | Method for evaluating capsule endoscope without reference image, electronic device, and medium | |
| CN111274865A (en) | Remote sensing image cloud detection method and device based on full convolution neural network | |
| CN108765359A (en) | A kind of fusion method of target in hyperspectral remotely sensed image and full-colour image based on JSKF models and NSCT technologies | |
| CN107146220B (en) | A kind of universal non-reference picture quality appraisement method | |
| CN111639587A (en) | Hyperspectral image classification method based on multi-scale spectrum space convolution neural network | |
| CN109740592B (en) | Image quality parameter-free assessment method based on memory | |
| CN113177548B (en) | Key area identification method for immune fixed electrophoresis | |
| CN108053396A (en) | A kind of more distorted image quality without with reference to evaluation method | |
| CN104103064A (en) | Reference-free noise image quality evaluation method based on gradient similarity | |
| CN112258464B (en) | Full-reference remote sensing image fusion quality evaluation method | |
| CN107169946A (en) | Image interfusion method based on non-negative sparse matrix Yu hypersphere color transformation | |
| CN114332534B (en) | Hyperspectral image small sample classification method | |
| KR102160687B1 (en) | Aviation image fusion method | |
| CN107180419B (en) | Median filtering detection method based on PCA (principal component analysis) network | |
| CN110796635B (en) | Light field image quality evaluation method based on shear wave transformation | |
| CN116758019A (en) | Multi-exposure fusion light field image quality evaluation method based on dynamic and static region division | |
| CN113450319B (en) | Super-resolution reconstruction image quality evaluation method based on KLT technology | |
| CN113192003B (en) | Spliced image quality evaluation method |
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 | ||
| TA01 | Transfer of patent application right | ||
| TA01 | Transfer of patent application right |
Effective date of registration: 20230925 Address after: South side of 3rd floor, Building D, Xixi Ginza, No. 780 Wen'er West Road, Xihu District, Hangzhou City, Zhejiang Province, 310000 Applicant after: Siwei Gaojing satellite remote sensing Co.,Ltd. Address before: 315211, Fenghua Road, Jiangbei District, Zhejiang, Ningbo 818 Applicant before: Ningbo University |
|
| GR01 | Patent grant | ||
| GR01 | Patent grant |