CN106997602A - SAR image registration method based on GPU and pyramid mutual information - Google Patents
SAR image registration method based on GPU and pyramid mutual information Download PDFInfo
- Publication number
- CN106997602A CN106997602A CN201710160992.8A CN201710160992A CN106997602A CN 106997602 A CN106997602 A CN 106997602A CN 201710160992 A CN201710160992 A CN 201710160992A CN 106997602 A CN106997602 A CN 106997602A
- Authority
- CN
- China
- Prior art keywords
- pyramid
- width
- image
- top layer
- gpu
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
- G06T2207/10044—Radar image
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20016—Hierarchical, coarse-to-fine, multiscale or multiresolution image processing; Pyramid transform
Abstract
The invention discloses a kind of SAR image registration method based on GPU and pyramid mutual information, the problem of mainly solving the slow registering speed of prior art image and low precision.Its implementation is：1.GPU generates two width pyramid top layer images through 1/4 down-sampling, while CPU calculates top layer images moving range；2. the shift position of top layer images is calculated in top layer images moving range；3.CPU calculates middle tomographic image moving range by top layer images shift position, while GPU generates tomographic image in the middle of two width through 1/2 down-sampling；4. tomographic image shift position in the middle of being calculated in middle tomographic image moving range；5. image moving range subject to registration is calculated according to middle tomographic image shift position；6. calculating final registration position in image moving range subject to registration, images after registration is generated.GPU is applied in pyramid mutual information registration algorithm by the present invention, the speed and precision of registration is improved, available for Real-time Remote Sensing data analysis.
Description
Technical field
The invention belongs to Radar Technology field, more particularly to a kind of method for registering of SAR image, available for Real-time Remote Sensing number
According to analysis.
Background technology
Pyramid mutual information image registration is one of current the most frequently used registration Algorithm, and mutual information is calculated and adopted with image by it
Sample combines, small with operand, the good advantage of registration effect.Pyramid mutual information image registration algorithm now exists more
Realized on CPU, the proposition such as Li Qiaoliang, Wang Guoyou, Liu Jianguo of the Central China University of Science and Technology based on batten pyramid and mutual information
Fast image registration algorithm, using batten pyramid reduction amount of calculation is built, algorithm single thread on CPU is run.Zhang Zanxia,
The layering remote sensing image processing method based on mutual information that Peng Jiaxiong, Wang Hongqun are proposed, using the golden word of wavelet transformation structural map picture
Tower carries out mutual information registration, is run also on CPU.
With the raising of SAR image resolution ratio, the image registration algorithm run time run on CPU is elongated, so that nothing
Method meets requirement of real-time, is occurred using the GPU registration Algorithms accelerated parallel；National University of Defense technology Zhao enters based on the distant of GPU
The GPU mutual information image registration algorithms of realization are given in sense image Parallel Processing algorithm and its paper of optimisation technique research,
But this method, for the GPU using single-precision floating point computing, registration accuracy is difficult to control.
The content of the invention
It is a kind of based on GPU and pyramid mutual information it is an object of the invention to for above-mentioned the deficiencies in the prior art, propose
SAR image registration method, to reduce the image registration time, and improve the precision of registration.
The technical scheme is that：
By the way of CPU+GPU is combined, carry out image sampling using GPU and mutual information is calculated, golden word is carried out using CPU
Tower has the calculating of regional extent and moving range per tomographic image.1/4 and 1/2 down-sampling is carried out to picture and builds 3 layers of pyramid,
By successively calculating registration position, final registering image is generated, implementation step includes as follows：
(1) two width SAR gray level images are obtained, and are sent it in GPU video memorys as two original SAR images, GPU pairs
This two width SAR image carries out 1/4 down-sampling respectively, obtains the first width pyramid top layer images D11With the second width pyramid top level diagram
As D12；The second width pyramid top layer images D is initialized by CPU simultaneously12Moving range；
(2) the second width pyramid top layer images D of pyramid top layer maximum mutual information correspondence is calculated12Shift position：
(2a) GPU chooses a mobile vector and mobile second width pyramid top layer images D in moving range12, generation
First width pyramid top layer has area image D11' and the shared area image D of the second width pyramid top layer12', then calculate this two
The mutual information MI of width image1；
(2b) traversal in moving range chooses new mobile vector, repeats (2a), obtains the set of association relationship；Find out
Maximum mutual information value in set, and preserve the corresponding second width pyramid top layer images D of the maximum mutual information value12Mobile position
Put；
(3) GPU carries out 1/2 down-sampling to two original SAR images, obtains the two images D in pyramid intermediate layer21With
D22；The second width pyramid top layer images D that CPU is preserved according to step (2b)12Shift position, is determined in the middle of the second width pyramid
Tomographic image D22Moving range；
(4) tomographic image D in the middle of maximum mutual information the second width pyramid of correspondence of pyramid intermediate layer is calculated22Shift position：
(4a) GPU is chosen in the moving range that step (3) is determined in a mobile vector and mobile second width pyramid
Between tomographic image D22, the first width pyramid intermediate layer of generation has area image D21' and the shared region in the second width pyramid intermediate layer
Image D22', then calculate the mutual information MI of this two images2；
(4b) traversal in the moving range that step (3) is determined chooses new mobile vector, repeats (4a), obtains mutual information
The set of value；Maximum mutual information value in set is found out, and preserves the maximum mutual information value the second width pyramid intermediate layer figure of correspondence
As D22Shift position；
(5) the second width pyramid centre tomographic image D that CPU is preserved according to step (4b)22Shift position, determines that the second width is former
The moving range of beginning SAR image, and a block space F is opened up for storing pyramid bottom association relationship and mobile arrow in internal memory
Amount, initialization space F association relationship for single-precision floating point bear it is infinite；
(6) shift position of pyramid bottom maximum mutual information second original SAR image of correspondence is calculated：
(6a) GPU chooses a mobile vector and mobile second original SAR figure in the moving range that step (5) is determined
Picture, the first width pyramid bottom of generation has area image D31' and the shared area image D of the second width pyramid bottom32', then count
Calculate the mutual information MI of this two images3；
(6b) CPU is by mutual information MI3It is compared with the association relationship that is stored in the F of space：If mutual information MI3Than space F
The association relationship of storage is big, then the association relationship for changing space F storages is MI3, while the mobile vector of space F storages is changed,
Otherwise, without any operation；
(6c) traversal in the moving range that step (5) is determined chooses new mobile vector, repeats (6a) and (6b)；Space
The association relationship of F storages is pyramid bottom maximum mutual information value MI, takes out the mobile vector that space F is preserved, GPU is according to this
Mobile vector moves second original SAR image again, obtains images after registration D.
The invention has the advantages that：
The present invention is realized due to applying to GPU in the mutual information calculating of the shared area image of every layer of pyramid with CPU
Traditional pyramid mutual information registration algorithm is compared, and reduces run time；Simultaneously because by generating every layer of consensus of pyramid
Area image, improves registration accuracy.
Test result shows that the present invention is compared with the pyramid mutual information registration algorithm that CPU is realized, the registering time significantly subtracts
It is small；The present invention is compared with traditional pyramid mutual information registration algorithm that GPU is realized, registration accuracy increases.
Brief description of the drawings
Fig. 1 is implementation process block diagram of the invention；
Fig. 2 is the original SAR image that emulation experiment of the present invention is used；
Fig. 3 is the images after registration that emulation experiment of the present invention is drawn；
Fig. 4 is the overlay chart picture of the images after registration and first original SAR image that are drawn with the present invention.
Embodiment
Reference picture 1, of the invention comprises the following steps that：
Step 1, GPU generates two width top layer images D11And D12, the second width pyramid top layer images D of CPU calculating12Movement
Scope.
(1a) GPU generates the first width pyramid top layer images D11With the second width pyramid top layer images D12：
(1a1) is according to the wide W of first original SAR image1With high H1, GPU opens up a width of in video memoryIt is a height ofThe first pyramid top layer two-dimensional array α1, whereinExpression is rounded downwards, ceil4Table
Show the multiple rounded up as 4；
(1a2) GPU is first the first width pyramid top layer images D11Distribute 32 × 32 thread blocks, each thread block reallocationIndividual thread, makes the first pyramid top layer two-dimensional array α1In all index (x1,y1) meetAnd meetElement respectively correspond to a GPU thread, whereinExpression takes upwards
It is whole；
(1a3) utilizes the first pyramid top layer two-dimensional array α1Middle element (x1,y1) corresponding thread, by the element assignment
For (4x in first original SAR image1,4y1) pixel value at coordinate, now the first pyramid top layer two-dimensional array α1Deposit
Store up the first width pyramid top layer images D11；
(1a4) is according to the wide W of second original SAR image2With high H2, opened up in GPU video memorys a width ofIt is a height ofThe second pyramid top layer two-dimensional array β1；
(1a5) GPU is first the second width pyramid top layer images D12Distribute 32 × 32 thread blocks, each thread block reallocationIndividual thread；Make the second pyramid top layer two-dimensional array β1In all index (x1',y1')
MeetWithElement respectively correspond to a GPU thread；
(1a6) utilizes the second pyramid top layer two-dimensional array β1Middle element (x1',y1') corresponding thread, the element is assigned
It is worth for second original SAR image (4x1',4y1') pixel value of coordinate, now the second pyramid top layer two-dimensional array β1Deposit
Store up the first width pyramid top layer images D12；
(1b) CPU calculates the second width pyramid top layer images D12Moving range：
Second original SAR image that CPU is inputted according to user different zones initial moving range, i.e.,：Level is moved
Dynamic interval [- w, w], vertically moves interval [- h, h], and anglec of rotation interval [- ang, ang] calculates the second width pyramid top level diagram
As D12It is in the moving range of different zones：Move horizontally intervalVertically move intervalAnglec of rotation area
Between [- ang, ang].
Step 2, GPU calculates the second width pyramid top layer images D12Shift position.
(2a) GPU generates the first width pyramid top layer and has area image D11' and the shared region of the second width pyramid top layer
Image D12'：
(2a1) translates the second width pyramid top layer images D12, the image Ds after being translated12；
The second width pyramid top layer images D that GPU is determined from step (1b)12A mobile vector is chosen in moving range
(x1,y1,ang1), wherein x1Represent the number of pixels moved horizontally, y1Represent the number of pixels vertically moved, ang1Represent to scheme
As D12Center is origin dextrorotation gyration, and according to the second width pyramid top layer images D12Wide W12With high H12, in video memory
In open up a width of W12+2x1, a height of H12+2y1Pyramid top layer translation array δ1, then by δ1Reset；
GPU is translation array δ1Distribute 32 × 32 thread blocks, each thread block reallocation
Individual thread, and make the second width pyramid top layer images D12Each pixel respectively correspond to a thread；
Thread is by the second width pyramid top layer images D12Middle index is (x12,y12) pixel copy it is flat to pyramid top layer
Move array δ1In (x12+x1,y12+y1) position, then translate array δ1Store the second width pyramid top layer images D12After translation
Image Ds12；
Image Ds after (2a2) rotation translation12, the image Dm after being moved1：
GPU opens up a width of (2x in video memory1+W12)cos(ang1)+W12sin(ang1), a height of (2x1+W12)sin(ang1)
+W12cos(ang1) pyramid top layer rotation array δ1', then pyramid top layer translation array δ is obtained respectively1With pyramid top layer
Rotate array δ1' first address, by the two first address and ang1It is input in the instruction of GPU affine transformation, the instruction is performed
Afterwards, image Ds12Postrotational image is stored in pyramid top layer rotation array δ1' in, that is, obtain the second width pyramid top layer
Image D12Image Dm after movement1；
(2a3) CPU calculates the image Dm after movement1With the first width pyramid top layer images D11The square of lap has
Regional extent：
(2a31) CPU is according to image D12Mobile vector (x1,y1,ang1), calculate the Breadth Maximum in the square consensus domain
L：
L=min (Wm1,W11,Hm1,H11),
Wherein, Wm1And Hm1For the image Dm after movement1Wide and height, W11And H11For the first width pyramid top layer images D11
Wide and height；
(2a32) is by the square consensus domain in image D11In be expressed as the first top layer have region S11, CPU zonings
S11Wide cx11With high cy11For：
CPU zonings S11Top left corner apex is in the first width pyramid top layer images D11In coordinate be：
((W11-cx11)/2+x1,(H11-cy11)/2-y1)；
(2a33) is by the square consensus domain in image D12In be expressed as the second top layer have region S12, CPU calculates the area
The wide cx in domain12With high cy12For：
cx12=cy12=cx11,
CPU zonings S12Top left corner apex is in the second width pyramid top layer images D12In coordinate be：
((Wm1-cx12)/2+x1,(Hm1-cy12)/2-y1)；
(2a4) GPU generates the first width pyramid top layer and has area image D11'：
GPU allocated sizes in video memory are cx11×cy11The first top layer have region two-dimensional array, the first top layer is total to
There is region S11In all pixels copy in the two-dimensional array, obtain the first width pyramid top layer and have area image D11'；
(2a5) GPU generates the second width pyramid top layer and has area image D12'：
GPU allocated sizes in video memory are cx12×cy12The second top layer have region two-dimensional array, the second top layer is total to
There is region S12In all pixels copy in the two-dimensional array, obtain the second width pyramid top layer and have area image D12'；
(2b) GPU calculates two width pyramid top layers and has area image D11' and D12' mutual information：
(2b1) GPU calculates the first width pyramid top layer and has area image D11' entropy H11：
Wherein P11Area image D is had for the first width pyramid top layer11' grey level histogram；
(2b2) GPU calculates the second width pyramid top layer and has area image D12' entropy H12：
Wherein P12Area image D is had for the second width pyramid top layer12' grey level histogram；
(2b3) calculates two width pyramid top layers and has area image D11' and D12' joint entropy H1：
Wherein P1Area image D is had for two width pyramid top layers11' and D12' joint grey level histogram；
(2b4) calculates two width top layers and has area image D11' and D12' mutual information MI1：
(2c) GPU preserves the second width pyramid top layer images D12Shift position：
GPU in the initial moving range of different zones, is opened according to second original SAR image in step (1b) in video memory
Ward off size be 32 × w × h × ang one-dimensional integer array as association relationship set ε；
The second width pyramid top layer images D that GPU is determined from step (1b)12Moving range in traversal choose new movement
Vector, and (2a) and (2b) is repeated, by obtained association relationship and the second width pyramid top layer images D chosen12Mobile vector
Add in association relationship set ε；
GPU opens up the top layer registration position array E that size is 3 in video memory1, then find out mutual information maximum in set ε
Value, by the corresponding second width top layer images D of the maximum mutual information value12Mobile vector deposit registration position array E1In, with level
Put array E1Store the second width pyramid top layer images D12Shift position.
Step 3, GPU generates tomographic image D in the middle of two width21And D22, CPU the second width pyramids of calculating centre tomographic image D22's
Moving range.
(3a) GPU generates tomographic image D in the middle of the first width pyramid21With tomographic image D in the middle of the second width pyramid22：
(3a1) is according to the wide W of first original SAR image1With high H1, GPU opens up a width of in video memory
It is a height ofThe first pyramid intermediate layer two-dimensional array α2；
(3a2) GPU is tomographic image D in the middle of the first width pyramid21Distribute 32 × 32 thread blocks；Each thread block reallocationIndividual thread, makes the first pyramid intermediate layer two-dimensional array α2In all index (x2,y2) full
FootAnd meetElement respectively correspond to a GPU thread；
(3a3) utilizes the first pyramid intermediate layer two-dimensional array α2Middle element (x2,y2) corresponding thread, the element is assigned
It is worth for (2x in first original SAR image2,2y2) pixel value at coordinate, now the first pyramid intermediate layer two-dimensional array α2
Tomographic image D in the middle of as the first width pyramid21；
(3a4) is according to the wide W of second original SAR image2With high H2, opened up in GPU video memorys a width of
It is a height ofThe second pyramid intermediate layer two-dimensional array β2；
(3a5) GPU is tomographic image D in the middle of the second width pyramid22Distribute 32 × 32 thread blocks, each thread block reallocationIndividual thread；Make the second pyramid intermediate layer two-dimensional array β2In all index (x2',
y2') meetWithElement respectively correspond to a GPU thread；
(3a6) utilizes the second pyramid intermediate layer two-dimensional array β2Middle element (x2',y2') corresponding thread, by the element
It is entered as second original SAR image (2x2',2y2') pixel value at coordinate, now the second pyramid intermediate layer two-dimensional array
β2Tomographic image D in the middle of as the second width pyramid22。
(3b) CPU calculates tomographic image D in the middle of the second width pyramid22Moving range：
The image D that CPU is preserved according to step (2c)12Shift position, first obtains the horizontal-shift a of the shift position2, it is vertical
Offset b2With anglec of rotation θ2；
Tomographic image D in the middle of pyramid is determined again22In the moving range of different zones, that is, it is [2a to move horizontally interval2-4,
2a2+ 4], it is [2b to vertically move interval2-4,2b2+ 4], anglec of rotation interval is [θ2-2,θ2+2]。
Step 4, tomographic image D in the middle of the second width pyramid is calculated22Shift position.
(4a) GPU generates two width pyramid intermediate layers and has area image D21' and D22'：
The second width pyramid centre tomographic image D that GPU is determined from step (3b)22A mobile vector is chosen in moving range
(x2,y2,ang2), wherein x2Represent the number of pixels moved horizontally, y2Represent the number of pixels vertically moved, ang2Represent second
Tomographic image D in the middle of width pyramid22Dextrorotation gyration, and tomographic image D in the middle of the second width pyramid is moved by the vector22, obtain
Image Dm after to the movement of pyramid intermediate layer2；
CPU calculates the image Dm after the movement2With tomographic image D in the middle of the first width pyramid21The square of lap has
Regional extent, and this is had into region tomographic image D in the middle of the first width pyramid21In be expressed as the first intermediate layer have region
S21, image Dm after being moved in pyramid intermediate layer2In be expressed as the second intermediate layer have region S22；
GPU is by region S21And S22Copy new video memory space to, form the first width pyramid intermediate layer and have area image
D21' and the shared area image D in the second width pyramid intermediate layer22'；
(4b) calculates two width pyramid intermediate layers and has area image D21' and D22' mutual information：
(4b1) GPU calculates the first width pyramid intermediate layer and has area image D21' entropy H21：
Wherein P21Area image D is had for the first width pyramid intermediate layer21' grey level histogram；
(4b2) GPU calculates the second width pyramid intermediate layer and has area image D22' entropy H22：
Wherein P22Area image D is had for the second width pyramid intermediate layer22' grey level histogram；
(4b3) GPU calculates tomographic image D in the middle of two width pyramids21' and D22' joint entropy H2：
Wherein P2Area image D is had for two width pyramid intermediate layers21' and D22' joint grey level histogram；
(4b4) GPU calculates two width pyramid intermediate layers and has area image D21' and D22' mutual information MI2：
(4c) GPU preserves tomographic image D in the middle of the second width pyramid22Shift position：
The association relationship set ε that GPU produces step (2c) is reset；In the second width pyramid determined again from step (3b)
Between tomographic image D22Traversal chooses tomographic image D in the middle of the second new width pyramid in moving range22Mobile vector, and this is moved
Vector is added in association relationship set ε, repeats (4a) and (4b), obtained association relationship is also added into association relationship set ε
In；
GPU opens up the intermediate layer registration position array E that size is 3 in video memory2, then find out in association relationship set ε
Maximum association relationship, by tomographic image D in the middle of the corresponding second width pyramid of the maximum mutual information value22In the middle of mobile vector deposit
Layer registration position array E2In, that is, obtain tomographic image D in the middle of pyramid22Shift position.
Step 5, CPU calculates the moving range of second original SAR image, and opens up registration position space F.
(5a) CPU calculates the moving range of second original SAR image：
The second width pyramid centre tomographic image D that CPU is preserved according to step (4c)22Shift position, first obtains the movement position
The horizontal-shift a put3, vertical shift b3With anglec of rotation θ3；Movement of second original SAR image in different zones is determined again
Scope, that is, it is [2a to move horizontally interval3-4,2a3+ 4], it is [2b to vertically move interval3-4,2b3+ 4], anglec of rotation interval is
[θ3-2,θ3+2]。
(5b) CPU opens up registration position space F：
CPU is used to store pyramid bottom association relationship and second original SAR image movement arrow in internal memory opening space F
Amount, and initialize space F association relationship born for single-precision floating point it is infinite.
Step 6, the registration position of second original SAR image is calculated.
(6a) GPU generates two width bottoms and has area image D31' and D32'：
GPU chooses a mobile vector (x out of step (5a) is determined second original SAR image moving range3,y3,
ang3), wherein x3Represent the number of pixels moved horizontally, y3Represent the number of pixels vertically moved, ang3Expression second is original
The dextrorotation gyration of SAR image, and second original SAR image is moved by the vector, obtain second original SAR image
Image Dm after movement3, obtain the image Dm after the movement3With the square consensus of first original SAR image lap
Domain；
The shared region is expressed as the first pyramid bottom in first original SAR image and has region S31, moving
Image Dm after dynamic3In be expressed as the second pyramid bottom have region S32；
GPU is by region S31And S32Copy new video memory space to, form the first width pyramid bottom and have area image D31'
Area image D is had with the second width pyramid bottom32'；
(6b) calculates two width bottoms and has area image D31' and D32' mutual information：
(6b1) GPU calculates the first width pyramid bottom and has area image D31' entropy H31：
Wherein P31Area image D is had for the first width pyramid bottom31' grey level histogram；
(6b2) calculates the second width bottom and has area image D32' entropy H32：
Wherein P32Area image D is had for the second width bottom32' grey level histogram；
(6b3) calculates two width bottoms and has area image D31' and D32' joint entropy H3：
Wherein P3Area image D is had for two width pyramid bottoms31' and D32' joint grey level histogram；
(6b4) calculates two width bottoms and has area image D31' and D32' mutual information MI3：
(6c) updates association relationship and mobile vector：
CPU is by mutual information MI3It is compared with the association relationship that is stored in the F of registration position space：If mutual information MI3Than
The association relationship of space F storages is big, then the association relationship for changing registration position space F storages is MI3, while changing registration position
The mobile vector of space F storages is (x3,y3,ang3), otherwise, without any operation；
(6d) GPU obtains the registration position of second original SAR image：
GPU traversals out of step (5) is determined moving range choose new mobile vector, repeat (6a), (6b) and (6c),
The association relationship of registration position space F storages is pyramid bottom maximum mutual information value MI, the mobile arrow stored according to space F
Measure the registration position to second original SAR image；
(6e) generates images after registration：
GPU obtains the mobile vector that space F is preserved, and moves second original SAR image again according to the mobile vector, obtains
To images after registration D.
The effect of the present invention can be further illustrated by following emulation experiment：
1st, experimental situation：
GPU：NVIDIA Tesla K40c,
Processor：Intel (R) Xeon (R) E5-2630v3 2.40GHz (2 processor),
Internal memory：64.0GB,
Hard disk：2T,
Operating system：Microsoft Windows 7,SP1 64,
Programming tool：Microsoft Visual Studio 2010, CUDA6.5.
2nd, experiment content and interpretation of result：
Experiment 1, with the present invention to second original SAR shown in the first original SAR image and Fig. 2 b shown in Fig. 2 a
Image carries out registration, as a result as shown in figure 3, the image again by Fig. 3 after registering is overlapping with Fig. 2 a progress, obtains overlay chart picture, such as
Shown in Fig. 4.
As seen from Figure 4, overlapping image clearly and profile show that the registration effect of the present invention is preferable without ghost image.
Experiment 2, the traditional pyramid mutual information registration method realized with the inventive method and existing CPU is to various sizes of
SAR image carry out registration, and measure obtain both approaches registration take, its result is as shown in table 1：
The inventive method of table 1 realizes the registering time-consuming contrast of traditional pyramid mutual information registration method with existing CPU
As shown in Table 1, the inventive method is compared with traditional pyramid mutual information registration method that existing CPU is realized, registration
It is time-consuming significantly to shorten；
Experiment 3, the traditional pyramid mutual information registration method realized with the inventive method and existing GPU known to several to matching somebody with somebody
Registration is to test registration accuracy again for the SAR image that level is put, and the registration position that both approaches are obtained is as shown in table 2.
The registration position of the inventive method of table 2 and the existing GPU traditional pyramid mutual information registration methods realized is contrasted
Picture size, horizontal-shift, vertical shift are in units of pixel in table 2, and second original SAR image is artificial
What mobile first original SAR image was obtained, second original SAR image is experiment picture.
As shown in Table 2, the present invention can correct registering all experiment pictures, its obtained registration position with it is correct registering
Position is identical；Though existing method correctly can test picture 1 by registration, fail correct registration and test picture 2, it obtains testing picture
The horizontal-shift of 2 registration positions differs 1 pixel, and the vertical shift of the experiment registration position of picture 2 with correct horizontal-shift
Also 1 pixel is differed with correct vertical shift；Existing method also fails to correct registration experiment picture 3, and it obtains testing picture 3
The horizontal-shift of registration position also differs 1 pixel with correct horizontal-shift.The tradition golden word of the invention realized with existing GPU
Tower mutual information registration method is compared, and registration accuracy is improved.
Above is to example of the present invention, not constituting any limitation of the invention, it is clear that the skill of this area
Art personnel can carry out various changes and modification in the case where not departing from this hair spirit and principle, but these still fall within the protection of the present invention
Within the scope of.
Claims (8)
1. the SAR image registration method based on GPU and pyramid mutual information, including：
(1) two width SAR gray level images are obtained, and are sent it in GPU video memorys as two original SAR images, GPU to this two
Width SAR image carries out 1/4 down-sampling respectively, obtains the first width pyramid top layer images D11With the second width pyramid top layer images
D12；The second width pyramid top layer images D is initialized by CPU simultaneously12Moving range；
(2) the second width pyramid top layer images D of pyramid top layer maximum mutual information correspondence is calculated12Shift position：
(2a) GPU chooses a mobile vector and mobile second width pyramid top layer images D in moving range12, generation first
Width pyramid top layer has area image D11' and the shared area image D of the second width pyramid top layer12', then calculate this two width figure
The mutual information MI of picture1；
(2b) traversal in moving range chooses new mobile vector, repeats (2a), obtains the set of association relationship；Find out set
Middle maximum mutual information value, and preserve the corresponding second width pyramid top layer images D of the maximum mutual information value12Shift position；
(3) GPU carries out 1/2 down-sampling to two original SAR images, obtains the two images D in pyramid intermediate layer21And D22；CPU
The the second width pyramid top layer images D preserved according to step (2b)12Shift position, determines tomographic image in the middle of the second width pyramid
D22Moving range；
(4) tomographic image D in the middle of maximum mutual information the second width pyramid of correspondence of pyramid intermediate layer is calculated22Shift position：
(4a) GPU chooses a mobile vector and mobile second width pyramid intermediate layer in the moving range that step (3) is determined
Image D22, the first width pyramid intermediate layer of generation has area image D21' and the shared area image in the second width pyramid intermediate layer
D22', then calculate the mutual information MI of this two images2；
(4b) traversal in the moving range that step (3) is determined chooses new mobile vector, repeats (4a), obtains association relationship
Set；Maximum mutual information value in set is found out, and preserves tomographic image D in the middle of the maximum mutual information value the second width pyramid of correspondence22
Shift position；
(5) the second width pyramid centre tomographic image D that CPU is preserved according to step (4b)22Shift position, determine second it is original
The moving range of SAR image, and a block space F is opened up for storing pyramid bottom association relationship and mobile vector in internal memory,
Initialization space F association relationship is born infinite for single-precision floating point.
(6) shift position of pyramid bottom maximum mutual information second original SAR image of correspondence is calculated：
(6a) GPU chooses a mobile vector and mobile second original SAR image in the moving range that step (5) is determined,
Generate the first width pyramid bottom and have area image D31' and the shared area image D of the second width pyramid bottom32', then calculate
The mutual information MI of this two images3；
(6b) CPU is by mutual information MI3It is compared with the association relationship that is stored in the F of space：If mutual information MI3Than space F storages
Association relationship it is big, then the association relationship for changing space F storages is MI3, while the mobile vector of space F storages is changed, otherwise,
Without any operation；
(6c) traversal in the moving range that step (5) is determined chooses new mobile vector, repeats (6a) and (6b)；Space F is deposited
The association relationship of storage is pyramid bottom maximum mutual information value MI, takes out the mobile vector that space F is preserved, GPU is according to the shifting
Dynamic vector moves second original SAR image again, obtains images after registration D.
2. GPU is adopted under carrying out 1/4 respectively to two original SAR images in the method as described in claim 1, wherein step (1)
Sample, is to have used GPU Unified Device computing architectures CUDA to carry out, step is as follows：
(1a) GPU obtains the first width pyramid top layer images D11：
(1a1) is according to the wide W of first original SAR image1With high H1, opened up in GPU video memorys a width ofIt is high
ForThe first pyramid top layer two-dimensional array α1, whereinExpression is rounded downwards, ceil4Expression rounds up as 4
Multiple；
(1a2) GPU distributes 32 × 32 thread blocks；Each thread block reallocationIndividual line
Journey, makes the first pyramid top layer two-dimensional array α1In all index (x1,y1) meetAnd meetElement respectively correspond to a GPU thread, whereinExpression rounds up；
(1a3) utilizes the first pyramid top layer two-dimensional array α1Middle element (x1,y1) corresponding thread, the element is entered as
(4x in one original SAR image1,4y1) pixel value at coordinate, now the first pyramid top layer two-dimensional array α1Store
First width pyramid top layer images D11；
(1b) GPU obtains the second width pyramid top layer images D12：
(1b1) is according to the wide W of second original SAR image2With high H2, opened up in GPU video memorys a width ofIt is high
ForThe second pyramid top layer two-dimensional array β1；
(1b2) GPU distributes 32 × 32 thread blocks, each thread block reallocationIndividual line
Journey；Make the second pyramid top layer two-dimensional array β1In all index (x1',y1') meetWithElement respectively correspond to a GPU thread；
(1b3) utilizes the second pyramid top layer two-dimensional array β1Middle element (x1',y1') corresponding thread, the element is entered as
Second original SAR image (4x1',4y1') pixel value of coordinate, now the second pyramid top layer two-dimensional array β1Store
First width pyramid top layer images D12。
3. GPU chooses a mobile vector and mobile the in moving range in the method as described in claim 1, step (2a)
Two width pyramid top layer images D12, it is the second width pyramid top layer images D determined from step (1)12One is chosen in moving range
Individual mobile vector (x1,y1,ang1), wherein x1Represent the number of pixels moved horizontally, y1The number of pixels vertically moved is represented,
ang1Represent with image D12Center is origin dextrorotation gyration, and the second width pyramid top layer images D is moved by the vector12,
Obtain the image and the first width pyramid top layer images D after movement11The square shared region of lap, and this is had into region
In D11In be expressed as S11, S is expressed as in image after movement12；GPU is by region S11And S12Copy new video memory space, shape to
Area image D is had into the first width pyramid top layer11' and the shared area image D of the second width pyramid top layer12'。
4. method as claimed in claim 2, wherein with mobile vector (x1,y1,ang1) mobile second width pyramid top layer images
D12, carry out as follows：
(2a1) translates the second width pyramid top layer images D12：
GPU is according to the second width pyramid top layer images D12Wide W12With high H12A width of W is opened up in video memory12+2x1, a height of H12+
2y1Pyramid top layer translation array δ1, and pyramid top layer is translated into array δ1Reset；
GPU 32 × 32 thread blocks of reallocation, each thread block distributionIndividual thread, and make the second width
Pyramid top layer images D12Each pixel respectively correspond to a thread；
Thread is by image D12Middle element (x12,y12) copy pyramid top layer translation array δ to1In (x12+x1,y12+y1) position,
Obtain the second width pyramid top layer images D12Image Ds after translation12；
Image Ds after (2a2) rotation translation12：
GPU opens up a width of (2x in video memory1+W12)cos(ang1)+W12sin(ang1), a height of (2x1+W12)sin(ang1)+
W12cos(ang1) pyramid top layer rotation array δ1', then pyramid top layer translation array δ is obtained respectively1With pyramid top layer
Rotate array δ1' first address, by the two first address and ang1It is input in the instruction of GPU affine transformation, the instruction is performed
Afterwards, image Ds12Postrotational image is stored in pyramid top layer rotation array δ1' in, that is, obtain the second width pyramid top layer
Image D12Image after movement.
5. the pyramid top layer images D that CPU is obtained according to step (2b) in method as claimed in claim 1, wherein step (3)12
Shift position calculates tomographic image D in the middle of pyramid22Moving range, be that image D is first obtained by CPU12The level of shift position is inclined
Move a1, vertical shift b1With anglec of rotation θ1；Tomographic image D in the middle of pyramid is determined again22In the moving range of different zones, i.e.,：
It is [2a to move horizontally interval1-4,2a1+ 4],
It is [2b to vertically move interval1-4,2b1+ 4],
Anglec of rotation interval is [θ1-2,θ1+2]。
6. GPU carries out 1/2 down-sampling to two original SAR images in method as claimed in claim 1, wherein step (3), obtain
Tomographic image D in the middle of first width pyramid21With tomographic image D in the middle of the second width pyramid22, carry out as follows：
(3a) GPU is according to the wide W of first original SAR image1With high H1, opened up in video memory a width ofIt is a height ofThe first pyramid intermediate layer two-dimensional array α2；
(3b) GPU first distributes 32 × 32 thread blocks, then is distributed for each thread blockNumber
The thread of amount, makes the first pyramid intermediate layer two-dimensional array α2In all index (x12,y12) meetWithElement respectively correspond to a GPU thread；
(3c) thread is by the first pyramid intermediate layer two-dimensional array α2Element (x12,y12) it is entered as first original SAR image
In (2x12,2y12) pixel value at coordinate, then the first pyramid intermediate layer two-dimensional array α2As the first width pyramid intermediate layer
Image D21；
(3d) GPU is according to the wide W of second original SAR image2With high H2, opened up in video memory a width ofIt is a height ofThe second pyramid intermediate layer two-dimensional array β2；
(3e) GPU distributes 32 × 32 thread blocks, each thread block reallocationIndividual thread；
Make the second pyramid intermediate layer two-dimensional array β2In all index (x22',y22') meetWithElement respectively correspond to a GPU thread；
(3f) thread is by the second pyramid intermediate layer two-dimensional array β2Element (x22',y22') it is entered as second original SAR figure
(the 2x as in22',2y22') pixel value at coordinate, then the second pyramid intermediate layer two-dimensional array β2Store the golden word of the second width
Tomographic image D in the middle of tower22。
7. GPU chooses a mobile vector and moved in moving range in the method as described in claim 1, wherein step (4a)
Tomographic image D in the middle of dynamic second width pyramid22, it is tomographic image D in the middle of the second width pyramid determined from step (3)22Moving range
One mobile vector (x of interior selection2,y2,ang2), wherein x2Represent the number of pixels moved horizontally, y2Represent the picture vertically moved
Plain number, ang2The dextrorotation gyration of representative image；Tomographic image D in the middle of second width pyramid is moved by the vector22, obtain
Image and the first width pyramid centre tomographic image D after this movement21The square shared region of lap, and this is had into region
In D21In be expressed as S21, S is expressed as in image after movement22；GPU is by the two regions S21And S22Copy new video memory to empty
Between, form the first width pyramid intermediate layer and have area image D21' and the shared area image D in the second width pyramid intermediate layer22'。
8. the pyramid intermediate layer figure that CPU is obtained according to step (4b) in the method as described in claim 1, wherein step (5)
As D22Shift position calculate pyramid second original SAR image of bottom moving range, be that image D is first obtained by CPU22
The horizontal-shift a of shift position2, vertical shift b2With anglec of rotation θ2；Determine second original SAR image in different zones again
Moving range, i.e.,：
It is [2a to move horizontally interval2-4,2a2+ 4],
It is [2b to vertically move interval2-4,2b2+ 4],
Anglec of rotation interval is [θ2-2,θ2+2]。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710160992.8A CN106997602B (en) | 2017-03-17 | 2017-03-17 | SAR image registration method based on GPU and pyramid mutual information |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710160992.8A CN106997602B (en) | 2017-03-17 | 2017-03-17 | SAR image registration method based on GPU and pyramid mutual information |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN106997602A true CN106997602A (en) | 2017-08-01 |
| CN106997602B CN106997602B (en) | 2019-08-06 |
Family
ID=59431801
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710160992.8A Active CN106997602B (en) | 2017-03-17 | 2017-03-17 | SAR image registration method based on GPU and pyramid mutual information |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN106997602B (en) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108742678A (en) * | 2018-06-01 | 2018-11-06 | 妙智科技（深圳）有限公司 | Method for registering images, device and computer readable storage medium |
| CN111445503A (en) * | 2020-03-25 | 2020-07-24 | 桂林电子科技大学 | Pyramid mutual information image registration method based on parallel programming model on GPU cluster |
| CN112862868A (en) * | 2021-01-31 | 2021-05-28 | 南京信息工程大学 | Motion sea wave image registration fusion method based on linear transformation and wavelet analysis |
| WO2022100065A1 (en) * | 2020-11-11 | 2022-05-19 | 广东拓斯达科技股份有限公司 | Image registration method and apparatus, electronic device and storage medium |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110222781A1 (en) * | 2010-03-15 | 2011-09-15 | U.S. Government As Represented By The Secretary Of The Army | Method and system for image registration and change detection |
| CN103839265A (en) * | 2014-02-26 | 2014-06-04 | 西安电子科技大学 | SAR image registration method based on SIFT and normalized mutual information |
| CN104751453A (en) * | 2015-03-11 | 2015-07-01 | 西安电子科技大学 | SAR (Synthetic Aperture Radar) image changing and detecting method based on CUDA (Compute Unified Device Architecture) and steady WT (Wavelet Transform) |
-
2017
- 2017-03-17 CN CN201710160992.8A patent/CN106997602B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110222781A1 (en) * | 2010-03-15 | 2011-09-15 | U.S. Government As Represented By The Secretary Of The Army | Method and system for image registration and change detection |
| CN103839265A (en) * | 2014-02-26 | 2014-06-04 | 西安电子科技大学 | SAR image registration method based on SIFT and normalized mutual information |
| CN104751453A (en) * | 2015-03-11 | 2015-07-01 | 西安电子科技大学 | SAR (Synthetic Aperture Radar) image changing and detecting method based on CUDA (Compute Unified Device Architecture) and steady WT (Wavelet Transform) |
Non-Patent Citations (6)
| Title |
|---|
| YANG HUANG等: "Research on CUDA-Based SIFT Registration of SAR Image", 《2011 FOURTH INTERNATIONAL SYMPOSIUM ON PARALLEL ARCHITECTURES, ALGORITHMS AND PROGRAMMING》 * |
| 刘晓慧等: "基于归一化互信息和金字塔分解优化的图像配准算法研究", 《微型电脑应用》 * |
| 强赞霞等: "基于互信息的分层遥感图像配准方法", 《计算机工程与应用》 * |
| 李乔亮等: "基于样条金字塔和互信息的快速图像配准", 《计算机应用研究》 * |
| 杨纪伟: "基于GPU的SAR图像处理并行算法及实现", 《中国优秀硕士学位论文全文数据库（电子期刊）》 * |
| 赵飞: "图像配准算法研究与实现", 《中国优秀硕士学位论文全文数据库（电子期刊）》 * |
Cited By (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108742678A (en) * | 2018-06-01 | 2018-11-06 | 妙智科技（深圳）有限公司 | Method for registering images, device and computer readable storage medium |
| CN108742678B (en) * | 2018-06-01 | 2022-02-18 | 妙智科技（深圳）有限公司 | Image registration method and device and computer-readable storage medium |
| CN111445503A (en) * | 2020-03-25 | 2020-07-24 | 桂林电子科技大学 | Pyramid mutual information image registration method based on parallel programming model on GPU cluster |
| CN111445503B (en) * | 2020-03-25 | 2023-04-25 | 桂林电子科技大学 | Pyramid mutual information image registration method based on parallel programming model on GPU cluster |
| WO2022100065A1 (en) * | 2020-11-11 | 2022-05-19 | 广东拓斯达科技股份有限公司 | Image registration method and apparatus, electronic device and storage medium |
| CN112862868A (en) * | 2021-01-31 | 2021-05-28 | 南京信息工程大学 | Motion sea wave image registration fusion method based on linear transformation and wavelet analysis |
| CN112862868B (en) * | 2021-01-31 | 2023-12-01 | 南京信息工程大学 | Motion sea wave image registration fusion method based on linear transformation and wavelet analysis |
Also Published As
| Publication number | Publication date |
|---|---|
| CN106997602B (en) | 2019-08-06 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106997602B (en) | SAR image registration method based on GPU and pyramid mutual information | |
| CN104616345B (en) | Octree forest compression based three-dimensional voxel access method | |
| CN105453139B (en) | Sparse GPU voxelization for 3D resurfacing | |
| CN106683173A (en) | Method of improving density of three-dimensional reconstructed point cloud based on neighborhood block matching | |
| Kropatsch et al. | Digital image analysis: selected techniques and applications | |
| CN110728671B (en) | Dense reconstruction method of texture-free scene based on vision | |
| CN104299228B (en) | A kind of remote sensing image dense Stereo Matching method based on Accurate Points position prediction model | |
| CN110838115B (en) | Ancient cultural relic three-dimensional model change detection method by contour line extraction and four-dimensional surface fitting | |
| CN1860522A (en) | System and method for real-time co-rendering of multiple attributes | |
| Liu et al. | Characterization of pore structure parameters of foam concrete by 3D reconstruction and image analysis | |
| CN106023147B (en) | The method and device of DSM in a kind of rapidly extracting linear array remote sensing image based on GPU | |
| CN103544710A (en) | Image registration method | |
| CN108010123A (en) | A kind of three-dimensional point cloud acquisition methods for retaining topology information | |
| CN107689079A (en) | The cloudland method for reconstructing that a kind of satellite cloud picture is combined with natural image | |
| CN110473233A (en) | Method for registering, computer equipment and storage medium | |
| CN107341837A (en) | Grid and vector data conversion and continuous yardstick expression based on image pyramid | |
| Novak et al. | Approximate registration of point clouds with large scale differences | |
| Peng et al. | Development of a framework for stereo image retrieval with both height and planar features | |
| CN107610161B (en) | A kind of InSAR fast image registration method based on Quadtree Partition | |
| CN116246017A (en) | Sea wave three-dimensional reconstruction method based on binocular image data | |
| CN113688816B (en) | Calculation method of visual odometer for improving ORB feature point extraction | |
| Wang et al. | Utilizing parallel networks to produce sub-pixel shifted images with multiscale spatio-spectral information for soft-then-hard sub-pixel mapping | |
| Han-dan et al. | Research on real-time visualization of large-scale 3D terrain | |
| Han et al. | Three-dimensional visualization of soil pore structure using computed tomography | |
| Bing et al. | Visual imaging method of 3D virtual scene based on VR technology |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |