gpsjdgs.cls

gps控制网设计 gps控制网设计 gps控制网设计

VERSION 1.0 CLASS
BEGIN
  MultiUse = -1  'True
END
Attribute VB_Name = "GPSjdgs"
Attribute VB_GlobalNameSpace = False
Attribute VB_Creatable = True
Attribute VB_PredeclaredId = False
Attribute VB_Exposed = False
Option Explicit
 
Private Type Kzdxx
    ID As Integer
    Name As String
    X As Single
    Y As Single
    Mx As Single
    My As Single
    Mp As Single
    E As Single
    F As Single
    Et As Single
 End Type
 Private Type JxxL
   Nk As Integer
   Mh As Integer
   P(1 To 3) As Single
   R(1 To 3) As Single
   detaLx As Single
   detaLy As Single
   detaFx As Single
   detaFy As Single
   Ms As Single
   Mt As Single
 End Type
 Private Type Bchxx
   H1 As Integer
   H2 As Integer
 End Type
 Private Type Fwjxx
   H1 As Integer
   H2 As Integer
 End Type
 
Private Point() As Kzdxx, BaseLine() As JxxL
Private S0() As Bchxx, T0() As Fwjxx, AM0 As Single
Private d() As Single, z() As Single, B() As Single, PP(2, 2) As Single, W() As Single
Private Mh%(), Nk%(), MHD%(), N%, N0%, Np%, NX%, ND%, NS%, NS2%, NT%, MH1%(), NK4%()
Private Mt As Single, D1 As Single, D2 As Single, E As Single, s As Single, t As Double
Private R As Single, aa As Single, Bb As Single, cc As Single, Dd As Single
Private Const rou = 206265#
Private Const Pi = 3.14159265
Dim Fp3$
Dim et1() As Single
Dim fwd%, fwf%, fwm#
Public Sub Main(filename As String)
  Dim i%, k%, NK1%, NK2%
  Dim Fp1 As String, Fp2 As String, Kzwmc As String, FP4$
  Fp1 = ""
  For i = 1 To Len(filename)
    If Mid(filename, i, 1) = "." Then
      Exit For
    Else
      Fp1 = Fp1 + Mid(filename, i, 1)
    End If
 Next i
  'Fp1 = filename + ".in"
  Fp2 = Fp1 + ".OUT"
  Fp3 = Fp1 + ".scr"
  Open filename For Input As #1
  Open Fp2 For Output As #2
  
  Input #1, Kzwmc
  Print #2, Tab(28); Kzwmc & "精度与可靠性估算"
  Print #2, ""
  Input #1, N0, Np, N
  Input #1, aa, Bb, cc, Dd
  Print #2, Tab(15); "起算点个数="; N0, "总点数="; Np, "基线个数="; N
  NX = 2 * (Np - N0)
  ReDim Point(Np), BaseLine(N), B(NX), W(NX), z(NX), MH1(N), Nk(Np), MHD(NX), Mh(2 * N)
  ReDim S0(NS), T0(NT), NK4(N)
  ReDim et1(Np)
  frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 20
  Read_Number
  frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 20
  FDZ
  frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 60
  ND = MHD(NX)
  ReDim d(ND)
  Adjust
   frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 50
  Precision
 frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 100
  Kkxjs
   frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 100
  Print_Number
  frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 50
  Close #1
  Close #2
  frmkzwsj.Presbar.Value = frmkzwsj.Presbar.Value + 100
  frmkzwsj.Presbar.Value = 0
  frmkzwsj.Presbar.Visible = False
  MsgBox "计算完毕！", vbOKOnly, "提示"
End Sub

Public Sub Read_Number()
  Dim i%, j%, k%, H1%, H2%, H As Single, NL%, NK1%, NK2%, KK%
  Dim Ms As Single, Ma As Single, A1 As Single, B1 As Single, C1 As Single, D1 As Single, AB As Single
  '读点号、点名与近似坐标
  For i = 1 To Np
    Input #1, Point(i).ID, Point(i).Name, Point(i).X, Point(i).Y
  Next i
  For i = 0 To Np
   Nk(i) = 0
  Next i
  '读基线向量的两端点号，并计算各基线的协方差阵
  For i = 1 To N
    Input #1, BaseLine(i).Nk, BaseLine(i).Mh
    MH1(i) = BaseLine(i).Mh
    NK4(i) = BaseLine(i).Nk
    Call PST(BaseLine(i).Nk, BaseLine(i).Mh)
    '根据标称精度估算边长中误差和方向中误差
    Ms = aa + Bb * s / 1000
    Ma = cc + Dd / (s / 1000)
 
    '方差以厘米的平方为单位
    A1 = Cos(t) ^ 2 * (Ms / 10) ^ 2 + (s * 100 / rou) ^ 2 * Ma ^ 2 * Sin(t) ^ 2
    B1 = Sin(t) * Cos(t) * (Ms / 10) ^ 2 - (s * 100 / rou) ^ 2 * Cos(t) * Sin(t) * Ma ^ 2
    C1 = Sin(t) ^ 2 * (Ms / 10) ^ 2 + (s * 100 / rou) ^ 2 * Cos(t) ^ 2 * Ma ^ 2
    '基线向量协方差阵求逆的权阵
    AB = A1 * C1 - B1 ^ 2
    D1 = A1
    A1 = C1 / AB
    C1 = D1 / AB
    B1 = -B1 / AB
    BaseLine(i).P(1) = A1
    BaseLine(i).P(2) = B1
    BaseLine(i).P(3) = C1
  Next i
  
  '确定每点的照准点个数NK与照准点点号MH
  For i = 1 To Np
    For j = 1 To N
      If i = NK4(j) Or i = MH1(j) Then
        Nk(i) = Nk(i) + 1
      End If
    Next j
  Next i
  Nk(0) = 0
  For i = 1 To Np
    Nk(i) = Nk(i) + Nk(i - 1)
  Next i
  For k = 1 To Np
    NK1 = Nk(k - 1) + 1
    NK2 = Nk(k)
    KK = 1
    For i = NK1 To NK2
      For j = KK To N
        If NK4(j) = k Then
          Mh(i) = BaseLine(j).Mh
          KK = j + 1
          Exit For
        Else
          If MH1(j) = k Then
            Mh(i) = BaseLine(j).Nk
            KK = j + 1
            Exit For
          End If
        End If
     Next j
   Next i
 Next k
End Sub

'计算法方程主对角元素的地址数组

Private Sub FDZ()
  Dim k%, K1%, j%, M1%, Mini%, N01%, NK1%, NK2%
  N01 = N0 + 1
  MHD(0) = 0
  For k = N01 To Np
    Mini = k
    NK1 = Nk(k - 1) + 1
    NK2 = Nk(k)
    For K1 = NK1 To NK2
      If Mh(K1) < Mini And Mh(K1) > N0 Then Mini = Mh(K1)
    Next K1
    M1 = 2 * (k - Mini + 1)
    j = 2 * (k - N0)
    MHD(j - 1) = MHD(j - 2) + M1 - 1
    MHD(j) = MHD(j - 1) + M1
 Next k
 
End Sub

'组成法方程系数数组
Private Sub Adjust()
 Dim VV1 As Single, VV2 As Single, NK1%, NK2%, D11 As Single, D22 As Single
 Dim i%, j%, k%, KI%, FS As Single, IM%, Lx As Single, Ly As Single, ID%, NTT%, K1%, K2%, K3%, K4%
 Call P00(W, NX)
 Call P00(d, ND)
 E = 0#
 For k = 1 To Np
   NK1 = Nk(k - 1) + 1
   NK2 = Nk(k)
   K1 = 2 * (k - N0) - 1
   For i = NK1 To NK2
     K2 = 2 * (Mh(i) - N0) - 1
     KI = 2 * (k - Mh(i))
     For j = 1 To N
       If k = BaseLine(j).Nk And Mh(i) = BaseLine(j).Mh Then
         If k > N0 Then
            d(MHD(K1)) = d(MHD(K1)) + BaseLine(j).P(1)
            d(MHD(K1 + 1) - 1) = d(MHD(K1 + 1) - 1) + BaseLine(j).P(2)
            d(MHD(K1 + 1)) = d(MHD(K1 + 1)) + BaseLine(j).P(3)
            If k > Mh(i) And Mh(i) > N0 Then
              d(MHD(K1) - KI) = d(MHD(K1) - KI) - BaseLine(j).P(1)
              d(MHD(K1) - KI + 1) = d(MHD(K1) - KI + 1) - BaseLine(j).P(2)
              d(MHD(K1 + 1) - KI - 1) = d(MHD(K1 + 1) - KI - 1) - BaseLine(j).P(2)
              d(MHD(K1 + 1) - KI) = d(MHD(K1 + 1) - KI) - BaseLine(j).P(3)
            End If
         End If
         If Mh(i) > N0 Then
            d(MHD(K2)) = d(MHD(K2)) + BaseLine(j).P(1)
            d(MHD(K2 + 1) - 1) = d(MHD(K2 + 1) - 1) + BaseLine(j).P(2)
            d(MHD(K2 + 1)) = d(MHD(K2 + 1)) + BaseLine(j).P(3)
            If Mh(i) > k And k > N0 Then
              KI = 2 * (Mh(i) - k)
              d(MHD(K2) - KI) = d(MHD(K2) - KI) - BaseLine(j).P(1)
              d(MHD(K2) - KI + 1) = d(MHD(K2) - KI + 1) - BaseLine(j).P(2)
              d(MHD(K2 + 1) - KI - 1) = d(MHD(K2 + 1) - KI - 1) - BaseLine(j).P(2)
              d(MHD(K2 + 1) - KI) = d(MHD(K2 + 1) - KI) - BaseLine(j).P(3)
            End If
         End If
         Exit For
       End If
     Next j
    Next i
  Next k
End Sub

  '评定精度
 Private Sub Precision()
  Dim k%, a As Single, Ame As Single, Amf As Single, QQQ As Single, Ha As Single, j%
  Dim N01%, i%, N1%, N2%, Q1 As Single, Q2 As Single, Q3 As Single, H1%, H2%
  For k = 1 To Np
     H1 = Point(k).ID
     If H1 > N0 Then
     H2 = 2 * (H1 - N0) - 1
      Call P00(W, NX)
      W(H2) = 1#
      Call PZZ
      Q1 = Abs(z(H2))
      Point(k).Mx = Sqr(Q1) '以厘米为单位
      Call P00(W, NX)
      W(H2 + 1) = 1#
      Call PZZ
      Q2 = z(H2)
      Q3 = Abs(z(H2 + 1))
      a = Pi / 2#
      If Q1 <> Q3 Then a = Atn(2 * Q2 / (Q1 - Q3)) / 2
      Ame = Sqr(Abs(Q1 + Q2 * Tan(a)))
      QQQ = Q1 + Q2 * Tan(a + Pi / 2)
      Amf = Sqr(Abs(QQQ))
      a = Degree(a)
      If a < 0# Then a = a + 179.596
      Point(k).My = Sqr(Q3) '以厘米为单位
      Point(k).Mp = Sqr(Q1 + Q3) '以厘米为单位
      If Ame > Amf Then
        Point(k).E = Ame '以厘米为单位
        Point(k).F = Amf '以厘米为单位
      Else
        Point(k).F = Ame '以厘米为单位
        Point(k).E = Amf '以厘米为单位
      End If
      Point(k).Et = a
   End If
  Next k
  N01 = N0 + 1
  For k = 1 To N
    N1 = BaseLine(k).Nk
    N2 = BaseLine(k).Mh
    Call PBB(N1, N2, 0)
      Call P00(W, NX)
      For j = 1 To NX
        W(j) = B(j)
      Next j
      Call PZZ
      D1 = 0#
      For j = 1 To NX
       D1 = D1 + z(j) * B(j)
      Next j
      BaseLine(k).Mt = Sqr(Abs(D1)) '以秒为单位
      Call PBB(N1, N2, 1)
      Call P00(W, NX)
      For j = 1 To NX
        W(j) = B(j)
      Next j
      Call PZZ
      D2 = 0#
      For j = 1 To NX
       D2 = D2 + z(j) * B(j)
      Next j
      BaseLine(k).Ms = Sqr(Abs(D2)) '以厘米为单位
 Next k
End Sub


'解方程
Private Sub PZZ()
  Dim k%, j%, K1%, KJ%, a As Single, c As Single, IJ%, K11%, i%
  For k = 1 To NX
    z(k) = W(k)
  Next k
  E = E + 1#
  For k = 1 To NX
    c = d(MHD(k))
    If k < NX Then
     K1 = k + 1
     For j = K1 To NX
      KJ = j - (MHD(j) - MHD(j - 1)) + 1
      If KJ <= k Then
         a = d(MHD(j) - j + k) / c
         If E < 1.1 Then
             For i = j To NX
                IJ = MHD(i) - i + j
                  If IJ > MHD(i - 1) And (IJ - j + k) > MHD(i - 1) Then
                     d(IJ) = d(IJ) - d(IJ - j + k) * a
                  End If
             Next i
          End If
          z(j) = z(j) - z(k) * a
      End If
    Next j
   End If
 Next k
 For k = 1 To NX
   K1 = NX + 1 - k
   If K1 <> NX Then
     K11 = K1 + 1
     For i = K11 To NX
       If MHD(i) - i + K1 > MHD(i - 1) Then
         z(K1) = -d(MHD(i) - i + K1) * z(i) + z(K1)
       End If
     Next i
    End If
    z(K1) = z(K1) / d(MHD(K1))
 Next k
End Sub

Private Sub P00(B() As Single, N%)
  Dim i%
  ReDim B(N)
  For i = 1 To N
    B(i) = 0
  Next i
End Sub

Private Function Radian(a As Single) As Single
  Dim Ra As Single, c As Single, FS As Single, Ib%, Ic%
  Ra = Pi / 180#
  Ib = Int(a)
  c = (a - Ib) * 100#
  Ic = Int(c)
  FS = (c - Ic) * 100#
  Radian = (Ib + Ic / 60# + FS / 3600#) * Ra
End Function

Private Function Degree(a As Single) As Single
  Dim B As Single, FS As Single, Ib%, Ia%
  a = 180# / Pi * a
  Ia = Int(a)
  B = (a - Ia) * 60
  Ib = Int(B)
  FS = (B - Ib) * 60
  Degree = a + Ib / 100# + FS / 10000#
End Function

'计算边长方位角
Private Sub PST(k%, Ka%)
  Dim D11 As Single, D22 As Single, D33 As Single, Ds2 As Single, Ds3 As Single
  D11 = Point(Ka).X - Point(k).X
  D22 = Point(Ka).Y - Point(k).Y
  s = Sqr(CDbl(D11 ^ 2 + D22 ^ 2 + D33 ^ 2))
  Ds2 = D22 / s
  Ds3 = Sqr(CDbl(1# - Ds2 ^ 2))
  If Abs(Ds3) <> 0 Then
    t = Atn(Ds2 / Ds3)