一维影像匹配运行主窗体.frm

摄影测量影象匹配

a(d9 + 1, 1) = 0
a(d9 + 1, 2) = 0
a(d9 + 1, 3) = 0
a(d9 + 1, 4) = 1
a(d9 + 1, 5) = xf1(i)
a(d9 + 1, 6) = yf1(i)
d9 = d9 + 2
Next i

d9 = 1
For i = 1 To nfs
l(d9, 1) = xf(i)
d9 = d9 + 1
l(d9, 1) = yf(i)
d9 = d9 + 1
Next i

ReDim mtxAt(6, 2 * nfs)
Call MTrans(2 * nfs, 6, a(), mtxAt()) '求转置距阵
Call MMul(6, 2 * nfs, 6, mtxAt(), a(), mtxC1())
Call MRinv(6, mtxC1())

ReDim chx(6, 2 * nfs)
Call MMul(6, 6, 2 * nfs, mtxC1(), mtxAt(), chx())

ReDim gz(6, 1)
Call MMul(6, 2 * nfs, 1, chx(), l(), gz())

Open "变换参数.txt" For Output As #1
Print #1, gz(1, 1)
Print #1, gz(2, 1)
Print #1, gz(3, 1)
Print #1, gz(4, 1)
Print #1, gz(5, 1)
Print #1, gz(6, 1)
Close #1
MsgBox ("已经转换完毕")
End Sub



Private Sub Command7_Click() '用于转换量测坐标到像平面坐标

Dim nf As Integer, xf() As Double, yf() As Double
Dim a(6, 1) As Double, kf As Integer, i As Integer
Dim xz() As Double, yz() As Double

If Option1.Value = True Then
nf = InputBox("输入需要转化的像点坐标的个数第一张")
ReDim xf(nf)
ReDim yf(nf)
Open "相片坐标1.txt " For Input As #1
For i = 1 To nf
Input #1, xf(i)
Input #1, yf(i)
Next i
Close #1
ElseIf Option2.Value = True Then
nf = InputBox("输入需要转化的像点坐标的个数第二张")
ReDim xf(nf)
ReDim yf(nf)
Open "相片坐标2.txt " For Input As #1
For i = 1 To nf
Input #1, xf(i)
Input #1, yf(i)
Next i
Close #1
Else
kf = 0
End If
ReDim xz(nf) As Double
ReDim yz(nf) As Double
Open "变换参数.txt" For Input As #1
For i = 1 To 6
Input #1, a(i, 1)
Next i
Close #1

For i = 1 To nf
xz(i) = a(1, 1) + a(2, 1) * xf(i) + a(3, 1) * yf(i)
yz(i) = a(4, 1) + a(5, 1) * xf(i) + a(6, 1) * yf(i)
Next i

If Option1.Value = True Then
Open "相片坐标1.txt " For Output As #1
For i = 1 To nf
Print #1, xz(i), yz(i)
Next i
Close #1
ElseIf Option2.Value = True Then
Open "相片坐标2.txt " For Output As #1
For i = 1 To nf
Print #1, xz(i), yz(i)
Next i
Close #1
Else
kf = 0
End If
MsgBox ("已经转换完毕")
End Sub

Private Sub Command4_Click()
Form1.Hide
Form4.Show
End Sub


Private Sub Command5_Click()

Call jiaohanshu(q1, w1, k1, a1, a2, a3, b1, b2, b3, c1, c2, c3)
Call jiaohanshu(q2, w2, k2, a11, a21, a31, b11, b21, b31, c11, c21, c31)
Command5.Enabled = False
Command1.Enabled = True

End Sub

Private Sub Form_Load()
Form8.Show
Form1.Hide
End Sub


Private Sub Command1_Click()
Form1.Hide
Form2.Show
End Sub


Private Sub Command2_Click()
Dim guodu1 As Integer, guodu2 As Integer, guodu3 As Integer
Dim nmz As Integer  '窗口的中心行列号
Dim i As Integer, j As Integer
Dim nbiaol As Integer, nbiaoh As Integer '每个水平像素单元在倾斜影象中相对于左上角的行列号

Dim caix1 As Single, caiy1 As Single   '采样点周围的四个像元素的中心坐标的数值
Dim caix2 As Single, caiy2 As Single
Dim caix3 As Single, caiy3 As Single
Dim caix4 As Single, caiy4 As Single

Dim hui1 As Single, hui2 As Single  '采样点周围的四个像素的灰度的数值
Dim hui3 As Single, hui4 As Single

Dim x1 As Single, y1 As Single     '采样点和周围像素相对平面位置关系

Dim k1 As Integer, k2 As Integer   '定义的中间过度性变量
Dim nqk As Integer, qkt As Single                 '用于对所计算的相关系数的个数以及用到的行参

Dim prgn As Integer                '用于控制进度条

Dim w As Integer                   '对相关系数的大小排序定义的有些变量
Dim imax As Single, max As Single

Dim ybjx1t As Single               '代替右相片的左边界x的数值
Dim qkjdu() As Single              '对qk协方差的值求绝对值
Dim qjn As Integer                '用于求绝对值循环的变量
Dim cun As Integer   '用于保存粗匹配求得的目标点的纵坐标的数值
Dim nprgn As Integer

'用于进行对求取的最大相关系数点进行抛物线的拟合来提高精度变量的定义
Dim fs As Double, fsa As Double, fsb As Double, fsc As Double
Dim fsa0 As Double, fsb0 As Double, fsc0 As Double
Dim fsl() As Double, fsx() As Double
Dim fsjian As Integer, bianti As Long   '过渡型变量
Dim fsxt() As Double, fsxtf() As Double, fsguo1() As Double
Dim vxg() As Double

nqk = 0
'''''求解右影象在其对应的水平影象中的x边界的坐标的数值
If bianbiaox = 1 Then
Call x_xt(ybjx1, ybjy1, a11, b11, c11, a21, b21, c21, a31, b31, c31)
Call x_xt(ybjx2, ybjy2, a11, b11, c11, a21, b21, c21, a31, b31, c31)
ybjx1t = ybjx1
Else
ybjx1 = xys + jiange / 2
ybjy1 = yys - jiange / 2
Call x_xt(ybjx1, ybjy1, a11, b11, c11, a21, b21, c21, a31, b31, c31)
ybjx2 = xyx - jiange / 2
ybjy2 = yyx + jiange / 2
Call x_xt(ybjx2, ybjy2, a11, b11, c11, a21, b21, c21, a31, b31, c31)
ybjx1t = ybjx1
End If

nprgn = Int((ybjx2 - (ybjx1 + (nm - 1) * jiange1)) / jiange1)
ProgressBar1.max = nprgn
ProgressBar1.Value = 0
'''''''''''''''''''''''''''''相关系数法搜索循环''''''''''''''''''''''''''''''''''''''''''''''''''''''
Do While ybjx1 < ybjx2 - (nm - 1) * jiange1
ReDim yq(nm, nm)
ReDim xyp(nm, nm)
ReDim yyp(nm, nm)
''''''计算右影象搜索窗口各个像素所对应的其水平象平面坐标
nqk = nqk + 1
ReDim Preserve qk(nqk)
guodu1 = 0
For i = 1 To nm
For j = 1 To nm
xyp(j, i) = ybjx1 + guodu1 * jiange1
Next j
guodu1 = guodu1 + 1
Next i

nmz = Int(nm / 2) + 1
guodu2 = 1
For i = nmz + 1 To nm
For j = 1 To nm
yyp(i, j) = zmyp - guodu2 * jiange1
Next j
guodu2 = guodu2 + 1
Next i

guodu3 = 1
For i = nmz - 1 To 1 Step -1
For j = 1 To nm
yyp(i, j) = zmyp + guodu3 * jiange1
Next j
guodu3 = guodu3 + 1
Next i

For i = 1 To nm
yyp(nmz, i) = zmyp
Next i
''''''对所建立的水平窗口中各个像素进行灰度重采样
ReDim yq(nm, nm)
For i = 1 To nm
For j = 1 To nm
Call xt_x(xyp(i, j), yyp(i, j), a11, b11, c11, a21, b21, c21, a31, b31, c31)
yyp(i, j) = yyp(i, j) + 1.863495
If Int((xyp(i, j) - (xys + jiange / 2)) / jiange) <> (xyp(i, j) - (xys + jiange / 2)) / jiange Or _
Int((yyp(i, j) - (yys + jiange / 2)) / jiange) <> (yyp(i, j) - (yys + jiange / 2)) / jiange Then

nbiaol = Int((xyp(i, j) - xys) / jiange) + 1
nbiaoh = Int((yys - yyp(i, j)) / jiange) + 1

caix2 = xys + nbiaol * jiange - jiange / 2
caiy2 = yys - nbiaoh * jiange + jiange / 2

caix1 = caix2 - jiange
caiy1 = caiy2

caix3 = caix2
caiy3 = caiy2 + jiange

caix4 = caix1
caiy4 = caiy3

huidul = Form3.Picture1.Point(caix1, caiy1)
Call huidujisuan(huidul, huidu)
hui1 = huidu

huidul = Form3.Picture1.Point(caix2, caiy2)
Call huidujisuan(huidul, huidu)
hui2 = huidu

huidul = Form3.Picture1.Point(caix3, caiy3)
Call huidujisuan(huidul, huidu)
hui3 = huidu

huidul = Form3.Picture1.Point(caix4, caiy4)
Call huidujisuan(huidul, huidu)
hui4 = huidu
x1 = xyp(i, j) - caix1
y1 = caiy2 - yyp(i, j)

Call chongcaiyang(x1, y1, hui1, hui2, hui3, hui4, yq(i, j))
Else

huidul = Form3.Picture1.Point(xyp(i, j), yyp(i, j))
Call huidujisuan(huidul, huidu)
yq(i, j) = huidu

End If
Next j
Next i
Call xiangguanxishu(zq(), yq(), qkt) '计算左右影象窗口的相关系数

qk(nqk) = qkt
ybjx1 = ybjx1 + jiange1
prgn = prgn + 1
If prgn < ProgressBar1.max Then
If prgn Mod 10 = 0 Then
ProgressBar1.Value = prgn
End If
End If
Loop
'''''''''''''''''''''''''''''''''''''''循环结束'''''''''''''''''''''''''''''''''''''''''''''''''''''
max = qk(1): imax = 1
For w = 2 To nqk    ''''''求取所求的相关系数最大的数值的下标的数值
If qk(w) > max Then
max = qk(w): imax = w
End If
Next w
'''''''''''''''''''''利用抛物线进行相关系数的拟合来精确的求取相关系数的最大值存在的地方''''''''''''''''''
If imax > 2 And imax < 5 And nqk > 7 Then
ReDim fsl(5, 1)
ReDim fsx(5, 3)
For i = -2 To 2 Step 1
fsjian = fsjian + 1
fsx(fsjian, 1) = 1
fsx(fsjian, 2) = i
fsx(fsjian, 3) = 2 * i
Next i

fsa0 = qkjdu(imax)
fsb0 = (qkjdu(imax + 1) - qkjdu(imax - 1)) / 2
fsc0 = (qkjdu(imax + 1) - 2 * qkjdu(imax) + qkjdu(imax - 1)) / 2

fsjian = 0
For i = -2 To 2 Step 1
fsjian = fsjian + 1
fsl(fsjian, 1) = qkjdu(i + imax) - (fsa0 + fsb0 * i + fsc0 * (i ^ 2))
Next i

ReDim fsxt(3, 5)
ReDim fsxtf(3, 3)
ReDim fsguo1(3, 5)
ReDim vxg(3, 1)

Call MTrans(5, 3, fsx(), fsxt())
Call MMul(3, 5, 3, fsxt(), fsx(), fsxtf())
Call MRinv(3, fsxtf())
Call MMul(3, 3, 5, fsxtf(), fsxt(), fsguo1())
Call MMul(3, 5, 1, fsguo1(), fsl(), vxg())

fsa0 = fsa0 + vxg(1, 1)
fsb0 = fsb0 + vxg(2, 1)
fsc0 = fsc0 + vxg(3, 1)

imax = imax - fsb0 / (2 * fsc0)
Else
bianti = 0
End If
'''''''''''''''''''''''''''''''''拟合结束'''''''''''''''''''''''''''''''''''''''''''''''
ymxp = ybjx1t + imax * jiange  ''''''将搜索的目标点在水平像平面中的坐标转换到其相应的倾斜像平面坐标系当中
ymyp = zmyp

Call xt_x(ymxp, ymyp, a11, b11, c11, a21, b21, c21, a31, b31, c31)
ymyp = ymyp + 1.863495
ymx = ymxp
ymy = ymyp
MsgBox "粗匹配已经完成"
Form1.Hide
Form2.Show
Form3.Show

Form3.Picture1.DrawWidth = 3 ''''''在右相片显示同名点
Form3.Picture1.Line (ymxp - 0.1, ymyp)-(ymxp + 0.1, ymyp), vbYellow
Form3.Picture1.Line (ymxp, ymyp - 0.1)-(ymxp, ymyp + 0.1), vbYellow

Form3.Picture1.ForeColor = vbMagenta