contou.bas

采用三角化的方法基于mapinfo的等值线算法例子。

    ReDim Xcontou(0 To N0), Ycontou(0 To N0), Zcontou(0 To N0)
    
    NContou = -1
    DX = (Xmax - Xmin) / (NX - 1)
    DY = (Ymax - Ymin) / (NY - 1)
    Yt = Ymin - DY
    Vmin = 10000000000#
    Vmax = -Vmin
    If (bSpace = False) Then
        For IY = 1 To NY
            Yt = Yt + DY
            Xt = Xmin - DX
            For IX = 1 To NX
                Xt = Xt + DX
                NContou = NContou + 1
                Input #1, Zcontou(NContou)
                Xcontou(NContou) = Xt
                Ycontou(NContou) = Yt
                Zgrid(IX, IY) = Zcontou(NContou)
                
                If (Zcontou(NContou) < Vmin) Then Vmin = Zcontou(NContou)
                If (Zcontou(NContou) > Vmax) Then Vmax = Zcontou(NContou)
            Next IX
        Next IY
    Else
        For IY = 1 To NY
            Yt = Yt + DY
            Xt = Xmin - DX
            For IX = 1 To NX
                Xt = Xt + DX
                NContou = NContou + 1
                Input #1, Zcontou(NContou)
                Xcontou(NContou) = Xt
                Ycontou(NContou) = Yt
                Zgrid(IX, IY) = Zcontou(NContou)
                
                If (Zcontou(NContou) < Vmin) Then Vmin = Zcontou(NContou)
                If (Zcontou(NContou) > Vmax) Then Vmax = Zcontou(NContou)
            Next IX
            Line Input #1, Temp
        Next IY
    End If
    Close (1)
ElseIf (DSAA_DSBB = "DSRB") Then
    DataType = 0
    Open TheContouPath For Binary Access Read Lock Read As #1
    Seek #1, 17
    Get #1, , lNY
    Get #1, , lNY
    Get #1, , lNX
    Get #1, , XminT
    Get #1, , YminT
    Get #1, , xStep
    Get #1, , yStep
    Get #1, , VminT
    Get #1, , VmaxT
        
    NX = lNX
    NY = lNY
    DX = xStep
    DY = yStep
    Xmin = XminT
    Xmax = Xmin + (NX - 1) * xStep
    Ymin = YminT
    Ymax = Ymin + (NY - 1) * yStep
    Vmin = VminT
    Vmax = VmaxT
    
    ReDim Zgrid(1 To NX, 1 To NY), bZgrid(1 To NX, 1 To NY)
    N0 = NX * NY
    ReDim Xcontou(0 To N0), Ycontou(0 To N0), Zcontou(0 To N0)
    
    Seek #1, 101
    NContou = -1
    Yt = Ymin - DY
    For IY = 1 To NY
        Yt = Yt + DY
        Xt = Xmin - DX
        For IX = 1 To NX
            Xt = Xt + DX
            NContou = NContou + 1
            Get #1, , DouValueTemp
            
            Xcontou(NContou) = Xt
            Ycontou(NContou) = Yt
            Zcontou(NContou) = DouValueTemp
            Zgrid(IX, IY) = DouValueTemp
        Next IX
    Next IY
    Close (1)
ElseIf (DSAA_DSBB = "DSBB") Then
    DataType = 0
    Open TheContouPath For Binary Access Read Lock Read As #1
    Seek #1, 5
    Get #1, , NX
    Get #1, , NY
    Get #1, , XminT
    Get #1, , XmaxT
    Get #1, , YminT
    Get #1, , YmaxT
    Get #1, , VminT
    Get #1, , VmaxT
    Xmin = XminT
    Xmax = XmaxT
    Ymin = YminT
    Ymax = YmaxT
    Vmin = VminT
    Vmax = VmaxT
    
    ReDim Zgrid(1 To NX, 1 To NY), bZgrid(1 To NX, 1 To NY)
    N0 = NX
    N0 = N0 * NY
    ReDim Xcontou(0 To N0), Ycontou(0 To N0), Zcontou(0 To N0)
    
    NContou = -1
    DX = (Xmax - Xmin) / (NX - 1)
    DY = (Ymax - Ymin) / (NY - 1)
    Yt = Ymin - DY
    For IY = 1 To NY
        Yt = Yt + DY
        Xt = Xmin - DX
        For IX = 1 To NX
            Xt = Xt + DX
            NContou = NContou + 1
            Get #1, , ValueTemp
            
            Xcontou(NContou) = Xt
            Ycontou(NContou) = Yt
            Zcontou(NContou) = ValueTemp
            Zgrid(IX, IY) = ValueTemp
        Next IX
    Next IY
    Close (1)
Else '一般格式Y,X,Z
    N0 = 1000
    ReDim Xcontou(0 To N0), Ycontou(0 To N0), Zcontou(0 To N0)
    N = -1
    Open TheContouPath For Input As #1
    Do While Not EOF(1)
        N = N + 1
        If (N > N0) Then
            N0 = N0 + 100
            ReDim Preserve Xcontou(0 To N0), Ycontou(0 To N0), Zcontou(0 To N0)
        End If
        Input #1, Ycontou(N), Xcontou(N), Zcontou(N)
    Loop
    Close (1)
    For I = N To 1 Step -1
        If (Ycontou(N) + Xcontou(N) < 1) Then
            N = N - 1
        Else
            Exit For
        End If
    Next I
    NContou = N
    Call SortXY(Xcontou, Ycontou, Zcontou, NContou)
    
    '判断是否是网格数据
    Call CheckGrid(Xcontou, Ycontou, NContou, NX, NY)
    
    If (DataType = 0) Then
        ReDim Zgrid(1 To NX, 1 To NY), bZgrid(1 To NX, 1 To NY)
        II = -1
        For IX = 1 To NX
            For IY = 1 To NY
                II = II + 1
                Zgrid(IX, IY) = Zcontou(II)
            Next IY
        Next IX
    End If
    
    Xmin = 10000000000#
    Xmax = -Xmin
    Ymin = 10000000000#
    Ymax = -Ymin
    Vmin = 10000000000#
    Vmax = -Vmin
        
    For I = 0 To NContou
        If (Xcontou(I) < Xmin) Then Xmin = Xcontou(I)
        If (Xcontou(I) > Xmax) Then Xmax = Xcontou(I)
        
        If (Ycontou(I) < Ymin) Then Ymin = Ycontou(I)
        If (Ycontou(I) > Ymax) Then Ymax = Ycontou(I)
    
        If (Zcontou(I) > Vmax) Then Vmax = Zcontou(I)
        If (Zcontou(I) < Vmin) Then Vmin = Zcontou(I)
    Next I
    
    'Begin判断经纬度是否颠倒
    If (Ymin >= -180 And Ymax <= 180 And Xmin >= -90 And Xmax <= 90) Then
        ReDim Xcontou(0 To NContou), Ycontou(0 To NContou), Zcontou(0 To NContou)
        Open TheContouPath For Input As #1
        For I = 0 To NContou
            Input #1, Xcontou(I), Ycontou(I), Zcontou(I)
        Next I
        Close (1)
        
        Call SortXY(Xcontou, Ycontou, Zcontou, NContou)
    
        '判断是否是网格数据
        Call CheckGrid(Xcontou, Ycontou, NContou, NX, NY)
    
        If (DataType = 0) Then
            ReDim Zgrid(1 To NX, 1 To NY), bZgrid(1 To NX, 1 To NY)
            II = -1
            For IX = 1 To NX
                For IY = 1 To NY
                    II = II + 1
                    Zgrid(IX, IY) = Zcontou(II)
                Next IY
            Next IX
        End If
    
        Xmin = 10000000000#
        Xmax = -Xmin
        Ymin = 10000000000#
        Ymax = -Ymin
        Vmin = 10000000000#
        Vmax = -Vmin
        
        For I = 0 To NContou
            If (Xcontou(I) < Xmin) Then Xmin = Xcontou(I)
            If (Xcontou(I) > Xmax) Then Xmax = Xcontou(I)
        
            If (Ycontou(I) < Ymin) Then Ymin = Ycontou(I)
            If (Ycontou(I) > Ymax) Then Ymax = Ycontou(I)
    
            If (Zcontou(I) > Vmax) Then Vmax = Zcontou(I)
            If (Zcontou(I) < Vmin) Then Vmin = Zcontou(I)
        Next I
    End If
    'End判断经纬度是否颠倒
End If
End Sub
'InNum-输入数,nDec-小数位数,nWidth-宽度,IC-左右对齐标志,OutStr-输出字符串
Private Function Formats(InNum As Variant, cFormat As String) As String
Dim I As Integer

I = Len(cFormat)
Formats = Format(InNum, cFormat)
If (Len(Formats) < I) Then Formats = Space(I - Len(Formats)) + Formats
End Function

'设置等值线值
Public Sub SetContouValue(ParValue() As Single)
Dim I As Integer, J As Integer, Y As Single
Dim StrValue As String

If (nDec = 0) Then
    FMT = "#####0"
Else
    FMT = "#####0."
    For I = 1 To nDec
        FMT = FMT + "0"
    Next I
End If
'Begin求最大标注长度
If (Len(Format(Vmin, FMT)) > Len(Format(Vmax, FMT))) Then
    StrValue = Trim(Format(Vmin, FMT))
Else
    StrValue = Trim(Format(Vmax, FMT))
End If
FillScaleWidth = Len(StrValue)
'End求最大标注长度

ValueMin = Vmin - Vd
J = -1
Do
    J = J + 1
    Y = Val(Format(Vmin + J * Vd, FMT))
    ParValue(J) = Y
    If (Y >= Vmax) Then Exit Do
Loop
Vn = J - 1

End Sub
'数据按X、Y排序
Public Sub SortXY(X() As Single, Y() As Single, Z() As Single, N As Long)
Dim N1 As Long, N2 As Long, X1 As Single
Dim I As Long, J As Long

'按X坐标排序
Call ShellSort(X, Y, Z, N, 0, N)
'按Y坐标排序
X1 = X(0)
N1 = 0
N2 = 0
For I = 1 To N
    If (X(I) = X1) Then
        N2 = N2 + 1
    Else
        If (N2 > N1) Then '相同X按Y坐标排序
            Call ShellSort(Y, X, Z, N, N1, N2)
        End If
        X1 = X(I)
        N1 = I
        N2 = I
    End If
Next I
If (N2 > N1) Then '相同X按Y坐标排序
    Call ShellSort(Y, X, Z, N, N1, N2)
End If

'平均重合点
J = -1
X1 = Z(0)
N1 = 1
For I = 1 To N
    If (Abs(X(I) - X(I - 1)) + Abs(Y(I) - Y(I - 1)) < 0.00001) Then
        X1 = X1 + Z(I)
        N1 = N1 + 1
    Else
        J = J + 1
        Z(J) = X1 / N1
        X(J) = X(I - 1)
        Y(J) = Y(I - 1)
        X1 = Z(I)
        N1 = 1
    End If
Next I
J = J + 1
Z(J) = X1 / N1
X(J) = X(N)
Y(J) = Y(N)
N = J
End Sub

'排序子程序
Private Sub ShellSort(X() As Single, Y() As Single, Z() As Single, N As Long, N1 As Long, N2 As Long)
Dim B As Long, M As Long, L As Long, I As Long, J As Long, K As Long
Dim Xtemp As Single, Ytemp As Single, Ztemp As Single, jk As Long

B = N2 - N1 + 1
B = Log(B) / Log(2#)
M = B
L = 2 ^ M
For I = 1 To M
    K = L - 1
    L = L / 2
    For J = K + N1 To N2
        Xtemp = X(J)
        Ytemp = Y(J)
        Ztemp = Z(J)
        jk = J - K
        Do While (jk > N1 - 1 And X(jk) > Xtemp)
            X(jk + K) = X(jk)
            Y(jk + K) = Y(jk)
            Z(jk + K) = Z(jk)
            jk = jk - K
            If (jk < 0) Then Exit Do
        Loop
        X(jk + K) = Xtemp
        Y(jk + K) = Ytemp
        Z(jk + K) = Ztemp
    Next J
Next I
End Sub
'判断是否为网格数据
Public Sub CheckGrid(X() As Single, Y() As Single, N As Long, NX As Integer, NY As Integer)
Dim IX As Long, IY As Long
Dim RR() As Single, VV() As Single, V As Single, C As Single
Dim Xtemp As Single, Ytemp As Single, Error1 As Integer
Dim I As Long

DataType = 0
'判断是否是网格数据
Xtemp = X(0)
Ytemp = Y(0)
NX = 0
NY = 0
DY = Y(2) - Y(1)
For I = 0 To N
    If (X(I) <> Xtemp) Then Exit For
    NY = NY + 1
Next I
If (NY < 3) Then 'X方向数据小于3个,肯定非网格数据
    DataType = 1
    GoTo Error1
End If
NX = Fix((N + 1) / NY)
If (NY * NX = N + 1) Then '有可能是网格数据
    DataType = 0
    '判断X坐标是否等间距
    DX = X(NY) - X(0)
    For IX = 2 To NX - 1
        Xtemp = X(IX * NY) - X((IX - 1) * NY)
        If (Abs(Xtemp - DX) > 0.000001) Then
            DataType = 1
            GoTo Error1
        End If
    Next IX
    '判断Y是否等间距
    DY = Y(1) - Y(0)
    For IY = 2 To NY - 1
        Ytemp = Y(IY) - Y(IY - 1)
        If (Abs(Ytemp - DY) > 0.000001) Then
            DataType = 1
            GoTo Error1
        End If
    Next IY
Else '肯定不是网格数据
    DataType = 1
End If
Error1:
End Sub
Private Sub OutMIFMIDHeader2D(Xmin As Single, Xmax As Single, Ymin As Single, Ymax As Single)
    Dim Columns() As String, ColumnsType() As String, ColumnsN As Integer
    Dim SymbolShape As Integer, SymbolColor As Long, SymbolSize As Integer
    Dim I As Integer, FileName As String
    Dim X0 As Double, Y0 As Double
    Dim MaxMarK As Integer
    Dim ValueMin As Double, Delta As Double, nDec As Integer, iColor As Integer
    Dim strTitle As String, strAxisT As String, strAxisY As String, StrMax As String
    Dim Yleng As Double

    PaperHeight = Screen.Height / Screen.TwipsPerPixelY
    PaperWidth = Screen.Width / Screen.TwipsPerPixelX

    Xleng = PaperWidth
    Yleng = (Ymax - Ymin) / (Xmax - Xmin) * PaperWidth
    Xminp = Xmin
    Xmaxp = Xmax
    Yminp = Ymin
    Ymaxp = Ymax
    XYFact = (Xleng / (Xmaxp - Xminp) + Yleng / (Ymaxp - Yminp)) / 2#
    WcsX0 = -Xminp * XYFact
    WcsY0 = PaperHeight + Yminp * XYFact

    I = InStr(TableName, ".")
    If (I > 0) Then
        TableNameT = Left(TableName, I - 1)
    Else
        TableNameT = TableName
    End If
    Call CheckTabName(TableNameT, "C")

    ColumnsN = 3
    ReDim Columns(1 To ColumnsN), ColumnsType(1 To ColumnsN)

    ColumnsN = 2
    Columns(1) = "等值线值"
    Columns(2) = "等值线值长度"
    
    TheMapInfoPath = App.Path + "\sssss"
    
    FileName = TheMapInfoPath + TableNameT + "_线条图"
   '' MsgBox "aa111" & FileName
    ColumnsType(1) = "Float"
    ColumnsType(2) = "Float"
    ColumnsType(3) = "SmallInt"
    
    Call MIFMID_Open(FileName, Columns, ColumnsType, ColumnsN, 0, 0, 0, 0)
    
End Sub
Private Function StringValue(value As Single) As String
Dim I As Integer, FMT As String

If (nDec = 0) Then
    FMT = "#####0"
Else
    FMT = "#####0."
    For I = 1 To nDec
        FMT = FMT + "0"
    Next I
End If
StringValue = Trim(Format(value, FMT))
End Function