A weight dictionary thermal infrared hyperspectral mixed pixel temperature emissivity solving method

Abstract

The present application belongs to the technical field of remote sensing, and relates to a weight dictionary thermal infrared hyperspectral mixed pixel temperature emissivity solving method. The method comprises: obtaining simulation data and performing inversion channel screening according to the simulation data; obtaining metadata and quality control marks and performing effective pixel screening; performing atmospheric correction; constructing an initial value dictionary and training to obtain an emissivity dictionary; calculating a ground surface emissivity spectral vector, reconstructing the emissivity spectral vector based on a ROMP algorithm, performing weight redistribution through a weight distribution algorithm, and iteratively solving to restore the ground surface temperature and the noise-free and unbiased emissivity spectrum of the mixed pixel. The present application combines the ROMP algorithm and the weight distribution algorithm into a WOMP algorithm, overcomes the greediness of the traditional algorithm, cooperates with the initial value dictionary and the emissivity dictionary, reconstructs the unbiased emissivity spectrum, realizes the synchronous inversion of the ground surface temperature and the emissivity spectrum, and significantly improves the stability, unbiasedness and decoupling capability of the mixed pixel inversion.

Description

Technical Field

[0001] This invention belongs to the field of remote sensing technology, specifically relating to a method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral hybrid pixel. Background Technology

[0002] Land surface temperature (LST) and land surface emissivity (LSE) are two key variables controlling terrestrial thermal infrared radiation. LST reflects the thermal state of the land surface and is an important parameter in global climate change research, surface energy balance analysis, urban heat island effect monitoring, agricultural drought early warning, and disaster assessment. LSE reflects the spectral characteristics of surface materials and is closely related to the composition, structure, and physical state of the land surface. It can be applied to land cover classification, mineral identification, soil moisture retrieval, and vegetation growth monitoring.

[0003] Because surface temperature and surface emissivity are highly coupled in thermal infrared radiance observations, temperature-emissivity separation (TES) has always been a core scientific problem in quantitative thermal infrared remote sensing. TES is essentially a solution process for an underdetermined problem: under ideal conditions with complete atmospheric correction, given the observations of N channels in the thermal infrared hyperspectral system, N radiative transfer equations can be established, but there are N+1 unknown parameters to solve for (including the emissivity of the N channels and the surface temperature). Therefore, physical constraints or parameterization methods must be introduced to reduce the number of unknowns in order to achieve simultaneous calculation of temperature and emissivity. Furthermore, atmospheric correction errors, surface spectral differences, instrument noise, and mixed pixel effects further increase the uncertainty of simultaneous LST and LSE inversion.

[0004] Existing temperature-emissivity separation methods can be mainly classified into the following categories: (1) The method based on the assumption of minimum emissivity is simple to calculate and has high operating efficiency, but it has a large error in regions with relatively flat emissivity spectra (such as water bodies and ice and snow) and is sensitive to noise; The method based on spectral smoothing constraints performs well on most natural surfaces, but it over-smooths on features with distinct spectral characteristics (such as minerals and rocks), losing important spectral details. (3) The method based on linear emissivity constraint has a performance that is highly dependent on the choice of basis function. In the mixed pixel scenario, since the emissivity spectrum is a mixture of the spectra of multiple ground objects, the predetermined basis function often cannot accurately describe it, resulting in constraint failure. (4) Signal processing-based methods have the potential to reduce effective spectral details while filtering out noise; temperature-emissivity separation methods based on dictionary sparse representation use a sparse dictionary to represent the emissivity spectrum as a linear combination of a few atoms in the dictionary, thereby reducing the number of unknowns. However, commonly used orthogonal matching pursuit algorithms are greedy algorithms, selecting the atom most relevant to the residual in each iteration, which can easily make the inversion results overly constrained by the training spectral library, especially unsuitable for mixed pixel scenarios. In mixed pixels, the emissivity spectrum is a linear mixture of the spectra of multiple ground objects. Traditional orthogonal matching pursuit algorithms tend to select spectral atoms of a single ground object, resulting in biased inversion results that cannot accurately reflect the true emissivity characteristics of mixed pixels.

[0005] In summary, existing temperature-emissivity separation methods still suffer from insufficient stability, biased inversion results, and poor decoupling capabilities in scenarios such as humid atmospheres, gray targets, and complex mixed pixels. These methods cannot meet the needs of large-scale and high-precision quantitative inversion of surface thermal parameters. Therefore, it is urgent to develop a temperature emissivity calculation method suitable for thermal infrared hyperspectral mixed pixels to improve the stability, unbiasedness, and decoupling capabilities of the inversion. Summary of the Invention

[0006] To address the aforementioned technical problems, the present invention proposes a solution. The present invention provides a method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel, comprising: Acquire simulation data and perform inversion channel selection based on the simulation data; the simulation data includes entrance pupil radiance, atmospheric transmittance, atmospheric upward radiation, atmospheric downward radiation, and surface radiance. Obtain metadata and quality control identifiers for thermal infrared hyperspectral payload observations, and perform effective pixel screening to obtain the effective entrance pupil radiance of each channel of the thermal infrared hyperspectral payload; Atmospheric correction is performed on the effective entrance pupil radiance to obtain the surface radiance and atmospheric downdraft radiation; An initial value dictionary is constructed, and the emissivity spectra of typical ground features are selected to train the initial value dictionary to obtain the emissivity dictionary; The largest surface brightness temperature value is selected as the initial value of surface temperature. The surface emissivity spectral vector is calculated and reconstructed based on the ROMP algorithm. The weights of the reconstructed emissivity spectral vector are then redistributed using a weight allocation algorithm. The surface temperature and the noise-free and offset-free emissivity spectra of the mixed pixels are then recovered through iterative solution.

[0007] In some optional embodiments, acquiring simulation data and performing inversion channel filtering based on the simulation data includes: Multiple atmospheric conditions and underlying surface parameters are used as inputs; Hyperspectral simulation data were obtained using a radiative transfer model. The spectral response function is convolved with the hyperspectral simulation data to obtain the simulation data; Based on simulation data, channels that are limited by convolution distortion effects and atmospheric transmittance less than a set threshold are eliminated.

[0008] In some optional embodiments, effective pixel screening is performed to obtain the effective entrance pupil radiance of each channel of the thermal infrared hyperspectral load, including: removing abnormal pixels based on metadata and quality control identifiers, and performing radiometric correction to obtain the entrance pupil radiance of each channel of the thermal infrared hyperspectral load.

[0009] In some optional embodiments, atmospheric correction is performed on the effective entrance pupil radiance to obtain the surface radiance and atmospheric downdraft radiation, including: Atmospheric parameters at the observation time and within the observation range are obtained. The atmospheric parameters and the observation angle are input into the radiative transfer model to calculate atmospheric upward radiation, atmospheric transmittance and atmospheric downward radiation. set up For real numbers, This represents the total number of channels in the thermal infrared hyperspectral image that can be used for inversion. Represents the surface radiance vector. The entrance pupil radiance vector, Represents the upward radiation vector of the atmosphere. This represents Hadema division, which involves dividing two vectors element by element. Let be the atmospheric transmittance vector, and let be the surface radiance of each channel: ; set up: This represents the surface emissivity spectral vector. This represents the Hadema product, which is the element-wise multiplication of two vectors. Represents the Planck function; The center wavelength vector of the channel. It refers to the surface temperature; Indicates length is A column vector whose elements are all 1; Here is the downward atmospheric radiation vector; the surface radiance satisfies the radiative transfer equation: .

[0010] In some optional embodiments, an initial value dictionary is constructed, including: Select several typical land cover emissivity spectra from the land cover spectral library; The emissivity spectra of typical ground features are convolved using the spectral response function of the thermal infrared hyperspectral channel to obtain the convolved spectral vector; Arrange the convolved spectral vectors in the column direction to form an initial value dictionary.

[0011] In some optional embodiments, the emissivity spectra of typical ground features are selected to train the initial value dictionary to obtain the emissivity dictionary, including: constructing a training data matrix, and using the orthogonal matching pursuit algorithm and the K-singular value decomposition algorithm alternately to train the initial value dictionary with the training data matrix as input.

[0012] In some optional embodiments, constructing a training data matrix includes: Let the spectral vector after each convolution be... ; Indicates the first One normalized spectral vector; Representing the An empirically constrained spectral vector; For mapping functions, it represents the operation of integrating normalized spectral vectors into empirically constrained spectral vectors; The empirical constraint matrix is ​​formed by concatenating empirical constraint spectral vectors along their column directions; the matrix and All belong to Through the empirical constraint matrix respectively By adding two different sets of white noise, we get: ; ; ; ; Training the initial value dictionary includes: set up: The number of atoms listed in the dictionary; This represents the total number of channels available for inversion in thermal infrared hyperspectral imaging. The dictionary to be trained; Represents the Frobenius norm; the training data matrix is ;No. indivual The sparse coefficient vector corresponding to the column vector sample in the data is represented as follows: ; The total number of spectra of ground features; The index of the spectral column vector is [index]; the sparse coefficient matrix is ​​[matrix]. , is represented as a matrix where sparse coefficient vectors are concatenated along their column directions, i.e. ; Represents the number of non-zero elements in a vector; Let be the sparsity of the dictionary training; then training on the initial dictionary can be represented as: .

[0013] In some optional embodiments, the emissivity spectral vector is reconstructed based on the ROMP algorithm, including: A compressed observation matrix is ​​constructed based on the reconstructed emissivity spectrum vector containing noise and bias. Generate a random Gaussian matrix and perform QR decomposition to obtain an orthogonal matrix; The empirically constrained spectral vector is obtained by normalizing and empirically constraining the emissivity spectral vector containing noise and bias. The compressed vector is calculated based on the orthogonal matrix, the compressed observation matrix, and the empirically constrained spectral vector; the compressed dictionary is calculated based on the orthogonal matrix, the compressed observation matrix, and the trained emissivity dictionary. Solve for the sparse coefficients using the OMP algorithm; Reconstruct the normalized spectral vector based on the trained emissivity dictionary and sparse coefficients; Based on the reconstructed normalized spectral vector, the emissivity spectral vector is reconstructed.

[0014] In some optional embodiments, the reconstructed emissivity spectral vector is weighted and reallocated to recover a noise-free and offset-free emissivity spectrum of the mixed pixels, including: Using an emissivity spectral vector that incorporates various ground features, noise, and biases as input, a constrained quadratic programming optimization method is used to solve for the preliminary weight estimates of the emissivity spectral vector. The candidate dictionary is constructed by selecting emissivity spectral vectors from the initial value dictionary, and calculating the residual vectors based on the emissivity spectral vectors that are mixed with various land cover information, noise and bias, and the emissivity spectral vectors that are re-estimated based on weights and candidate dictionaries. The candidate dictionary after residual allocation is then constructed. The ROMP algorithm is used to reconstruct the emissivity spectrum of each candidate in the candidate dictionary after residual allocation, and the weights of the reconstructed emissivity spectrum are redistributed to restore the noise-free and offset-free emissivity spectrum of the mixed pixels.

[0015] In some optional embodiments, iteratively solving to recover the surface temperature and the noise-free and offset-free emissivity spectrum of the mixed pixels includes: Based on the Planck function, the atmospherically corrected surface radiance is converted into surface brightness temperature, and the largest surface brightness temperature value among several channels is selected as the initial value of surface temperature. The surface emissivity spectral vector is directly calculated based on the current estimated surface temperature and the radiative transfer equation. Using the surface emissivity spectral vector as input, the emissivity spectral vector is reconstructed based on the ROMP algorithm, and the weights of the reconstructed emissivity spectral vector are redistributed using a weight allocation algorithm to recover the unbiased emissivity spectral vector. Establish a cost function and use the golden section method to find the temperature value corresponding to the minimum cost function, which is used as the optimized surface temperature value for this iteration. Linear fitting is performed on the cost function values ​​in the last few iterations. If the absolute value of the slope of the fitted line is less than or the number of iterations exceeds the set number, the iteration terminates, and the optimized surface temperature value is output as the final inversion value of the surface temperature, and the unbiased emissivity spectral vector is output as the final inversion value of the surface emissivity spectrum; otherwise, the iteration is repeated.

[0016] The beneficial effects of this invention are: (1) This invention addresses the challenge of simultaneously retrieving surface temperature and emissivity from thermal infrared hyperspectral remote sensing data. It utilizes the compressed sensing theory of sparse dictionaries to sparsely represent emissivity, thereby reducing the number of unknowns in the inversion equation and enabling the simultaneous calculation of surface temperature and surface emissivity. (2) This invention improves the training method of the emissivity dictionary. It not only directly integrates the initial convolutional emissivity spectrum into the initial value dictionary to provide a reference for subsequent weight estimation, but also uses mean normalization, empirical constraints and Gaussian noise augmentation techniques, combined with OMP and KSVD algorithms, to train an emissivity dictionary with spectral shape and value range information. (3) The method of the present invention uses the ROMP algorithm and the weight allocation algorithm, and combines the two into the WOMP algorithm to overcome the greediness of the traditional OMP algorithm. It also uses the initial value dictionary and emission dictionary to reconstruct the unbiased emissivity spectrum to meet the needs of simultaneous inversion of surface temperature and surface emissivity spectrum of mixed pixels. The WOMP algorithm can provide the category weight information of the emissivity spectrum, that is, the probability that the emissivity is biased towards a certain type of land cover, thereby providing direct discrimination information for subsequent decoupling of mixed pixels. Attached Figure Description

[0017] Figure 1 This is a schematic diagram of a weighted dictionary thermal infrared hyperspectral hybrid pixel temperature emissivity calculation method provided in an embodiment of the present invention; Figure 2 This is a flowchart of the inversion channel selection based on simulated data provided in an embodiment of the present invention; Figure 3 This is a flowchart provided by an embodiment of the present invention for obtaining entrance pupil radiance by filtering effective pixels based on metadata; Figure 4 This is a flowchart of atmospheric correction for obtaining surface radiance and atmospheric downdraft radiation provided in an embodiment of the present invention; Figure 5 This is a flowchart of constructing an initial value dictionary and training an emissivity dictionary provided in an embodiment of the present invention; Figure 6 This is a flowchart of reconstructing the emissivity spectral vector based on the ROMP algorithm provided in an embodiment of the present invention; Figure 7 This is a flowchart of weight reallocation for the reconstructed emissivity spectral vector provided in an embodiment of the present invention; Figure 8 The flowchart illustrates the iterative solution used in embodiments of the present invention to recover the surface temperature and the noise-free and offset-free emissivity spectrum of the mixed pixels. Detailed Implementation

[0018] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. The components of the embodiments of the present invention described and shown in the accompanying drawings can generally be arranged and designed in various different configurations.

[0019] Example 1 As an example, to address the problems existing in the prior art, this embodiment provides a method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral hybrid pixel.

[0020] The implementation details of the method in this embodiment are described below. The following content is only for the convenience of understanding and is not necessary for implementing this solution.

[0021] The weighted dictionary method for calculating the temperature emissivity of mixed thermal infrared and hyperspectral pixels, as described in this embodiment, can be applied to electronic devices with communication, computing, and data storage capabilities. (See attached...) Figure 1 As shown, the method provided in this embodiment includes steps 110-150.

[0022] Step 110: Obtain simulation data and perform inversion channel selection based on the simulation data; the simulation data includes entrance pupil radiance, atmospheric transmittance, atmospheric upward radiation, atmospheric downward radiation and surface radiance.

[0023] For details, see attached. Figure 2 As shown, the process of acquiring simulation data and performing inversion channel filtering based on the simulation data includes: Multiple atmospheric conditions and underlying surface parameters are used as inputs; Hyperspectral simulation data were obtained using a radiative transfer model. The spectral response function is convolved with the hyperspectral simulation data to obtain the simulation data; Based on simulation data, channels that are limited by convolution distortion effects and have atmospheric transmittance below a set threshold are eliminated. Preferably, channels with atmospheric transmittance less than 0.5 are eliminated.

[0024] Channel availability analysis is achieved through channel inversion based on simulated data. The simulated data is obtained through simulation using mainstream radiative transfer models. These models require various atmospheric conditions and underlying surface parameters as input to comprehensively characterize the radiation observation scenario of the payload. In the thermal infrared band, underlying surface parameters generally refer to surface temperature and surface emissivity. Then, the spectral performance parameters of the thermal infrared hyperspectral payload are obtained, including the center wavelength and full width at half maximum (FWHM) of each channel. Based on these spectral performance parameters, a spectral response function (SRF) is established. The SRF is then used to convolve the thermal infrared simulated data to obtain the entrance pupil radiance, atmospheric transmittance, upward atmospheric radiation, downward atmospheric radiation, and surface radiance of each channel. Inversion based on the convolved simulated data identifies channels affected by convolution distortion (significant differences exist between spectral multiplication followed by convolution and convolution followed by multiplication), and these channels are then eliminated. Additionally, channels with atmospheric transmittance less than 0.5 also need to be eliminated.

[0025] Step 120: Obtain metadata and quality control identifiers for thermal infrared hyperspectral payload observations, and perform effective pixel screening to obtain the effective entrance pupil radiance of each channel of the thermal infrared hyperspectral payload.

[0026] As attached Figure 3 The diagram shows the flowchart for filtering valid pixels and radiometric correction based on metadata. Filtering valid pixels based on metadata first requires obtaining the metadata and quality control identifiers for a single observation. If the cloud cover identifier of a pixel is higher than 5% or the sensor signal-to-noise ratio identifier is less than the design specification of 95%, that pixel needs to be removed. Secondly, radiometric correction is performed on the remaining reliable pixels based on the radiometric correction parameters to obtain the entrance pupil radiance of each channel of the thermal infrared hyperspectral payload.

[0027] Specifically, effective pixel screening is performed to obtain the effective entrance pupil radiance of each channel of the thermal infrared hyperspectral payload. This includes: removing abnormal pixels based on metadata and quality control identifiers, and performing radiometric correction to obtain the entrance pupil radiance of each channel of the thermal infrared hyperspectral payload. Abnormal pixels are defined as pixels with a cloud coverage identifier higher than 5% or a sensor signal-to-noise ratio identifier less than 95% of the design specification.

[0028] Step 130: Perform atmospheric correction on the effective entrance pupil radiance to obtain the surface radiance and atmospheric downdraft radiation.

[0029] Specifically, atmospheric correction is performed on the effective entrance pupil radiance to obtain the surface radiance and downward atmospheric radiation, including: Atmospheric parameters at the observation time and within the observation range are obtained. The atmospheric parameters and the observation angle are input into the radiative transfer model to calculate atmospheric upward radiation, atmospheric transmittance and atmospheric downward radiation. In the embodiments of the present invention, atmospheric correction is carried out based on an atmospheric radiative transfer model, as shown in the appendix. Figure 4 The flowchart shown in this embodiment illustrates atmospheric correction of entrance pupil radiance. First, the coverage area, observation time, and observation angle provided by the metadata are obtained. Then, atmospheric profiles and other parameters within the observation range at that time are queried and downloaded. These atmospheric parameters and observation angle are used as inputs to the radiative transfer model to calculate the atmospheric upward radiation, atmospheric transmittance, and atmospheric downward radiation for each pixel and channel, thereby obtaining the surface radiance for each channel.

[0030] Based on the radiative transfer equation, surface radiance can be expressed as a coupled term of surface temperature, surface emissivity, and downward atmospheric radiation.

[0031] set up For real numbers, This represents the total number of channels in the thermal infrared hyperspectral image that can be used for inversion. Represents the surface radiance vector. The entrance pupil radiance vector, Represents the upward radiation vector of the atmosphere. This represents Hadema division, which involves dividing two vectors element by element. Let be the atmospheric transmittance vector, and let be the surface radiance of each channel: ; set up: This represents the surface emissivity spectral vector. This represents the Hadema product, which is the element-wise multiplication of two vectors. Represents the Planck function; The center wavelength vector of the channel. It refers to the surface temperature; Indicates length is A column vector whose elements are all 1; Here is the downward atmospheric radiation vector; the surface radiance satisfies the radiative transfer equation: .

[0032] Step 140: Construct an initial value dictionary. Select typical ground features and their emissivity spectra to train the initial value dictionary to obtain the emissivity dictionary.

[0033] As attached Figure 5 The flowchart shown illustrates the construction of the initial value dictionary and the training emission rate dictionary.

[0034] Specifically, constructing an initial value dictionary includes: Select several typical land cover emissivity spectra from the land cover spectral library; The emissivity spectra of typical ground features are convolved using the spectral response function of the thermal infrared hyperspectral channel to obtain the convolved spectral vector; Arrange the convolved spectral vectors in the column direction to form an initial value dictionary.

[0035] In practical applications, the spectral types of emissivity spectra of typical ground features include water bodies, snow and ice, vegetation, sandy areas, soil, rocks, and man-made objects. The spectral response function of the thermal infrared hyperspectral channel is used to... Convolving the spectra of typical ground features, the spectral vector of each convolved feature can be represented as follows: ,in, Index representing the spectrum, Represents the total number of spectra. The spectral vectors after convolution are arranged in the column direction to form an initial value dictionary. The emissivity spectrum after convolution is normalized by mean, subject to empirical constraints, and augmented with Gaussian white noise to construct the training data matrix. .

[0036] Specifically, the emissivity spectra of typical ground features are selected to train the initial value dictionary to obtain the emissivity dictionary. This includes: constructing a training data matrix, and alternately using the orthogonal matching pursuit algorithm and the K-singular value decomposition (KSVD) algorithm, with the training data matrix as input, to train the initial value dictionary.

[0037] In some optional embodiments, constructing a training data matrix includes: Let the spectral vector after each convolution be... ; Indicates the first One normalized spectral vector; Representing the An empirically constrained spectral vector; For mapping functions, it represents the operation of integrating normalized spectral vectors into empirically constrained spectral vectors; The empirical constraint matrix is ​​formed by concatenating empirical constraint spectral vectors along their column directions; the matrix and All belong to Through the empirical constraint matrix respectively By adding two different sets of white noise, we get: ; ; ; ; Specifically, the white noise has a mean of 0 and a variance of 10. -5 .

[0038] Training the initial value dictionary includes: set up: The number of atoms listed in the dictionary; This represents the total number of channels available for inversion in thermal infrared hyperspectral imaging. The dictionary to be trained; Represents the Frobenius norm; the training data matrix is ;No. indivual The sparse coefficient vector corresponding to the column vector sample in the data is represented as follows: ; The total number of spectra of ground features; The index of the spectral column vector is [index]; the sparse coefficient matrix is ​​[matrix]. , is represented as a matrix where sparse coefficient vectors are concatenated along their column directions, i.e. ; Represents the number of non-zero elements in a vector; Let be the sparsity of the dictionary training; then training on the initial dictionary can be represented as: .

[0039] Under normal circumstances, Set as ,in This indicates a rounding down operation. During the dictionary training iterations, when the dictionary to be trained is fixed at a certain moment, the OMP (Orthogonal Matching Pursuit) algorithm is used to estimate the sparse coefficients corresponding to the training samples; when the sparse coefficients are fixed at a certain moment, the KSVD algorithm is used to update the dictionary atoms; when the above optimization error is less than the default threshold or the number of iterations is greater than a set number (e.g., 15), the iteration stops, and the trained emissivity dictionary is output. .

[0040] This invention improves the training method of the emissivity dictionary. It not only directly integrates the initial convolutional emissivity spectrum into the initial value dictionary to provide a reference for subsequent weight estimation, but also uses mean normalization, empirical constraints and Gaussian noise augmentation techniques, combined with OMP and KSVD algorithms, to train an emissivity dictionary with spectral shape and value range information.

[0041] Step 150: Select the largest surface brightness temperature value as the initial value of surface temperature, calculate the surface emissivity spectral vector, and reconstruct the emissivity spectral vector based on the ROMP (Ratio Orthogonal Matching Pursuit) algorithm. Then, redistribute the weights of the reconstructed emissivity spectral vector through a weight allocation algorithm, and iteratively solve to recover the surface temperature and the noise-free and offset-free emissivity spectrum of the mixed pixels.

[0042] As attached Figure 6 As shown, the emissivity spectral vector is reconstructed based on the ROMP algorithm, including: A compressed observation matrix is ​​constructed based on the reconstructed emissivity spectrum vector containing noise and bias. Generate a random Gaussian matrix and perform QR decomposition to obtain an orthogonal matrix; The empirically constrained spectral vector is obtained by normalizing and empirically constraining the emissivity spectral vector containing noise and bias. The compressed vector is calculated based on the orthogonal matrix, the compressed observation matrix, and the empirically constrained spectral vector; the compressed dictionary is calculated based on the orthogonal matrix, the compressed observation matrix, and the trained emissivity dictionary. Solve for the sparse coefficients using the OMP algorithm; Reconstruct the normalized spectral vector based on the trained emissivity dictionary and sparse coefficients; Based on the reconstructed normalized spectral vector, the emissivity spectral vector is reconstructed.

[0043] Specifically, let: This represents the emissivity spectral vector that needs to be reconstructed and contains noise and bias. To compress the observation matrix, where, Indicates column direction Normalization; Indicates length is A row vector whose elements are all 1; It is a random Gaussian matrix. represent The orthogonal matrix obtained after QR decomposition To The empirically constrained spectral vector obtained after normalization and empirical constraints; for The compressed vector obtained after compressed observation; for A compressed dictionary is obtained after compressed observation; Let be the sparse coefficients obtained after OMP optimization, and Indicates the sparsity of the inversion; The reconstructed normalized spectral vector; This is the reconstructed emissivity spectrum vector; for The reverse operation; and They are respectively The maximum and minimum values ​​of are then used to solve for the sparse coefficients using the OMP algorithm, expressed as: ; ; ; The reconstructed emissivity spectral vector is represented as: .

[0044] This invention addresses the challenge of simultaneously retrieving surface temperature and emissivity from thermal infrared hyperspectral remote sensing data. It utilizes the compressed sensing theory of sparse dictionaries to sparsely represent emissivity, reducing the number of unknowns in the inversion equations. This satisfies the requirement that the number of unknowns is less than the number of equations, thereby enabling the simultaneous calculation of surface temperature and surface emissivity.

[0045] As attached Figure 7 As shown, the reconstructed emissivity spectral vector is weighted and redistributed to recover the noise-free and offset-free emissivity spectrum of the mixed pixels, including: Using an emissivity spectral vector that incorporates various ground features, noise, and biases as input, a constrained quadratic programming optimization method is used to solve for the preliminary weight estimates of the emissivity spectral vector. The candidate dictionary is constructed by selecting emissivity spectral vectors from the initial value dictionary, and calculating the residual vectors based on the emissivity spectral vectors that are mixed with various land cover information, noise and bias, and the emissivity spectral vectors that are re-estimated based on weights and candidate dictionaries. The candidate dictionary after residual allocation is then constructed. The ROMP algorithm is used to reconstruct the emissivity spectrum of each candidate in the candidate dictionary after residual allocation, and the weights of the reconstructed emissivity spectrum are redistributed to restore the noise-free and offset-free emissivity spectrum of the mixed pixels.

[0046] The method of this invention combines the ROMP algorithm and the weight allocation algorithm into the WOMP (Weighted Orthogonal Matching Pursuit) algorithm. set up These are preliminary estimates of the weights. For the weight items to be optimized, for The square of the norm, It is a parameter used to adjust the strength of the centralization regularization term; Defined as The maximum value of the candidate spectrum, generally speaking, Set to 7, representing the seven land cover types: water, ice and snow, vegetation, soil, rocks, sand, and man-made structures; Represents the first part of a vector The index of the maximum value; This represents a list of candidate spectral indices, where the weight value of the candidate spectrum corresponding to each index in the list is greater than 10. -2 And placed in front One; order This represents the number of actual candidate spectra. for Preliminary estimate of the weights corresponding to the index. The normalized weight vector has the following preliminary weight estimates: ; ; ; .

[0047] set up Indicates in Based on The emissivity spectral vector selected by the index. express The candidate dictionary formed; This represents the re-estimated emissivity spectral vector based on weights and candidate dictionaries. For the residual vector, The candidate dictionary after residual allocation is then The constructed candidate dictionary is represented as follows: ; The residual vector is represented as: ; The candidate dictionary after residual allocation is represented as follows: .

[0048] Let: the candidate emissivity spectrum vector reconstructed by the ROMP algorithm; It is A set arranged in column direction; This represents the emissivity spectrum finally recovered using the WOMP algorithm; for For each candidate spectrum in the dataset, the ROMP algorithm can independently reconstruct it. All reconstructed spectra can then be weighted and redistributed to recover the noise-free and offset emissivity spectra of the mixed pixels. This process is represented as follows: ; .

[0049] ; This represents the inverse operation of WOMP to recover the emissivity spectrum. In the solution, the normalized weight vector... The variables are unknowns, and the coefficient matrix is ​​sparse. The non-zero elements in the variable are also unknown variables, among which, Each column of sparse coefficient vectors in the equation can be reconstructed into a candidate emissivity spectrum. To ensure that the number of unknowns is less than the number of equations, the inversion sparsity is generally set to 1. Since the analytical solutions to the above unknowns cannot be obtained directly, iterative solutions are needed to solve for the unknowns and simultaneously recover the surface temperature and surface emissivity spectra.

[0050] The radiative transfer equation is expressed as: ; As attached Figure 8 As shown, the iterative solution recovers the surface temperature and the noise-free and offset-free emissivity spectrum of the mixed pixels, including: Based on the Planck function, the atmospherically corrected surface radiance is converted into surface brightness temperature, and the largest surface brightness temperature value among several channels is selected as the initial value of surface temperature. The surface emissivity spectral vector is directly calculated based on the current estimated surface temperature and the radiative transfer equation. Using the surface emissivity spectral vector as input, the emissivity spectral vector is reconstructed based on the ROMP algorithm, and the weights of the reconstructed emissivity spectral vector are redistributed using a weight allocation algorithm to recover the unbiased emissivity spectral vector. Establish a cost function and use the golden ratio method in Find the temperature value that minimizes the cost function within the interval, and use it as the optimized surface temperature value for this iteration; This indicates the lower limit of the temperature search below the initial surface temperature of 1K. This indicates the upper limit for temperature search that is 1K higher than the initial surface temperature.

[0051] Linear fitting is performed on the cost function values ​​in the last few iterations. If the absolute value of the slope of the fitted line is less than or the number of iterations exceeds the set number, the iteration terminates, and the optimized surface temperature value is output as the final inversion value of the surface temperature, and the unbiased emissivity spectral vector is output as the final inversion value of the surface emissivity spectrum; otherwise, the iteration is repeated.

[0052] Let the cost function be... The cost function is expressed as: .

[0053] For example, if a linear fit is performed on the cost function in the last three iterations, and the slope of the fitted line is less than 10... -5 If the number of iterations exceeds 15, the iteration terminates and the output is... As the final inversion value of surface temperature, the output is... As the final inversion value of the surface emissivity spectrum, otherwise, let The iteration is then repeated. This embodiment of the invention utilizes the local golden section method to optimize the surface temperature inversion value and improves the iteration convergence condition. Specifically, it determines the iteration to proceed by judging whether the slope of the last three cost function values ​​is gentle, thereby avoiding the occurrence of local false convergence.

[0054] Traditional OMP algorithms exhibit significant greediness in reconstructing emissivity spectra using sparse dictionaries, meaning the retrieved emissivity spectra tend to favor emissivity values ​​of a single class in the spectral library, making them unsuitable for surface inversion of mixed pixels. This invention overcomes the greediness of traditional OMP algorithms by combining the ROMP algorithm and a weight allocation algorithm into the WOMP algorithm. Furthermore, by employing an initial value dictionary and an emission dictionary, it reconstructs unbiased emissivity spectra, thus meeting the requirement of simultaneous inversion of surface temperature and surface emissivity spectra for mixed pixels. The WOMP algorithm provides class weight information for the emissivity spectrum, indicating the probability that the emissivity is biased towards a certain land cover class, thereby providing direct discriminative information for subsequent decoupling of mixed pixels.

[0055] The above are merely preferred embodiments of the present invention and are not intended to limit the present invention. Various modifications and variations can be made to the present invention by those skilled in the art. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral hybrid pixel, characterized in that, include: Acquire simulation data and perform inversion channel filtering based on the simulation data; The simulation data includes entrance pupil radiance, atmospheric transmittance, upward atmospheric radiation, downward atmospheric radiation, and surface radiance. Obtain metadata and quality control identifiers for thermal infrared hyperspectral payload observations, and perform effective pixel screening to obtain the effective entrance pupil radiance of each channel of the thermal infrared hyperspectral payload; Atmospheric correction is performed on the effective entrance pupil radiance to obtain the surface radiance and atmospheric downdraft radiation; An initial value dictionary is constructed, and the emissivity spectra of typical ground features are selected to train the initial value dictionary to obtain the emissivity dictionary; The largest surface brightness temperature value is selected as the initial value of surface temperature. The surface emissivity spectral vector is calculated and reconstructed based on the ROMP algorithm. The weights of the reconstructed emissivity spectral vector are then redistributed using a weight allocation algorithm. The surface temperature and the noise-free and offset-free emissivity spectra of the mixed pixels are then recovered through iterative solution.

2. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, Acquire simulation data and perform inversion channel filtering based on the simulation data, including: Multiple atmospheric conditions and underlying surface parameters are used as inputs; Hyperspectral simulation data were obtained using a radiative transfer model. The spectral response function is convolved with the hyperspectral simulation data to obtain the simulation data; Based on simulation data, channels that are limited by convolution distortion effects and atmospheric transmittance less than a set threshold are eliminated.

3. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, Effective pixel screening is performed to obtain the effective entrance pupil radiance of each channel of the thermal infrared hyperspectral load. This includes: removing abnormal pixels based on metadata and quality control identifiers, and performing radiometric correction to obtain the entrance pupil radiance of each channel of the thermal infrared hyperspectral load.

4. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, Atmospheric correction is performed on the effective entrance pupil radiance to obtain the surface radiance and downward atmospheric radiation, including: Atmospheric parameters at the observation time and within the observation range are obtained. The atmospheric parameters and the observation angle are input into the radiative transfer model to calculate atmospheric upward radiation, atmospheric transmittance and atmospheric downward radiation. set up For real numbers, This represents the total number of channels in the thermal infrared hyperspectral image that can be used for inversion. Represents the surface radiance vector. The entrance pupil radiance vector, Represents the upward radiation vector of the atmosphere. This represents Hadema division, which involves dividing two vectors element by element. Given the atmospheric transmittance vector, the surface radiance of each channel is expressed as: ； set up: This represents the surface emissivity spectral vector. This represents the Hadema product, which is the element-wise multiplication of two vectors. Represents the Planck function; The center wavelength vector of the channel. It refers to the surface temperature; Indicates length is A column vector whose elements are all 1; Here is the downward atmospheric radiation vector; the surface radiance satisfies the radiative transfer equation: 。 5. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, Construct an initial value dictionary, including: Select several typical land cover emissivity spectra from the land cover spectral library; The emissivity spectra of typical ground features are convolved using the spectral response function of the thermal infrared hyperspectral channel to obtain the convolved spectral vector; Arrange the convolved spectral vectors in the column direction to form an initial value dictionary.

6. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, The initial value dictionary is trained by selecting the emissivity spectra of typical ground features to obtain the emissivity dictionary. This process includes: constructing a training data matrix, and alternately using the orthogonal matching pursuit algorithm and the K-singular value decomposition algorithm, with the training data matrix as input, to train the initial value dictionary.

7. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 6, characterized in that, Construct the training data matrix, including: Let the spectral vector after each convolution be... ; Indicates the first One normalized spectral vector; Representing the An empirically constrained spectral vector; The mapping function represents the operation of integrating normalized spectral vectors into empirically constrained spectral vectors; The empirical constraint matrix is ​​formed by concatenating empirical constraint spectral vectors along their column directions; the matrix and All belong to Through the empirical constraint matrix respectively By adding two different sets of white noise, we get: ； ； ； ； Training the initial value dictionary includes: set up: The number of atoms listed in the dictionary; This represents the total number of channels available for inversion in thermal infrared hyperspectral imaging. This represents the dictionary to be trained; Represents the Frobenius norm; the training data matrix is ;No. indivual The sparse coefficient vector corresponding to the column vector sample in the data is represented as follows: ; The total number of spectra of ground features; The index of the spectral column vector is [index]; the sparse coefficient matrix is ​​[matrix]. , is represented as a matrix where sparse coefficient vectors are concatenated along their column directions, i.e. ; Represents the number of non-zero elements in a vector; Let be the sparsity of the dictionary training; then training the initial dictionary can be represented as: 。 8. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, The emissivity spectral vector is reconstructed based on the ROMP algorithm, including: A compressed observation matrix is ​​constructed based on the reconstructed emissivity spectrum vector containing noise and bias. Generate a random Gaussian matrix and perform QR decomposition to obtain an orthogonal matrix; The empirically constrained spectral vector is obtained by normalizing and empirically constraining the emissivity spectral vector containing noise and bias. The compressed vector is calculated based on the orthogonal matrix, the compressed observation matrix, and the empirically constrained spectral vector; the compressed dictionary is calculated based on the orthogonal matrix, the compressed observation matrix, and the trained emissivity dictionary. Solve for the sparse coefficients using the OMP algorithm; Reconstruct the normalized spectral vector based on the trained emissivity dictionary and sparse coefficients; Based on the reconstructed normalized spectral vector, the emissivity spectral vector is reconstructed.

9. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, The reconstructed emissivity spectral vector is weighted and redistributed to recover the noise-free and offset-free emissivity spectrum of the mixed pixels, including: Using an emissivity spectral vector that incorporates various ground features, noise, and biases as input, a constrained quadratic programming optimization method is used to solve for the preliminary weight estimates of the emissivity spectral vector. The candidate dictionary is constructed by selecting emissivity spectral vectors from the initial value dictionary, and calculating the residual vectors based on the emissivity spectral vectors that are mixed with various land cover information, noise and bias, and the emissivity spectral vectors that are re-estimated based on weights and candidate dictionaries. The candidate dictionary after residual allocation is then constructed. The ROMP algorithm is used to reconstruct the emissivity spectrum of each candidate in the candidate dictionary after residual allocation, and the weights of the reconstructed emissivity spectrum are redistributed to restore the noise-free and offset-free emissivity spectrum of the mixed pixels.

10. The method for calculating the temperature emissivity of a weighted dictionary thermal infrared hyperspectral mixed pixel according to claim 1, characterized in that, Iterative solutions recover the surface temperature and the noise-free and offset-free emissivity spectra of the mixed pixels, including: Based on the Planck function, the atmospherically corrected surface radiance is converted into surface brightness temperature, and the largest surface brightness temperature value among several channels is selected as the initial value of surface temperature. The surface emissivity spectral vector is directly calculated based on the current estimated surface temperature and the radiative transfer equation. Using the surface emissivity spectral vector as input, the emissivity spectral vector is reconstructed based on the ROMP algorithm, and the weights of the reconstructed emissivity spectral vector are redistributed using a weight allocation algorithm to recover the unbiased emissivity spectral vector. Establish a cost function and use the golden section method to find the temperature value corresponding to the minimum cost function, which is used as the optimized surface temperature value for this iteration. Linear fitting is performed on the cost function values ​​in the last few iterations. If the absolute value of the slope of the fitted line is less than or the number of iterations exceeds the set number, the iteration terminates, and the optimized surface temperature value is output as the final inversion value of the surface temperature, and the unbiased emissivity spectral vector is output as the final inversion value of the surface emissivity spectrum; otherwise, the iteration is repeated.

Images