End Sub

Private Sub CFM_Click()

   SetHourglass
    
   'THIS PROGRAM IS FOR Applied General C-F MODEL PROCEDURE
   'ORIGINATED BY : CHEN YONGLIANG
   'DATE:           OCTOBER 10,2000
    
    Dim xmin As Double, xmax As Double, ymin As Double, xo As Double, yo As Double, ymax As Double, xx As Double, yy As Double, yg As Double
    Dim dx As Double, dy As Double
    Dim m As Integer, np As Integer, k As Integer, l As Integer
    Dim grid_no  As Integer, i As Integer, ii As Integer, ij As Integer, j As Integer, deposit_num As Integer, ixx As Integer, iyy As Integer
    Dim grid_co() As Double, deposit_co(1000, 2) As Double, evidence_co(10000, 2) As Double, joint_co(100) As Double
    Dim grid_dep() As Integer, point_num(1000) As Integer
    Dim cf() As Double, believe() As Double
    Dim filename As String * 80, filename1 As String * 80, filename2 As String * 80
    Dim region_num As Integer, kk As Integer, i1 As Integer, kkks As Integer, iii As Integer, jjj As Integer
    Dim kks As Integer, kkd As Integer, kke1 As Integer, kke2 As Integer, kke As Integer, kkef As Integer
    Dim xs As Double, xx1 As Double, xx2 As Double, xx3 As Double, yy1 As Double, yy2 As Double, yy3 As Double
    Dim ccc As Double, py As Double, pyz As Double, pyzf As Double, scf As Double, scf1 As Double, mincf As Double
    Dim file() As String * 80
    
    On Error GoTo errprint
    
   'The number of map layers
    m = Val(Text1.Text)
    
   'The number of scanning lines
    k = Val(Text2.Text)
    
    ReDim file(2 * m + 1)
    
    'Open *.txt file and input map layer file names
    filename = Text3.Text
    Open filename For Input As #1
    For i = 1 To 2 * m + 1
      Input #1, file(i)
    Next i
    Close #1
    
   'Open the *.mif file and input the data of deposit map layer
    filename = file(1)
    Open filename For Input As #1
    Call InputDep(deposit_co, deposit_num, 1, xmin, ymin, xmax, ymax)
    Close #1
  
   'The ratio between the width and height of the map area
    ccc = (xmax - xmin) / (ymax - ymin)
    
   'The number of columns
    l = CInt(k * ccc)

    dx = (xmax - xmin) / CDbl(l)
    dy = (ymax - ymin) / CDbl(k)
    
    ReDim grid_co(k * l, 12), grid_dep(k * l, m + 1), believe(k * l), cf(m, 3)
 
   'Generate the uniform grid cells
    For i = 1 To k
      For j = 1 To l
      
        grid_no = (i - 1) * l + j
        xx = xmin + j * dx - dx / 2#
        yy = ymin + i * dy - dy / 2#
        xx1 = xx - dx / 2#
        yy1 = yy - dy / 2#
        xx2 = xx + dx / 2#
        yy2 = yy + dy / 2#
        
       'The X_Y coordinates of the center of the cell
        grid_co(grid_no, 1) = xx
        grid_co(grid_no, 2) = yy
        
       'The four boundaries of the cell
        grid_co(grid_no, 3) = xx1
        grid_co(grid_no, 4) = yy1
        grid_co(grid_no, 5) = xx2
        grid_co(grid_no, 6) = yy1
        grid_co(grid_no, 7) = xx2
        grid_co(grid_no, 8) = yy2
        grid_co(grid_no, 9) = xx1
        grid_co(grid_no, 10) = yy2
        grid_co(grid_no, 11) = xx1
        grid_co(grid_no, 12) = yy1
        
      Next j
    Next i
  
   'Determining grid cells whether containning a deposit
    For i = 1 To deposit_num
      xx = deposit_co(i, 1)
      yy = deposit_co(i, 2)
      ixx = Fix((xx - xmin) / dx) + 1
      iyy = Fix((yy - ymin) / dy) + 1
      ii = (iyy - 1) * l + ixx
      grid_dep(ii, 1) = 1
    Next i
   
   'Initialize the grid_dep
    For i = 1 To k * l
      For j = 1 To m
        grid_dep(i, j + 1) = -1
      Next j
    Next i
 
   'Determining map patterns whether exist in grid cells
    For j = 1 To m
    
     'Open the map pattern file and input the data
      filename = file(2 * j)
      Open filename For Input As #2
      Call InputEvi(evidence_co, point_num, region_num, 2)
       
      For ij = 1 To k
         yg = grid_co((ij - 1) * l + 1, 2)
         kk = 0
         For i1 = 1 To region_num
            kkks = 0
            For ii = kk + 1 To kk + point_num(i1) - 1
              xx1 = evidence_co(ii, 1)
              yy1 = evidence_co(ii, 2)
              xx2 = evidence_co(ii + 1, 1)
              yy2 = evidence_co(ii + 1, 2)
              If ii < kk + point_num(i1) - 1 Then
                xx3 = evidence_co(ii + 2, 1)
                yy3 = evidence_co(ii + 2, 2)
              ElseIf ii = kk + point_num(i1) - 1 Then
                xx3 = evidence_co(kk + 2, 1)
                yy3 = evidence_co(kk + 2, 2)
              End If
              
             'Determining the intersecting points of the scanning line and the polygon boundaries
              If yy1 <> yy2 Then
                If (yg > yy1 And yg < yy2) Or (yg < yy1 And yg > yy2) Then
                  kkks = kkks + 1
                  joint_co(kkks) = (xx2 - xx1) * (yg - yy1) / (yy2 - yy1) + xx1
                ElseIf yg = yy2 Then
                  If (yy3 > yg And yy1 < yg) Or (yy3 < yg And yy1 > yg) Then
                    kkks = kkks + 1
                    joint_co(kkks) = xx2
                  ElseIf yy3 > yg And yy1 > yg Then
                    kkks = kkks + 2
                    joint_co(kkks - 1) = xx2
                    joint_co(kkks) = xx2
                  ElseIf yy3 < yg And yy1 < yg Then
                    ixx = Fix((xx2 - xmin) / dx) + 1
                    iyy = Fix((yy2 - ymin) / dy) + 1
                    grid_dep((iyy - 1) * l + ixx) = 1
                  End If
                End If
              ElseIf yy1 = yg And yy2 = yg Then
                kkks = kkks + 2
                joint_co(kkks - 1) = xx1
                joint_co(kkks) = xx2
              End If
            Next ii
       
           'Sort the X_coordinates
            For iii = 1 To kkks - 1
              For jjj = iii + 1 To kkks
                If joint_co(iii) > joint_co(jjj) Then
                   xs = joint_co(iii)
                   joint_co(iii) = joint_co(jjj)
                   joint_co(jjj) = xs
                End If
              Next jjj
            Next iii
            
           'Determining the grid_dep() values
            kks = kkks \ 2
            For iii = 1 To kks
              For jjj = 1 To l
                If grid_co((ij - 1) * l + jjj, 1) >= joint_co(2 * iii - 1) And grid_co((ij - 1) * l + jjj, 1) <= joint_co(2 * iii) Then
                  grid_dep((ij - 1) * l + jjj, j + 1) = 1
                End If
              Next jjj
            Next iii
            kk = kk + point_num(i1)
         Next i1
      Next ij
      Close #2
      
     'Open the unmapped region file and input the data
      filename = file(2 * j + 1)
      Open filename For Input As #2
      Call InputEvi(evidence_co, point_num, region_num, 2)
       
      For ij = 1 To k
         yg = grid_co((ij - 1) * l + 1, 2)
         kk = 0
         For i1 = 1 To region_num
            kkks = 0
            For ii = kk + 1 To kk + point_num(i1) - 1
              xx1 = evidence_co(ii, 1)
              yy1 = evidence_co(ii, 2)
              xx2 = evidence_co(ii + 1, 1)
              yy2 = evidence_co(ii + 1, 2)
              If ii < kk + point_num(i1) - 1 Then
                xx3 = evidence_co(ii + 2, 1)
                yy3 = evidence_co(ii + 2, 2)
              ElseIf ii = kk + point_num(i1) - 1 Then
                xx3 = evidence_co(kk + 2, 1)
                yy3 = evidence_co(kk + 2, 2)
              End If
              
             'Determining the intersecting points of the scanning line and the polygon boundaries
              If yy1 <> yy2 Then
                If (yg > yy1 And yg < yy2) Or (yg < yy1 And yg > yy2) Then
                  kkks = kkks + 1
                  joint_co(kkks) = (xx2 - xx1) * (yg - yy1) / (yy2 - yy1) + xx1
                ElseIf yg = yy2 Then
                  If (yy3 > yg And yy1 < yg) Or (yy3 < yg And yy1 > yg) Then
                    kkks = kkks + 1
                    joint_co(kkks) = xx2
                  ElseIf yy3 > yg And yy1 > yg Then
                    kkks = kkks + 2
                    joint_co(kkks - 1) = xx2
                    joint_co(kkks) = xx2
                  ElseIf yy3 < yg And yy1 < yg Then
                    ixx = Fix((xx2 - xmin) / dx) + 1
                    iyy = Fix((yy2 - ymin) / dy) + 1
                    grid_dep((iyy - 1) * l + ixx) = 0
                  End If
                End If
              ElseIf yy1 = yg And yy2 = yg Then
                kkks = kkks + 2
                joint_co(kkks - 1) = xx1
                joint_co(kkks) = xx2
              End If
            Next ii
       
           'Sort the X_coordinates
            For iii = 1 To kkks - 1
              For jjj = iii + 1 To kkks
                If joint_co(iii) > joint_co(jjj) Then
                   xs = joint_co(iii)
                   joint_co(iii) = joint_co(jjj)
                   joint_co(jjj) = xs
                End If
              Next jjj
            Next iii
            
           'Determining the grid_dep() values
            kks = kkks \ 2
            For iii = 1 To kks
              For jjj = 1 To l
                If grid_co((ij - 1) * l + jjj, 1) >= joint_co(2 * iii - 1) And grid_co((ij - 1) * l + jjj, 1) <= joint_co(2 * iii) Then
                  grid_dep((ij - 1) * l + jjj, j + 1) = 0
                End If
              Next jjj
            Next iii
            kk = kk + point_num(i1)
         Next i1
      Next ij
      Close #2
    Next j
     

  'Calculating the certainty factors for all map patterns
   kkd = 0
   For i = 1 To k * l
     If grid_dep(i, 1) = 1 Then
       kkd = kkd + 1
     End If
   Next i
   For j = 1 To m
       kke1 = 0
       kke2 = 0
       kke = 0
       kkef = 0
       For i = 1 To k * l
          If grid_dep(i, j + 1) = 1 Then
             kke1 = kke1 + grid_dep(i, 1) * grid_dep(i, j + 1)
             kke = kke + 1
          ElseIf grid_dep(i, j + 1) = -1 Then
             kke2 = kke2 + grid_dep(i, 1) * (-grid_dep(i, j + 1))
             kkef = kkef + 1
          End If
       Next i
       py = CDbl(kkd) / CDbl(k * l)
       pyz = CDbl(kke1) / CDbl(kke)
       pyzf = CDbl(kke2) / CDbl(kkef)
       If pyz >= py Then
          cf(j, 1) = (pyz - py) / ((1 - py) * pyz)
       Else
          cf(j, 1) = (pyz - py) / (py * (1 - pyz))
       End If
       If pyzf >= py Then
          cf(j, 2) = (pyzf - py) / ((1 - py) * pyzf)
       Else
          cf(j, 2) = (pyzf - py) / (py * (1 - pyzf))
       End If
   Next j
   
   'Output the certainty factors of all map patterns
    Open "c:\cfmweights.txt" For Output As #20
    Print #20, "The two estimated certainty factors:"
    For j = 1 To m
       Print #20, cf(j, 1), "  ", cf(j, 2)
    Next j
    Close #20
   
  'Synthesizing the certainty factors for all grids
   If m >= 2 Then
     For i = 1 To k * l
       If grid_dep(i, 2) = 1 And grid_dep(i, 3) = 1 Then
         If cf(1, 1) >= 0 And cf(2, 1) >= 0 Then
           scf = cf(1, 1) + cf(2, 1) - cf(1, 1) * cf(2, 1)
         ElseIf (cf(1, 1) > 0 And cf(2, 1) < 0) Or (cf(1, 1) < 0 And cf(2, 1) > 0) Then
           If Abs(cf(1, 1)) <= Abs(cf(2, 1)) Then
              mincf = Abs(cf(1, 1))
           ElseIf Abs(cf(1, 1)) > Abs(cf(2, 1)) Then
              mincf = Abs(cf(2, 1))
           End If
           scf = (cf(1, 1) + cf(2, 1)) / (1 - mincf)
         ElseIf cf(1, 1) < 0 And cf(2, 1) < 0 Then
           scf = cf(1, 1) + cf(2, 1) + cf(1, 1) * cf(2, 1)
         End If
       ElseIf grid_dep(i, 2) = 1 And grid_dep(i, 3) = -1 Then
         If cf(1, 1) >= 0 And cf(2, 2) >= 0 Then
           scf = cf(1, 1) + cf(2, 2) - cf(1, 1) * cf(2, 2)
         ElseIf (cf(1, 1) > 0 And cf(2, 2) < 0) Or (cf(1, 1) < 0 And cf(2, 2) > 0) Then
           If Abs(cf(1, 1)) <= Abs(cf(2, 2)) Then
             mincf = Abs(cf(1, 1))
           ElseIf Abs(cf(1, 1)) > Abs(cf(2, 2)) Then
             mincf = Abs(cf(2, 2))
           End If
           scf = (cf(1, 1) + cf(2, 2)) / (1 - mincf)
         ElseIf cf(1, 1) < 0 And cf(2, 2) < 0 Then
           scf = cf(1, 1) + cf(2, 2) + cf(1, 1) * cf(2, 2)
         End If
       ElseIf grid_dep(i, 2) = 1 And grid_dep(i, 3) = 0 Then
         scf = cf(1, 1)
       ElseIf grid_dep(i, 2) = -1 And grid_dep(i, 3) = 1 Then
         If cf(1, 2) >= 0 And cf(2, 1) >= 0 Then
           scf = cf(1, 2) + cf(2, 1) - cf(1, 2) * cf(2, 1)
         ElseIf (cf(1, 2) > 0 And cf(2, 1) < 0) Or (cf(1, 2) < 0 And cf(2, 1) > 0) Then
           If Abs(cf(1, 2)) <= Abs(cf(2, 1)) Then
             mincf = Abs(cf(1, 2))
           ElseIf Abs(cf(1, 2)) > Abs(cf(2, 1)) Then
             mincf = Abs(cf(2, 1))
           End If
           scf = (cf(1, 2) + cf(2, 1)) / (1 - mincf)
         ElseIf cf(1, 2) < 0 And cf(2, 1) < 0 Then
           scf = cf(1, 2) + cf(2, 1) + cf(1, 2) * cf(2, 1)
         End If
       ElseIf grid_dep(i, 2) = -1 And grid_dep(i, 3) = -1 Then
         If cf(1, 2) >= 0 And cf(2, 2) >= 0 Then
           scf = cf(1, 2) + cf(2, 2) - cf(1, 2) * cf(2, 2)
         ElseIf (cf(1, 2) > 0 And cf(2, 2) < 0) Or (cf(1, 2) < 0 And cf(2, 2) > 0) Then
           If Abs(cf(1, 2)) <= Abs(cf(2, 2)) Then
             mincf = Abs(cf(1, 2))
           ElseIf Abs(cf(1, 2)) > Abs(cf(2, 2)) Then
             mincf = Abs(cf(2, 2))
           End If
           scf = (cf(1, 2) + cf(2, 2)) / (1 - mincf)
         ElseIf cf(1, 2) < 0 And cf(2, 2) < 0 Then
           scf = cf(1, 2) + cf(2, 2) + cf(1, 2) * cf(2, 2)
         End If
       ElseIf grid_dep(i, 2) = -1 And grid_dep(i, 3) = 0 Then
           scf = cf(1, 2)
       ElseIf grid_dep(i, 2) = 0 And grid_dep(i, 3) = 1 Then
           scf = cf(2, 1)
       ElseIf grid_dep(i, 2) = 0 And grid_dep(i, 3) = -1 Then
           scf = cf(2, 2)
       ElseIf grid_dep(i, 2) = 0 And grid_dep(i, 3) = 0 Then
           scf = 0
       End If
       For j = 3 To m
         If grid_dep(i, j + 1) = 1 Then
           If scf >= 0 And cf(j, 1) >= 0 Then
             scf1 = scf + cf(j, 1) - scf * cf(j, 1)
           ElseIf (scf > 0 And cf(j, 1) < 0) Or (scf < 0 And cf(j, 1) > 0) Then