'********************************************************************************
'*   模塊名稱：cbxfxs
'*   功能描述：計算二階車貝雪夫系數
'*   設    計：江文華、陳道勁
'*   說    明：GridValue（1 To MY, 1 To MX）
'*   返回說明：cbxfData(1 to 9) cbxfData(1)=A00
'                               cbxfData(2)=A01
'                               cbxfData(3)=A10
'                               cbxfData(4)=A11
'                               cbxfData(5)=A02
'                               cbxfData(6)=A20
'                               cbxfData(7)=A12
'                               cbxfData(8)=A21
'                               cbxfData(9)=A22
'********************************************************************************


Sub CalculateCBXF(cFilein As String) 'MX 為緯向格點數,MY為經向格點數

'cFilein為輸入數據文件名
'MX 為緯向格點數,MY為經向格點數
'cbxfData 數組用于存放車貝雪夫二階系數A00、A01、A10、A11、A02、A20、A12、A21、A22
'Dim ZX() As Single, ZY() As Single
'Dim GridValue() As Single
''從數據文件中讀取數據

Dim GridValue()
Dim ZX()
Dim ZY()


Open cFilein For Input As #11
Input #11, MY, MX
ReDim GridValue(1 To MY, 1 To MX)
ReDim ZX(1 To MY)
ReDim ZY(1 To MX)
For J = 1 To MY
    For I = 1 To MX
     Input #11, GridValue(J, I)
    Next
Next
Close #11
'


Dim CI As Single, CJ As Single ', I As Long, J As Long
CI = 0
For I = 1 To MX
CI = CI + I
Next
CI = CI / MX
'
CJ = 0
For J = 1 To MY
   CJ = CJ + J
Next
CJ = CJ / MY
'
Dim M0 As Single, N0 As Single, S01 As Single, S02 As Single, S10 As Single, S20 As Single
M0 = (MX ^ 2 - 1) / 12
N0 = (MY ^ 2 - 1) / 12
'計算經向平均值
For I = 1 To MX
  ZY(I) = 0
    For J = 1 To MY
     ZY(I) = ZY(I) + GridValue(J, I)
    Next
  ZY(I) = ZY(I) / MY
Next
'計算緯向平均值
For J = 1 To MY
  ZX(J) = 0
    For I = 1 To MX
      ZX(J) = ZX(J) + GridValue(J, I)
    Next
  ZX(J) = ZX(J) / MX
Next
'計算系數A00(A00=cbxfData(1))
cbxfData(1) = 0
    For J = 1 To MY
       cbxfData(1) = cbxfData(1) + ZX(J)
    Next
cbxfData(1) = cbxfData(1) / MY
'計算系數A01(A01=cbxfData(2))
 cbxfData(2) = 0
 S01 = 0
    For J = 1 To MY
        S01 = S01 + (J - CJ) ^ 2
        cbxfData(2) = cbxfData(2) + ZX(J) * (J - CJ)
    Next
cbxfData(2) = cbxfData(2) / S01
'計算系數A10(A10=cbxfData(3))
 cbxfData(3) = 0
 S10 = 0
    For I = 1 To MX
        S10 = S10 + (I - CI) ^ 2
        cbxfData(3) = cbxfData(3) + ZY(I) * (I - CI)
    Next
cbxfData(3) = cbxfData(3) / S10
'計算系數A11(A11=cbxfData(4))
 cbxfData(4) = 0
For I = 1 To MX
   For J = 1 To MY
        cbxfData(4) = cbxfData(4) + GridValue(J, I) * (I - CI) * (J - CJ)
   Next
Next
   cbxfData(4) = cbxfData(4) / (S10 * S01)
'計算系數A02(A02=cbxfData(5))
 S02 = 0
 cbxfData(5) = 0
For J = 1 To MY
    S02 = S02 + ((J - CJ) ^ 2 - N0) ^ 2
    cbxfData(5) = cbxfData(5) + ZX(J) * ((J - CJ) ^ 2 - N0)
Next
cbxfData(5) = cbxfData(5) / S02
'計算系數A20(A20=cbxfData(6))
 cbxfData(6) = 0
 S20 = 0
For I = 1 To MX
  S20 = S20 + ((I - CI) ^ 2 - M0) ^ 2
  cbxfData(6) = cbxfData(6) + ZY(I) * ((I - CI) ^ 2 - M0)
Next
 cbxfData(6) = cbxfData(6) / S20
 '計算系數A12(A12=cbxfData(7))
  cbxfData(7) = 0
 For I = 1 To MX
   For J = 1 To MY
    cbxfData(7) = cbxfData(7) + GridValue(J, I) * (I - CI) * ((J - CJ) ^ 2 - N0)
   Next
 Next
 cbxfData(7) = cbxfData(7) / (S10 * S02)
  '計算系數A21(A21=cbxfData(8))
      cbxfData(8) = 0
 For I = 1 To MX
   For J = 1 To MY
    cbxfData(8) = cbxfData(8) + GridValue(J, I) * ((I - CI) ^ 2 - M0) * (J - CJ)
   Next
 Next
 cbxfData(8) = cbxfData(8) / (S20 * S01)
   '計算系數A22(A22=cbxfData(9))
    cbxfData(9) = 0
 For I = 1 To MX
   For J = 1 To MY
     cbxfData(9) = cbxfData(9) + GridValue(J, I) * ((I - CI) ^ 2 - M0) * ((J - CJ) ^ 2 - N0)
   Next
 Next
   cbxfData(9) = cbxfData(9) / (S20 * S02)
 Erase GridValue, ZX, ZY
End Sub