Three-dimensional medical volume data zero watermark method based on sixteen-element discrete cosine transform

Abstract

The invention discloses a three-dimensional medical volume data zero-watermarking method based on sixteen-element discrete cosine transform, which comprises the following steps of: slicing and grouping original three-dimensional medical volume data, and then performing Tucker decomposition on each medical image slice group to obtain a first energy characteristic graph; performing two-dimensional non-downsampling discrete wavelet transform on the first energy feature maps obtained in each group to obtain corresponding low-frequency sub-bands, performing non-overlapping partitioning on the low-frequency sub-bands in each group and constructing sixteen-element image blocks, calculating direct-current components of the image blocks after sixteen-element discrete cosine transform in a space domain, and performing discrete cosine transform on the direct-current components to obtain a first energy feature map; and finally, obtaining a binary image feature matrix according to the L2 norm of each direct-current component and applying the binary image feature matrix to zero-watermark authentication. The method not only has strong robustness for various conventional attacks in digital watermark attacks, but also meets the requirement of high fidelity during image data copyright protection in the medical field, and has high security at the same time.

Description

Technical Field

[0001] The embodiments of the present invention relate to the field of image processing technology, and in particular to a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform. Technical Background

[0002] With the rapid development of multimedia technology and network applications in recent years, digitization has become popular in the medical field. Applying computer information network technology to the entire medical process has gradually become a new model of modern medical care, and is also the development direction and management goal of public medical care. Some new medical technologies are also developing rapidly and are widely used in medical practice, such as digital medical testing, remote diagnosis and treatment, and medical digital imaging technology, which have accelerated the process of medical digitization. The medical field is rapidly transforming from traditional medicine to telemedicine, and a large amount of medical information can be quickly and efficiently transmitted on the public network, including a large number of medical images.

[0003] At present, medical devices such as computed tomography (CT) and magnetic resonance imaging (MRI) output three-dimensional medical volume data, which inevitably face the risk of tampering and theft during transmission. Moreover, in the medical field, any slight change may indicate a certain lesion, so any change to the original image may cause distortion of the medical image, thus affecting the doctor's diagnosis. In order to solve these problems, we use zero watermark technology to protect the copyright information of medical volume data. Zero watermark technology is usually a method of generating a unique key by performing an XOR operation on the extracted original image features and the binary watermark image to achieve the purpose of watermark embedding, so that the pixel information in the original image will not be changed. This technology not only balances the imperceptibility and robustness of traditional embedded watermarks, but also does not cause any loss to the original image, so it will not affect the normal use value of the original image. Zero watermark technology is one of the important research directions in the field of information hiding.

[0004] However, there are relatively few zero-watermark algorithms currently used in medical volume data, and the algorithms are not very robust to conventional image attacks. With the increase in the transmission of medical information on the Internet, it has become particularly important and urgent to improve the robustness of zero-watermark algorithms for medical volume data to reduce the risk of tampering and theft of its content. Summary of the invention

[0005] The embodiment of the present invention provides a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform, which can effectively resist various conventional image processing attacks to protect the copyright information of three-dimensional medical volume data and solve the problem of poor robustness of existing three-dimensional medical volume data zero watermark algorithms.

[0006] The embodiment of the present invention provides a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform, comprising:

[0007] Embed zero watermark in original three-dimensional medical volume data and extract zero watermark from the three-dimensional medical volume data to be authenticated;

[0008] The embedding of a zero watermark in the original three-dimensional medical volume data comprises:

[0009] The Fibonacci Q matrix is ​​used as the scrambling matrix to scramble the original binary watermark image W of size N×N to obtain the scrambled binary watermark image W A ;

[0010] The original three-dimensional medical volume data V of size M×M×P is sliced ​​and grouped, and the processed medical image slice group is recorded as U i , i = 1, 2, 3, ..., 15; for each medical image slice group U i Tucker decomposition is performed on all of them, and the first energy feature graph after Tucker decomposition of each group is extracted, marked as E i , i=1,2,3,...,15; for E i Each component in is subjected to k-layer two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i , i=1,2,3,……,15;

[0011] For the low frequency subband L i Each component in is divided into non-overlapping sub-blocks, and then L i The sub-blocks at the corresponding positions in each component in the hexadecimal image block are used to construct the hexadecimal image block and the DC coefficient D of the image block after the hexadecimal discrete cosine transform is calculated in the spatial domain. st , D st The magnitudes of the real and imaginary parts of st Indicates that, finally using O st The L2 norm of is used to obtain the binary robust feature matrix T;

[0012] The binary image feature matrix J is obtained by scrambling the binary robust feature matrix T using the two-dimensional Sinusoidal-Chebyshev chaotic system and the Spatial Pyramid Matching chaotic system.

[0013] The scrambled binary watermark image W corresponding to the original binary watermark image W A Perform an XOR operation with the binary image feature matrix J to obtain the authenticated zero watermark image H, save H to the watermark database of the third-party registration agency, and save the control parameters and process parameters used in the binary zero watermark embedding process as the key;

[0014] The zero-watermark extraction of the three-dimensional medical volume data to be authenticated comprises:

[0015] The three-dimensional medical volume data V′ with a size of M×M×P and to be authenticated is sliced ​​and grouped, and the processed medical image slice group is recorded as U i ′, i = 1, 2, 3, ..., 15; for each medical image slice group U i ′ are all subjected to Tucker decomposition, and the first energy feature graph after Tucker decomposition of each group is extracted, marked as E i ′, i=1,2,3,...,15; for E i Each component in ′ is transformed by k layers of two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i ′, i = 1, 2, 3, ..., 15;

[0016] For the low frequency subband L i Each component in ′ is divided into non-overlapping sub-blocks, and then L i The sub-blocks at the corresponding positions in each component of ′ construct the hexadecimal image block and calculate the DC coefficient D of the image block after the hexadecimal discrete cosine transform in the spatial domain st ′, D st The magnitudes of the real and imaginary parts of ′ are expressed by a one-dimensional vector O st ' indicates that finally O st The L2 norm of ′ is used to obtain the binary robust feature matrix T′;

[0017] The binary image feature matrix J′ is obtained by scrambling the binary robust feature matrix T′ using the two-dimensional Sinusoidal-Chebyshev chaotic system and the Spatial Pyramid Matching chaotic system.

[0018] Take out the certified zero watermark image H stored in the watermark database of the third-party registration agency, perform an XOR operation with the binary image feature matrix J′, and obtain the binary watermark image W to be descrambled. A ′;

[0019] The Fibonacci Q inverse matrix is ​​used as the anti-scrambling matrix to treat the binary watermark image W for the anti-scrambling operation. A ′ is descrambled to extract the final watermark information W′, and finally the copyright of the three-dimensional medical volume data V′ to be authenticated is determined according to the content information displayed by W′.

[0020] The embodiment of the present invention provides a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform, which uses the Fibonacci Q matrix as the scrambling matrix to scramble the original watermark image, and uses the two-dimensional Sinusoidal-Chebyshev chaotic system and the Spatial Pyramid Matching chaotic system to encrypt and scramble the feature matrix, thereby improving the security of the zero watermark algorithm. Tucker decomposition and hexadecimal discrete cosine transform are used to construct a zero watermark by utilizing the correlation between each slice of three-dimensional medical volume data and its low-frequency components, and by directly calculating the DC coefficient of the hexadecimal image block after the hexadecimal discrete cosine transform in the spatial domain, the operation cost is effectively reduced and the timeliness of the algorithm is enhanced. At the same time, before the zero watermark is constructed, the low-frequency component is extracted by using the two-dimensional non-subsampled discrete wavelet transform, which effectively improves the robustness and watermark embedding capacity of the method. BRIEF DESCRIPTION OF THE DRAWINGS

[0021] In order to more clearly explain the technical solutions in the embodiments of the present invention, the drawings required in the description of the embodiments are briefly described below. Obviously, the drawings described below are only part of the embodiments of the present invention, and ordinary technicians in this field can obtain other drawings based on these drawings without creative work.

[0022] Figure 1 It is a flow chart of an embedding method in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention;

[0023] Figure 2 It is a flow chart of an extraction method in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention;

[0024] Figure 3 It is a schematic diagram of a single slice of MRI and CT brain three-dimensional medical volume data in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention;

[0025] Figure 4 It is a schematic diagram of an original binary watermark image in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention;

[0026] Figure 5 It is a test NC value table of MRI and CT brain three-dimensional medical volume data under various non-geometric attacks in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention;

[0027] Figure 6It is a test NC value table of MRI and CT brain three-dimensional medical volume data under multiple geometric attacks in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention. DETAILED DESCRIPTION

[0028] In order to enable those skilled in the art to better understand the scheme of the present invention, the technical scheme in the embodiment of the present invention will be described in detail below in conjunction with the drawings in the embodiment of the present invention. Obviously, the described embodiment is only a part of the embodiment of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work should fall within the scope of protection of the present invention.

[0029] Figure 1 is a flow chart of an embedding method in a three-dimensional medical volume data zero watermarking method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention, Figure 2 The present invention provides a flowchart of a method for extracting a zero watermark of three-dimensional medical volume data based on hexadecimal discrete cosine transform provided by an embodiment of the present invention. This embodiment is applicable to the case where a zero watermark is used to protect the copyright of three-dimensional medical volume data.

[0030] like Figure 1 As shown, the present embodiment provides an embedding method in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform, which specifically includes the following steps:

[0031] S110, using the Fibonacci Q matrix as the scrambling matrix to perform a scrambling operation on the original binary watermark image W of size N×N, to obtain a scrambled binary watermark image W A .

[0032] Figure 4 Schematic diagram of an original binary watermark image in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention, with a size of N×N. In this embodiment, the original binary watermark image provided above can be used to implement binary watermark scrambling processing.

[0033] The Fibonacci Q matrix is ​​used as the scrambling matrix to perform the scrambling operation on the original binary watermark image W, which is expressed by the following formula:

[0034]

[0035] In the formula, Q p is the p-order Fibonacci Q matrix, i.e., the scrambling matrix, N is the length and width of the original binary watermark image W, x and y are the position coordinates of the pixels in the original binary watermark image W before scrambling, and x′ and y′ are the scrambled binary watermark image W after scrambling.A The position coordinates of the pixels in , mod(·) is the remainder function.

[0036] The p-order Fibonacci Q matrix is ​​implemented by the following formula:

[0037]

[0038] In the formula, F p+1 , F p , F p-1 are the number of terms in the Fibonacci sequence. The order p is recorded as the key Key1 and saved.

[0039] In this embodiment, the Fibonacci Q matrix is ​​used as the scrambling matrix to perform a scrambling operation on the original binary watermark image W, which has low computational complexity and high chaos degree, and can enhance the timeliness and security of the overall watermark algorithm.

[0040] S120, slicing and grouping the original three-dimensional medical volume data V of size M×M×P, and performing grouping on each medical image slice group U i Tucker decomposition is performed on all groups, and the first energy feature graph E after Tucker decomposition is extracted. i , for E i Each component in is subjected to k-layer two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i .

[0041] For example, the original three-dimensional medical volume data V of size M×M×P is sliced ​​and divided into 15 medical image slice groups, denoted as U i , i=1,2,3,……,15, and is implemented as follows:

[0042]

[0043] Where P is the total number of slices after slicing of the three-dimensional medical volume data, floor(·) is the function of rounding down to an integer, and mod(·) is the function of taking the remainder.

[0044] Will U i Converted into tensor form and subjected to Tucker decomposition, the medical image slice group U i It can be expressed as:

[0045] U i =G i ×1V i (1) ×2V i (2) ×3V i (3)

[0046] In the formula, i=1,2,3,...,15, U i is the tensor representation of each medical image slice group, G i is the core tensor of each group, V (1)i , V (2)i , V (3)i is the factor matrix of each group of tensors in different directions.

[0047] Let tensor X i =G i ×1V i (1) ×2V i (2) is the tensor U i , extract each group of sub-tensors X i The first sub-image of the medical image slice group U i The first energy characteristic graph after Tucker decomposition is denoted as E i , i=1,2,3,……,15.

[0048] For E i Perform k-layer two-dimensional non-subsampled discrete wavelet transform to obtain E i The low-frequency subband L corresponding to each component in i , i=1,2,3,……,15.

[0049] Three-dimensional medical volume data is a typical high-dimensional medical image data, which usually needs to be sliced ​​when processed. By slicing and grouping the three-dimensional medical volume data, and then performing Tucker decomposition and extracting the first energy feature map, the correlation between each slice in each medical image slice group can be effectively preserved, and most of the energy of the three-dimensional medical volume data is retained as a whole, which helps to enhance the overall robustness of the medical volume data zero watermark algorithm. At the same time, the use of grouping operations combined with Tucker decomposition can improve the timeliness of the algorithm without affecting the robustness, and avoid processing all slices. The low-frequency subband image obtained by using two-dimensional non-subsampled discrete wavelet transform can preserve most of the information in the first energy feature map after Tucker decomposition without reducing the size. Based on this low-frequency subband, the robust feature matrix is ​​constructed, which can enhance the robustness of the algorithm while ensuring the watermark embedding capacity.

[0050] S130, for low frequency sub-band L i Each component in is divided into non-overlapping sub-blocks, and then L i The sub-blocks corresponding to the positions of each component in the image block are used to construct the hexadecimal image block and the DC coefficient D of the image block after the hexadecimal discrete cosine transform is calculated in the spatial domain. st , D st The magnitudes of the real and imaginary parts ofst Indicates that, finally using O st The L2 norm of is used to obtain the binary robust feature matrix T.

[0051] The pair L i Divide the non-overlapping sub-blocks into r×r blocks, and then use L i The sub-blocks corresponding to the positions of each component in the image block are used to construct the hexadecimal image block and the DC coefficient D of the image block after the hexadecimal discrete cosine transform is calculated in the spatial domain. st , which is implemented as follows:

[0052]

[0053] In the formula, u g is a unit pure hexadecimal number that satisfies

[0054] L st (x, y) is the L after block division i The hexadecimal image block constructed by the sub-blocks corresponding to the positions of each component in is expressed as:

[0055]

[0056] C st (u,v) is L st The spectrum of (x,y) in the hexadecimal discrete cosine transform domain can be expressed as:

[0057]

[0058] Where, s = 1, 2, 3, ..., M / r, t = 1, 2, 3, ..., M / r, C st Distribution of (u,v) in the hexadecimal space.

[0059] From the above formula, we can see that when u and v are both equal to 0, that is, the hexadecimal image block L st (x,y) DC coefficient D after hexadecimal discrete cosine transform st It can be expressed in the airspace as:

[0060]

[0061] Then by L st The definition of the (x,y) hexadecand is:

[0062]

[0063]

[0064] Other relevant parameters are:

[0065]

[0066] After the above calculation, we can get the sixteen-dimensional image block L st (x,y) DC coefficient D after hexadecimal discrete cosine transform st The real part and the imaginary part, s = 1, 2, 3, ..., M / r, t = 1, 2, 3, ..., M / r.

[0067] Using the one-dimensional vector O st Indicates D st The magnitude of the real and imaginary parts of :

[0068]

[0069] The last st The L2 norm of the coefficient feature matrix S is obtained, which is implemented using the following formula:

[0070] S(s,t)=||O st ||2

[0071] Wherein, ||·||2 is the calculation symbol of L2 norm, s=1,2,3,……,M / r, t=1,2,3,……,M / r.

[0072] The mean value of the coefficient feature matrix S is used to binarize the matrix to obtain the binary robust feature matrix T, which is implemented using the following formula:

[0073]

[0074] Using the low frequency subband L i The construction of the hexadecimal image block by the sub-blocks at the corresponding positions of each component after non-overlapping division can retain the correlation between the low-frequency components to the greatest extent. At the same time, directly calculating the DC coefficients of the hexadecimal image block after the hexadecimal discrete cosine transform in the spatial domain can effectively reduce the computational cost and improve the timeliness of the algorithm. The binary robust feature matrix is ​​constructed based on the robustness of the relationship between the L2 norm and the overall mean of the DC coefficients of each hexadecimal image block after the hexadecimal discrete cosine transform, which also improves the robustness of the algorithm.

[0075] S140. Using a two-dimensional Sinusoidal-Chebyshev chaotic system and a Spatial Pyramid Matching chaotic system, a scrambling operation is performed on the binary robust feature matrix T to obtain a binary image feature matrix J.

[0076] The binary robust feature matrix T is expanded into a binary extended feature matrix R of size N×N, which is implemented as follows:

[0077]

[0078] Where N is the length and width of the original binary watermark image W, M is the length and width of a single slice in the original three-dimensional medical volume data V, r is the length and width of the sub-block, and repmat(·) is the repeated array copy function.

[0079] A two-dimensional Sinusoidal-Chebyshev system is used to generate a one-dimensional chaotic sequence of length N×N and the mean of the chaotic sequence is used for binarization to obtain a one-dimensional binary sequence Y. The reshape function is then used to reconstruct Y into a two-dimensional matrix and XORed with the binary extended feature matrix R to obtain the encrypted binary feature matrix B. The implementation method is as follows:

[0080] B = xor(R, reshape(Y, [NN]))

[0081] Where xor(·) is the function for performing XOR calculation, and reshape(·) is the function for reconstructing array.

[0082] The two-dimensional Sinusoidal-Chebyshev system is expressed by the following formula:

[0083]

[0084] Among them, a, b, c, d are system control parameters, m i , n i is the intermediate variable in the iterative process of the two-dimensional Sinusoidal-Chebyshev chaotic system, m i+1 , n i+1 m i , n i The next state of , mod(·) is the remainder function.

[0085] The Spatial Pyramid Matching chaotic system is used to generate a one-dimensional chaotic sequence Z of length N×N. The index vector of the one-dimensional chaotic sequence Z in descending order is used to scramble the encrypted binary feature matrix B to obtain the binary image feature matrix J. The following method is used:

[0086]

[0087] Where sort(·) is the sorting function, and index is the index vector obtained by sorting the one-dimensional chaotic sequence Z in descending order.

[0088] The Spatial Pyramid Matching chaotic system is expressed by the following formula:

[0089]

[0090] Where η and μ are system control parameters, h is a fixed parameter, q(t) is the intermediate variable in the iterative process of the Spatial Pyramid Matching chaotic system, and q(t+1) is the next state of q(t).

[0091] When using the two-dimensional Sinusoidal-Chebyshev chaotic system to obtain a one-dimensional binary sequence Y with a length of N×N, the system control parameters a, b, c, d and the system initial values ​​m0, n0 are used as the key Key2. When using the Spatial PyramidMatching chaotic system to generate a one-dimensional chaotic sequence Z with a length of N×N, the fixed parameter h, the system control parameters η and μ and the system initial value q(0) are used as the key Key3. Key2 and Key3 are saved.

[0092] The two-dimensional Sinusoidal-Chebyshev chaotic system and Spatial Pyramid Matching chaotic system are used to encrypt and scramble the binary robust feature matrix T, which can not only improve the security of the algorithm, but also improve the robustness of the algorithm through the random characteristics of the generated binary sequence. The binary robust feature matrix T is expanded by the repmat function, which solves the problem of size mismatch between the feature matrix and the scrambled watermark when XORing, and improves the overall stability of the algorithm.

[0093] S150: The scrambled binary watermark image W corresponding to the original binary watermark image W A Perform an XOR operation with the binary image feature matrix J to obtain the authenticated zero watermark image H, save H to the watermark database of the third-party registration agency, and save the control parameters and process parameters used in the binary zero watermark embedding process as the key, thus completing the zero watermark embedding process.

[0094] H = xor(W A ,J)

[0095] Where xor(·) is a function that performs XOR calculation.

[0096] Through the above steps, the final authentication zero-watermark image H and the related keys Key1, Key2 and Key3 used in the zero-watermark embedding process are obtained for the final judgment of the copyright ownership of the three-dimensional medical volume data to be authenticated.

[0097] In addition, this embodiment also provides an extraction method in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform, such as Figure 2 As shown, the specific steps include:

[0098] S210, slicing and grouping the three-dimensional medical volume data V′ of size M×M×P to be authenticated, and performing grouping on each medical image slice group U i ′ are all subjected to Tucker decomposition, and the first energy feature graph E after Tucker decomposition of each group is extracted i ′, for E i Each component in ′ is subjected to k-layer two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i ′.

[0099] The 3D medical volume data V′ with a size of M×M×P and to be authenticated is sliced ​​and divided into 15 medical image slice groups, denoted as U i ′, i = 1, 2, 3, ..., 15, and is implemented as follows:

[0100]

[0101] Will U i ′ is converted into tensor form and subjected to Tucker decomposition, then the medical image slice group U i ′ can be expressed as:

[0102]

[0103] Let tensor is the tensor U i ′, extract each group of sub-tensors X i ′, namely the medical image slice group U i ′The first energy characteristic graph after Tucker decomposition, denoted as E i ′, i=1,2,3,……,15.

[0104] For E i ′ Perform k-layer two-dimensional non-subsampled discrete wavelet transform to obtain E i The low-frequency subband L corresponding to each component in ′ i ′, i=1,2,3,……,15.

[0105] S220, for low frequency sub-band L i 'Divide into non-overlapping sub-blocks, and then use the L i The sub-blocks corresponding to the positions of each component in ′ construct a hexadecimal image block and calculate the DC coefficient D′ of the image block after the hexadecimal discrete cosine transform in the spatial domain st , D′ st The magnitude of the real part and each imaginary part is expressed by a one-dimensional vector O′ st Indicates that, finally using O′ st The L2 norm of is used to obtain the binary robust feature matrix T′.

[0106] For L i ′ divides the non-overlapping sub-blocks into r×r blocks, and then uses L i The sub-blocks corresponding to the positions of each component in ′ construct a hexadecimal image block and calculate the DC coefficient D′ of the image block after the hexadecimal discrete cosine transform in the spatial domain st , D′ st The magnitude of the real part and each imaginary part is expressed by a one-dimensional vector O′ st Indicates that, finally using O′ st The L2 norm of the coefficient characteristic matrix S' is obtained, which is implemented using the following formula:

[0107] S'(s,t)=||O' st ||2

[0108] The mean value of the coefficient feature matrix S′ is used to binarize the matrix to obtain the binary robust feature matrix T′, which is implemented using the following formula:

[0109]

[0110] S230, using a two-dimensional Sinusoidal-Chebyshev chaotic system and a Spatial Pyramid Matching chaotic system to perform a scrambling operation on the binary robust feature matrix T′ to obtain a binary image feature matrix J′.

[0111] The binary robust feature matrix T′ is expanded into a binary extended feature matrix R′ of size N×N, which is implemented as follows:

[0112]

[0113] Using the key Key2, a two-dimensional Sinusoidal-Chebyshev system is used to generate a one-dimensional chaotic sequence of length N×N and the mean of the chaotic sequence is used for binarization to obtain a one-dimensional binary sequence Y. Then, the reshape function is used to reconstruct Y into a two-dimensional matrix and perform an XOR operation with the binary extended feature matrix R′ to obtain the encrypted binary feature matrix B′. The implementation method is as follows:

[0114] B'=xor(R',reshape(Y,[NN]))

[0115] Using the key Key3, the Spatial Pyramid Matching chaotic system is used to generate a one-dimensional chaotic sequence Z, and the index vector of the one-dimensional chaotic sequence Z in descending order is used to scramble the encrypted binary feature matrix B′ to obtain the binary image feature matrix J′. The formula is as follows:

[0116]

[0117] S240, take out the authenticated zero watermark image H stored in the watermark database of the third-party registration agency, perform an XOR operation with the binary image feature matrix J′, and obtain the binary watermark image W to be descrambled. A ′.

[0118] W A '=xor(H,J')

[0119] S250, using the Fibonacci Q inverse matrix as the anti-scrambling matrix to perform an anti-scrambling operation on the binary watermark image W A ′ is descrambled to extract the final watermark information W′, and finally the copyright of the three-dimensional medical volume data V′ to be authenticated is determined according to the content information displayed by W′.

[0120] The Fibonacci Q inverse matrix is ​​used as the anti-scrambling matrix to perform the anti-scrambling operation on the binary watermark image W A ′ is descrambled to extract the final watermark information W′. Finally, the copyright of the three-dimensional medical volume data V′ to be authenticated is determined based on the content information displayed by W′, including:

[0121] Using the key Key1, the Fibonacci Q inverse matrix is ​​used as the anti-scrambling matrix to perform the anti-scrambling operation on the binary watermark image W A ′ is descrambled to obtain the descrambled watermark image, and the final watermark image W′ can be extracted. The implementation method is as follows:

[0122]

[0123] In the formula, Q -p is the p-order Fibonacci Q inverse matrix, i.e., the inverse scrambling matrix, x1 and y1 are the binary watermark images W to be inversely scrambled. A ′, x1′ and y1′ are the position coordinates of the pixels in the final watermark image W′ extracted after descrambling.

[0124] The p-order Fibonacci Q inverse matrix is ​​implemented by the following formula:

[0125]

[0126] The copyright ownership of the three-dimensional medical volume data V′ to be authenticated is determined according to the content displayed by the finally extracted watermark image W′.

[0127] The present embodiment provides a three-dimensional medical volume data zero watermark method based on sixteen-dimensional discrete cosine transform. By grouping and processing three-dimensional medical volume data slices and performing Tucker decomposition, the first energy feature map of each group after Tucker decomposition is extracted, which can preserve the correlation between the slices to the greatest extent, which is conducive to improving the overall robustness of the algorithm; then, the first energy feature map of each group after Tucker decomposition is subjected to a two-dimensional non-subsampled discrete wavelet transform to obtain the corresponding low-frequency subband of the same size, which is helpful to improve the embedding capacity of the watermark. The sixteen-dimensional image block is constructed by using the sub-blocks corresponding to the position of each component after non-overlapping division of the low-frequency subband, and the DC coefficient of the image block after sixteen-dimensional discrete cosine transform is directly calculated in the spatial domain, so as to obtain a binary robust feature matrix, which can not only preserve the correlation between the low-frequency components, but also effectively reduce the calculation cost. At the same time, the Fibonacci Q matrix and the dual chaotic system are used to scramble and encrypt the original watermark image and the binary robust feature matrix respectively, which enhances the overall security of the algorithm. It is robust against a variety of image processing attacks, such as noise attacks, filtering attacks, JPEG compression attacks, shearing attacks, rotation attacks, scaling attacks, and translation attacks.

[0128] The process and effect of a three-dimensional medical volume data zero watermarking method based on hexadecimal discrete cosine transform provided by the present invention are specifically described below with reference to specific examples.

[0129] In order to verify the effectiveness of the present invention, the simulation experiment uses MRI and CT brain three-dimensional medical volume data with a size of 128×128×99 for testing. The schematic diagram of a single slice of MRI and CT brain three-dimensional medical volume data is shown in FIG. Figure 3 The original binary watermark image uses a 64×64 "CAUC" image, see Figure 4 .

[0130] The present invention uses normalized cross-correlation (NC) to evaluate the similarity between the final extracted binary watermark image W′ and the original binary watermark image W. NC is defined as follows:

[0131]

[0132] Where N is the size of the original binary watermark image, W(i,j) represents the pixel value of the original binary watermark image at point (i,j), and W′(i,j) represents the pixel value of the final extracted binary watermark image at point (i,j). The closer the NC value is to 1, the more similar the final extracted binary watermark image W′ is to the original binary watermark image W, which means that the method is more robust.

[0133] Figure 5 This is a table of NC values ​​for testing MRI and CT brain three-dimensional medical volume data under several non-geometric attacks in a three-dimensional medical volume data zero watermarking method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention. Figure 5 It can be seen that when facing noise, filtering, and compression attacks, the NC values ​​of the extracted watermarks of MRI and CT brain three-dimensional medical volume data are close to 1, indicating that this method has a very high ability to resist non-geometric attacks and has strong robustness.

[0134] Figure 6 This is a test NC value table of MRI and CT brain three-dimensional medical volume data under geometric attacks such as shearing, scaling, translation and rotation in a three-dimensional medical volume data zero watermark method based on hexadecimal discrete cosine transform provided by an embodiment of the present invention. Figure 6 It can be seen that for MRI and CT brain three-dimensional medical volume data, the three-dimensional medical volume data zero-watermark method provided by this embodiment shows good robustness against common geometric attacks.

[0135] The above specific implementations do not constitute a limitation on the protection scope of the present invention. It should be understood by those skilled in the art that various modifications, combinations, sub-combinations and substitutions can be made according to design requirements and other factors. Any modification, equivalent substitution and improvement made within the spirit and principle of the present invention should be included in the protection scope of the present invention.

Claims

1. A zero watermarking method for three-dimensional medical volume data based on hexadecimal discrete cosine transform, characterized in that: include: Embed zero watermark in original three-dimensional medical volume data and extract zero watermark from the three-dimensional medical volume data to be authenticated; The embedding of a zero watermark in the original three-dimensional medical volume data comprises: The Fibonacci Q matrix is ​​used as the scrambling matrix to scramble the original binary watermark image W of size N×N to obtain the scrambled binary watermark image W A ; The original three-dimensional medical volume data V of size M×M×P is sliced ​​and grouped, and the processed medical image slice group is recorded as U i , i = 1, 2, 3, ..., 15; for each medical image slice group U i Tucker decomposition is performed on all of them, and the first energy feature graph after Tucker decomposition of each group is extracted, marked as E i , i=1,2,3,...,15; for E i Each component in is subjected to k-layer two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i , i=1,2,3,……,15; For the low frequency subband L i Each component in is divided into non-overlapping sub-blocks, and then L i The sub-blocks at the corresponding positions in each component in the hexadecimal image block are used to construct the hexadecimal image block and the DC coefficient D of the image block after the hexadecimal discrete cosine transform is calculated in the spatial domain. st , D st The magnitudes of the real and imaginary parts of st Indicates that, finally using O st The L2 norm of gets the binary robust feature matrix T; The binary image feature matrix J is obtained by scrambling the binary robust feature matrix T using the two-dimensional Sinusoidal-Chebyshev chaotic system and the Spatial Pyramid Matching chaotic system. The scrambled binary watermark image W corresponding to the original binary watermark image W A Perform an XOR operation with the binary image feature matrix J to obtain the authenticated zero watermark image H, save H to the watermark database of the third-party registration agency, and save the control parameters and process parameters used in the binary zero watermark embedding process as the key; The zero-watermark extraction of the three-dimensional medical volume data to be authenticated comprises: The three-dimensional medical volume data V′ with a size of M×M×P and to be authenticated is sliced ​​and grouped, and the processed medical image slice group is recorded as U i ′, i = 1, 2, 3, ..., 15; for each medical image slice group U i ′ are all subjected to Tucker decomposition, and the first energy feature graph after Tucker decomposition of each group is extracted, marked as E i ′, i=1,2,3,...,15; for E i Each component in ′ is transformed by k layers of two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i ′, i = 1, 2, 3, ..., 15; For the low frequency subband L i Each component in ′ is divided into non-overlapping sub-blocks, and then L i The sub-blocks at the corresponding positions in each component of ′ construct the hexadecimal image block and calculate the DC coefficient D of the image block after the hexadecimal discrete cosine transform in the spatial domain st ′, D st The magnitudes of the real and imaginary parts of ′ are expressed by a one-dimensional vector O st ' indicates that finally O st The L2 norm of ′ is used to obtain the binary robust feature matrix T′; The binary image feature matrix J′ is obtained by scrambling the binary robust feature matrix T′ using the two-dimensional Sinusoidal-Chebyshev chaotic system and the Spatial Pyramid Matching chaotic system. Take out the certified zero watermark image H stored in the watermark database of the third-party registration agency, perform an XOR operation with the binary image feature matrix J′, and obtain the binary watermark image W to be descrambled. A ′; The Fibonacci Q inverse matrix is ​​used as the anti-scrambling matrix to treat the binary watermark image W for the anti-scrambling operation. A ′ is descrambled to extract the final watermark information W′, and finally the copyright of the three-dimensional medical volume data V′ to be authenticated is determined according to the content information displayed by W′.

2. The method according to claim 1, characterized in that The Fibonacci Q matrix is ​​used as the scrambling matrix to perform a scrambling operation on the original binary watermark image W of size N×N to obtain a scrambled binary watermark image W A ,include: The Fibonacci Q matrix is ​​used as the scrambling matrix to perform the scrambling operation on the original binary watermark image W, which is expressed by the following formula: In the formula, Q p is the p-order Fibonacci Q matrix, i.e., the scrambling matrix, N is the length and width of the original binary watermark image W, x and y are the position coordinates of the pixels in the original binary watermark image W before scrambling, and x′ and y′ are the scrambled binary watermark image W after scrambling. A The position coordinates of the pixels in the image, mod(·) is the remainder function; The p-order Fibonacci Q matrix is ​​implemented by the following formula: In the formula, F p+1 , F p , F p-1 are the number of terms in the Fibonacci sequence. The order p is regarded as the key Key1 and saved.

3. The method according to claim 1, characterized in that The original three-dimensional medical volume data V of size M×M×P is sliced ​​and grouped, and each medical image slice group U i Tucker decomposition is performed on all groups, and the first energy feature graph E after Tucker decomposition is extracted. i , for E i Each component in is subjected to k-layer two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i ,include: The original 3D medical volume data V of size M×M×P is sliced ​​and divided into 15 medical image slice groups, denoted as U i , i=1,2,3,……,15, and is implemented as follows: Where, P is the total number of slices after slicing of 3D medical volume data, floor(·) is the function of rounding down to an integer, and mod(·) is the function of taking the remainder; Will U i Converted into tensor form and subjected to Tucker decomposition, the medical image slice group U i It can be expressed as: U i =G i ×1V i (1) ×2V i (2) ×3V i (3) In the formula, i=1,2,3,...,15, U i is the tensor representation of each medical image slice group, G i is the core tensor in each group, V (1)i , V (2)i , V (3)i is the factor matrix of each group of tensors in different directions; Let tensor X i =G i ×1V i (1) ×2V i (2) is the tensor U i , extract each group of sub-tensors X i The first sub-image of the medical image slice group U i The first energy characteristic graph after Tucker decomposition is denoted as E i , i=1,2,3,……,15; For E i Perform k-layer two-dimensional non-subsampled discrete wavelet transform to obtain E i The low-frequency subband L corresponding to each component in i , i=1,2,3,……,15.

4. The method according to claim 1, characterized in that The low frequency sub-band L i Divide into non-overlapping sub-blocks, and then use the L i The sub-blocks corresponding to the positions of each component in the image block are used to construct the hexadecimal image block and the DC coefficient D of the image block after the hexadecimal discrete cosine transform is calculated in the spatial domain. st , D st The magnitudes of the real and imaginary parts of st Indicates that, finally using O st The L2 norm of gets the binary robust feature matrix T, including: For L i Divide the non-overlapping sub-blocks into r×r blocks, and then use the L i The sub-blocks corresponding to the positions of each component in the image block are used to construct the hexadecimal image block and the DC coefficient D of the image block after the hexadecimal discrete cosine transform is calculated in the spatial domain. st , D st The magnitudes of the real and imaginary parts of st Indicates that, finally using O st The L2 norm of the coefficient feature matrix S is obtained, which is implemented using the following formula: S(s,t)=||O st ||2 Where, ||·||2 is the calculation symbol of L2 norm, s=1,2,3,……,M / r, t=1,2,3,……,M / r; The mean value of the coefficient feature matrix S is used to binarize the matrix to obtain the binary robust feature matrix T, which is implemented using the following formula:

5. The method according to claim 1, characterized in that The two-dimensional Sinusoidal-Chebyshev chaotic system and the Spatial Pyramid Matching chaotic system are used to perform a scrambling operation on the binary robust feature matrix T to obtain a binary image feature matrix J, including: The binary robust feature matrix T is expanded into a binary extended feature matrix R of size N×N, which is implemented as follows: Where N is the length and width of the original binary watermark image W, M is the length and width of a single slice in the original three-dimensional medical volume data V, r is the length and width of the sub-block, and repmat(·) is the repeat array copy function; A two-dimensional Sinusoidal-Chebyshev system is used to generate a one-dimensional chaotic sequence of length N×N and the mean of the chaotic sequence is used for binarization to obtain a one-dimensional binary sequence Y. The reshape function is then used to reconstruct Y into a two-dimensional matrix and XORed with the binary extended feature matrix R to obtain the encrypted binary feature matrix B. The implementation method is as follows: B = xor(R, reshape(Y, [NN])) In the formula, xor(·) is the function to perform XOR calculation, and reshape(·) is the function to reshape the array. The two-dimensional Sinusoidal-Chebyshev system is expressed by the following formula: In the formula, a, b, c, d are system control parameters, m i , n i is the intermediate variable in the iterative process of the two-dimensional Sinusoidal-Chebyshev chaotic system, m i+1 , n i+1 m i , n i The next state of , mod(·) is the remainder function; The Spatial Pyramid Matching chaotic system is used to generate a one-dimensional chaotic sequence Z of length N×N. The index vector of the one-dimensional chaotic sequence Z in descending order is used to scramble the encrypted binary feature matrix B to obtain the binary image feature matrix J. The following method is used: In the formula, sort(·) is the sorting function, and index is the index vector obtained by sorting the one-dimensional chaotic sequence Z in descending order; The Spatial Pyramid Matching chaotic system is expressed by the following formula: In the formula, η and μ are system control parameters, h is a fixed parameter, q(t) is the intermediate variable in the iterative process of the Spatial Pyramid Matching chaotic system, and q(t+1) is the next state of q(t); When using the two-dimensional Sinusoidal-Chebyshev chaotic system to obtain a one-dimensional binary sequence Y with a length of N×N, the system control parameters a, b, c, d and the system initial values ​​m0, n0 are used as the key Key2. When using the Spatial Pyramid Matching chaotic system to generate a one-dimensional chaotic sequence Z with a length of N×N, the fixed parameter h, the system control parameters η and μ and the system initial value q(0) are used as the key Key3. Key2 and Key3 are saved.

6. The method according to claim 1, characterized in that The three-dimensional medical volume data V′ of size M×M×P to be authenticated is sliced ​​and grouped, and each medical image slice group U i ′ are all subjected to Tucker decomposition, and the first energy feature graph E after Tucker decomposition of each group is extracted i ′, for E i Each component in ′ is subjected to k-layer two-dimensional non-subsampled discrete wavelet transform to obtain the low-frequency subband L corresponding to each component i ',include: The 3D medical volume data V′ with a size of M×M×P and to be authenticated is sliced ​​and divided into 15 medical image slice groups, denoted as U i ′, i = 1, 2, 3, ..., 15, and is implemented as follows: Will U i ′ is converted into tensor form and subjected to Tucker decomposition, then the medical image slice group U i ′ can be expressed as: Let tensor is the tensor U i ′, extract each group of sub-tensors X i ′, namely the medical image slice group U i ′The first energy characteristic graph after Tucker decomposition, denoted as E i ′, i = 1, 2, 3, ..., 15; For E i ′ Perform k-layer two-dimensional non-subsampled discrete wavelet transform to obtain E i The low-frequency subband L corresponding to each component in ′ i ′, i=1,2,3,……,15.

7. The method according to claim 1, characterized in that The low frequency sub-band L i 'Divide into non-overlapping sub-blocks, and then use the L i The sub-blocks corresponding to the positions of each component in ′ construct a hexadecimal image block and calculate the DC coefficient D′ of the image block after the hexadecimal discrete cosine transform in the spatial domain st , D′ st The magnitude of the real part and each imaginary part is expressed by a one-dimensional vector O′ st Indicates that, finally using O′ st The L2 norm of gets the binary robust feature matrix T′, including: For L i ′ divides the non-overlapping sub-blocks into r×r blocks, and then uses L i The sub-blocks corresponding to the positions of each component in ′ construct a hexadecimal image block and calculate the DC coefficient D′ of the image block after the hexadecimal discrete cosine transform in the spatial domain st , D′ st The magnitude of the real part and each imaginary part is expressed by a one-dimensional vector O′ st Indicates that, finally using O′ st The L2 norm of the coefficient characteristic matrix S' is obtained, which is implemented using the following formula: S'(s,t)=||O' st ||2 The mean value of the coefficient feature matrix S′ is used to binarize the matrix to obtain the binary robust feature matrix T′, which is implemented using the following formula:

8. The method according to claim 1, characterized in that The binary watermark image W to be subjected to the descrambling operation A ′ is descrambled to extract the final watermark image W′, including: The binary robust feature matrix T′ is expanded into a binary extended feature matrix R′ of size N×N, which is implemented as follows: Using the key Key2, a two-dimensional Sinusoidal-Chebyshev system is used to generate a one-dimensional chaotic sequence of length N×N and the mean of the chaotic sequence is used for binarization to obtain a one-dimensional binary sequence Y. Then, the reshape function is used to reconstruct Y into a two-dimensional matrix and perform an XOR operation with the binary extended feature matrix R′ to obtain the encrypted binary feature matrix B′. The implementation method is as follows: B'=xor(R',reshape(Y,[NN])) Using the key Key3, the Spatial Pyramid Matching chaotic system is used to generate a one-dimensional chaotic sequence Z, and the index vector of the one-dimensional chaotic sequence Z in descending order is used to scramble the encrypted binary feature matrix B′ to obtain the binary image feature matrix J′. The formula is as follows: Extract the authenticated zero watermark image H stored in the watermark database of the third-party registration agency, perform an XOR operation with the binary image feature matrix J′, and obtain the binary watermark image W to be descrambled. A ′; Using the key Key1, the Fibonacci Q inverse matrix is ​​used as the anti-scrambling matrix to perform the anti-scrambling operation on the binary watermark image W A ′ performs an inverse scrambling operation to obtain an inverse scrambled watermark image, that is, extract the final watermark image W′. The implementation method is as follows: In the formula, Q -p is the p-order Fibonacci Q inverse matrix, i.e., the inverse scrambling matrix, x1 and y1 are the binary watermark images W to be inversely scrambled. A ′, x1′ and y1′ are the position coordinates of the pixels in the final watermark image W′ extracted after descrambling; The p-order Fibonacci Q inverse matrix is ​​implemented by the following formula: The copyright ownership of the three-dimensional medical volume data V′ to be authenticated is determined according to the content displayed by the finally extracted watermark image W′.

Images