Option Explicit On
Option Strict On
Public Class Form1

    Private Sub ToolStripButton1_Click(ByVal sender As System.Object, ByVal e As System.EventArgs) Handles ToolStripButton1.Click
        Call main()
    End Sub

    Private Sub main()
        '作図用変数宣言
        Dim kpt As Integer    '拡大率（画面表示：1，ファイル保存：2）
        Dim HSIZE As Integer  '画像横寸法（px）
        Dim DHL As Integer    '画像左余白寸法：縦軸描画用（px）
        Dim DHR As Integer    '画像右余白寸法（px）
        Dim VSIZE As Integer  '画像縦寸法（px）
        Dim DVL As Integer    '画像下余白寸法：横軸描画用（px）
        Dim DVU As Integer    '画像上余白寸法（px）
        Dim sxjiku As String  'ｘ軸名（単位含む）
        Dim xmin As Double    'ｘ軸最小値
        Dim xmax As Double    'ｘ軸最大値
        Dim dx As Double      'ｘ軸増分
        Dim slog As String    'ｘ軸の種類（L：対数，N：普通）
        Dim thePicBox As PictureBox = PictureBox1
        Dim bmp As Bitmap
        Dim g As Graphics
        Dim kxi0 As Integer
        Dim kxf0 As Integer
        Dim kyi0 As Integer
        Dim kyf0 As Integer

        'データ格納用変数宣言
        Dim ndata As Integer
        Dim datax() As Double
        Dim datay() As Double
        Dim alpha As Double = 0.5 'Hazen公式

        'データ入出力用変数宣言
        Dim sr As System.IO.StreamReader
        Dim dat As String
        Dim sbuf() As String
        Dim delim() As Char = {","c}
        Dim fname1 As String = ""
        Dim fname2 As String = ""

        Dim i As Integer
        Dim j As Integer
        Dim k As Integer
        Dim work As Double

        'データ入力
        OpenFileDialog1.InitialDirectory = My.Computer.FileSystem.CurrentDirectory
        If OpenFileDialog1.ShowDialog() = Windows.Forms.DialogResult.OK Then fname1 = OpenFileDialog1.FileName
        If SaveFileDialog1.ShowDialog() = Windows.Forms.DialogResult.OK Then fname2 = SaveFileDialog1.FileName
        sr = New System.IO.StreamReader(fname1, System.Text.Encoding.Default)
        dat = sr.ReadLine()
        dat = sr.ReadLine() : sbuf = dat.Split(delim) : sxjiku = sbuf(0) 'x軸名
        dat = sr.ReadLine() : sbuf = dat.Split(delim) : HSIZE = CInt(sbuf(0)) : DHL = CInt(sbuf(1)) : DHR = CInt(sbuf(2))
        dat = sr.ReadLine() : sbuf = dat.Split(delim) : VSIZE = CInt(sbuf(0)) : DVL = CInt(sbuf(1)) : DVU = CInt(sbuf(2))
        dat = sr.ReadLine() : sbuf = dat.Split(delim)
        xmin = CDbl(sbuf(0)) : xmax = CDbl(sbuf(1)) : dx = CDbl(sbuf(2)) : slog = sbuf(3)
        dat = sr.ReadLine() : sbuf = dat.Split(delim) : ndata = CInt(sbuf(0)) - 1
        ReDim datax(ndata)
        ReDim datay(ndata)
        k = 0
        Do While sr.Peek() >= 0
            dat = sr.ReadLine() : sbuf = dat.Split(delim)
            datax(k) = CDbl(sbuf(0))
            k = k + 1
        Loop
        sr.Close()
        ndata = k - 1

        'プロット領域定義
        kxi0 = DHL
        kxf0 = HSIZE - DHR
        kyi0 = DVU
        kyf0 = VSIZE - DVL
        'ソート(大きい順)
        For i = 0 To ndata
            For j = i + 1 To ndata
                If datax(i) < datax(j) Then
                    work = datax(i)
                    datax(i) = datax(j)
                    datax(j) = work
                End If
        Next j, i

        'ｘ軸対数値の場合0データ処理
        If slog = "L" Then
            For i = 0 To ndata
                If datax(i) <= 0 Then datax(i) = 0.00001
                datax(i) = Math.Log10(datax(i))
                If datax(i) < xmin Then datax(i) = xmin
            Next i
        End If

        '超過確率計算とｙ軸データ処理
        For i = 0 To ndata
            datay(i) = (CDbl(i + 1) - alpha) / (CDbl(ndata + 1) + 1.0 - 2.0 * alpha)
            work = datay(i)
            datay(i) = TODA(work)
        Next i

        For i = 0 To 1
            Select Case i
                Case 0 '画面出力
                    kpt = 1
                    thePicBox = PictureBox1
                    thePicBox.Visible = True
                Case 1 '画像ファイル出力
                    kpt = 2
                    thePicBox = PictureBox2
                    thePicBox.Visible = False
            End Select
            thePicBox.Size = New Size(kpt * HSIZE, kpt * VSIZE)
            bmp = New Bitmap(thePicBox.Width, thePicBox.Height)
            thePicBox.Image = bmp
            g = Graphics.FromImage(thePicBox.Image)
            g.FillRectangle(Brushes.White, 0, 0, thePicBox.Width, thePicBox.Height)
            Call PLOT(g, kpt, sxjiku, xmin, xmax, dx, kxi0, kxf0, kyi0, kyf0, ndata, datax, datay, slog)
            g.Dispose()
        Next i
        PictureBox1.Left = 0
        PictureBox1.Top = ToolStrip1.Height
        Me.ClientSize = New Size(PictureBox1.Width, PictureBox1.Height + ToolStrip1.Height)

        PictureBox2.Image.Save(fname2, System.Drawing.Imaging.ImageFormat.Png)

    End Sub


    Private Sub PLOT(ByVal g As Graphics, ByVal kpt As Integer, ByVal sxjiku As String, _
                 ByVal xmin As Double, ByVal xmax As Double, ByVal dx As Double, _
                 ByVal kxi0 As Integer, ByVal kxf0 As Integer, ByVal kyi0 As Integer, ByVal kyf0 As Integer, _
                 ByVal ndata As Integer, ByRef datax() As Double, ByRef datay() As Double, ByVal slog As String)

        Dim i As Integer
        Dim ymin As Double : Dim ymax As Double
        Dim syjiku As String = "超過確率 (%)"
        Dim xx As Double : Dim yy As Double
        Dim kxx1 As Integer : Dim kyy1 As Integer
        Dim kxx2 As Integer : Dim kyy2 As Integer
        Dim kxx As Integer : Dim kyy As Integer
        Dim str As String
        Dim aa As Double
        Dim bb As Double
        Dim rr As Double
        Dim dc As Integer
        Dim f As New Font("ＭＳ ゴシック", kpt * 10)
        Dim TextSize1 As New System.Drawing.SizeF
        Dim TextSize2 As New System.Drawing.SizeF
        Dim kxi As Integer : Dim kxf As Integer : Dim kyi As Integer : Dim kyf As Integer
        Dim LPen As New System.Drawing.Pen(System.Drawing.Color.Black)

        kxi = kpt * kxi0 : kxf = kpt * kxf0 : kyi = kpt * kyi0 : kyf = kpt * kyf0

        'ｘ−ｙ軸描画
        LPen.DashStyle = Drawing2D.DashStyle.Dash
        'ｘ軸描画
        If slog = "N" Then Call DRN_XJIKU(g, f, LPen, kpt, xmin, xmax, dx, kxi, kxf, kyi, kyf)
        If slog = "L" Then Call DRL_XJIKU(g, f, LPen, kpt, xmin, xmax, kxi, kxf, kyi, kyf)
        'ｙ軸描画
        ymax = TODA(0.999)
        ymin = TODA(0.001)
        Call DRH_YJIKU(g, f, LPen, kpt, ymin, ymax, kxi, kxf, kyi, kyf)
        '枠線の描画
        g.DrawRectangle(New Pen(Color.Black, 2), kxi, kyi, kxf - kxi, kyf - kyi)
        'ｙ軸名描画
        str = syjiku
        TextSize2 = g.MeasureString(str, f)
        TextSize1 = g.MeasureString("-1.0", f)
        kxx = CInt(kxi - kpt * 10 - TextSize1.Width - TextSize2.Height)
        kyy = CInt((kyi + kyf) / 2 + TextSize2.Width / 2)
        Call INC_STR(g, f, str, kxx, kyy, -90)
        'ｘ軸名描画
        str = sxjiku
        TextSize2 = g.MeasureString(str, f)
        kxx = CInt((kxi + kxf) / 2 - TextSize2.Width / 2)
        kyy = CInt(kyf + kpt * 10 + TextSize1.Height)
        Call INC_STR(g, f, str, kxx, kyy, 0)
        '丸印描画
        'データプロット
        dc = kpt * 7
        For i = 0 To ndata
            xx = datax(i)
            yy = datay(i)
            kxx = kxi + CInt((xx - xmin) * (kxf - kxi) / (xmax - xmin)) - CInt(dc / 2)
            kyy = kyf - CInt((yy - ymin) * (kyf - kyi) / (ymax - ymin)) - CInt(dc / 2)
            g.DrawEllipse(New Pen(Color.Red, kpt * 1), kxx, kyy, dc, dc)
        Next i
        '回帰式(y=aa*x+bb,rr)
        Call KAIKI(ndata, datax, datay, aa, bb, rr)
        '回帰曲線描画
        xx = xmin : yy = aa * xx + bb
        If yy >= ymax Then
            yy = ymax : xx = (yy - bb) / aa
        End If
        kxx1 = kxi + CInt((xx - xmin) * (kxf - kxi) / (xmax - xmin))
        kyy1 = kyf - CInt((yy - ymin) * (kyf - kyi) / (ymax - ymin))
        xx = xmax : yy = aa * xx + bb
        If yy < ymin Then
            yy = ymin : xx = (yy - bb) / aa
        End If
        kxx2 = kxi + CInt((xx - xmin) * (kxf - kxi) / (xmax - xmin))
        kyy2 = kyf - CInt((yy - ymin) * (kyf - kyi) / (ymax - ymin))
        g.DrawLine(New Pen(Color.Blue, 1), kxx1, kyy1, kxx2, kyy2)
        '回帰係数・相関係数表示
        str = "00000000000"
        TextSize1 = g.MeasureString(str, f)
        g.FillRectangle(Brushes.White, kxf - TextSize1.Width - kpt * 10, kyi + kpt * 10, TextSize1.Width, TextSize1.Height * 4 + kpt * 10)
        g.DrawRectangle(Pens.Black, kxf - TextSize1.Width - kpt * 10, kyi + kpt * 10, TextSize1.Width, TextSize1.Height * 4 + kpt * 10)
        str = "y=a*x+b"
        g.DrawString(str, f, Brushes.Black, kxf - TextSize1.Width - kpt * 5, kyi + kpt * 10 + kpt * 5)
        str = "a=" & aa.ToString("0.00000")
        g.DrawString(str, f, Brushes.Black, kxf - TextSize1.Width - kpt * 5, kyi + kpt * 10 + kpt * 5 + TextSize1.Height * 1)
        str = "b=" & bb.ToString("0.00000")
        g.DrawString(str, f, Brushes.Black, kxf - TextSize1.Width - kpt * 5, kyi + kpt * 10 + kpt * 5 + TextSize1.Height * 2)
        str = "r=" & rr.ToString("0.00000")
        g.DrawString(str, f, Brushes.Black, kxf - TextSize1.Width - kpt * 5, kyi + kpt * 10 + kpt * 5 + TextSize1.Height * 3)
        f.Dispose()

    End Sub


    Private Sub INC_STR(ByVal g As Graphics, ByVal f As System.Drawing.Font, ByVal str As String, _
                    ByVal kxx As Integer, ByVal kyy As Integer, ByVal ang As Single)
        '軸名描画
        g.ScaleTransform(1.0, 1.0)                '横・縦の表示比率を設定
        g.TranslateTransform(kxx, kyy)            '表示位置の設定（表示位置を原点とする座標移動）
        g.RotateTransform(ang)                    '表示角度を指定
        g.DrawString(str, f, Brushes.Black, 0, 0) '描画実行
        g.ResetTransform()                        '単位行列にリセット
    End Sub


    Private Sub DRN_XJIKU(ByVal g As Graphics, ByVal f As System.Drawing.Font, ByVal LPen As System.Drawing.Pen, _
                      ByVal kpt As Integer, _
                      ByVal xmin As Double, ByVal xmax As Double, ByVal dx As Double, _
                      ByVal kxi As Integer, ByVal kxf As Integer, _
                      ByVal kyi As Integer, ByVal kyf As Integer)
        Dim i As Integer
        Dim ix As Integer
        Dim xx As Double
        Dim wv As Double
        Dim str As String
        Dim kxx As Integer
        Dim TextSize1 As New System.Drawing.SizeF

        '普通ｘ軸描画
        ix = CInt((xmax - xmin) / dx)
        For i = 0 To ix
            xx = xmin + CDbl(i) * dx
            kxx = kxi + CInt((xx - xmin) * (kxf - kxi) / (xmax - xmin))
            wv = xmin + CDbl(i) * dx
            str = wv.ToString("0")
            If CInt(dx * 1000.0) Mod 10 <> 0 Then str = wv.ToString("0.000")
            If CInt(dx * 1000.0) Mod 10 = 0 And CInt(dx * 100.0) Mod 10 <> 0 Then str = wv.ToString("0.00")
            If (CInt(dx * 1000.0) Mod 10 = 0 And CInt(dx * 100.0) Mod 10 = 0) And CInt(dx * 10.0) Mod 10 <> 0 Then str = wv.ToString("0.0")
            TextSize1 = g.MeasureString(str, f)
            g.DrawLine(LPen, kxx, kyi, kxx, kyf)
            g.DrawString(str, f, Brushes.Black, kxx - TextSize1.Width / 2, kyf + kpt * 3)
        Next i
    End Sub


    Private Sub DRL_XJIKU(ByVal g As Graphics, ByVal f As System.Drawing.Font, ByVal LPen As System.Drawing.Pen, _
                          ByVal kpt As Integer, _
                          ByVal xmin As Double, ByVal xmax As Double, _
                          ByVal kxi As Integer, ByVal kxf As Integer, _
                          ByVal kyi As Integer, ByVal kyf As Integer)
        Dim i As Integer
        Dim j As Integer
        Dim xx As Double
        Dim str As String
        Dim kxx As Integer
        Dim TextSize1 As New System.Drawing.SizeF

        '対数ｘ軸描画
        For i = CInt(xmin) To CInt(xmax) - 1
            For j = 1 To 9
                xx = System.Math.Log(CDbl(j) * System.Math.Pow(10.0, CDbl(i))) / System.Math.Log(10.0)
                kxx = kxi + CInt((xx - xmin) * (kxf - kxi) / (xmax - xmin))
                g.DrawLine(LPen, kxx, kyi, kxx, kyf)
            Next j
        Next i
        For i = CInt(xmin) To CInt(xmax)
            xx = CDbl(i)
            kxx = kxi + CInt((xx - xmin) * (kxf - kxi) / (xmax - xmin))
            str = CStr(System.Math.Pow(10.0, CDbl(i)))
            TextSize1 = g.MeasureString(str, f)
            g.DrawString(str, f, Brushes.Black, kxx - TextSize1.Width / 2, kyf + kpt * 3)
        Next i
    End Sub


    Private Sub DRH_YJIKU(ByVal g As Graphics, ByVal f As System.Drawing.Font, ByVal LPen As System.Drawing.Pen, _
                          ByVal kpt As Integer, _
                          ByVal ymin As Double, ByVal ymax As Double, _
                          ByVal kxi As Integer, ByVal kxf As Integer, _
                          ByVal kyi As Integer, ByVal kyf As Integer)
        Dim i As Integer
        Dim yy As Double
        Dim str As String = ""
        Dim kyy As Integer
        Dim TextSize1 As New System.Drawing.SizeF
        Dim work As Double

        For i = 1 To 11
            If i = 1 Then str = "99.9"
            If i = 2 Then str = "99"
            If i = 3 Then str = "95"
            If i = 4 Then str = "90"
            If i = 5 Then str = "80"
            If i = 6 Then str = "50"
            If i = 7 Then str = "20"
            If i = 8 Then str = "10"
            If i = 9 Then str = "5"
            If i = 10 Then str = "1"
            If i = 11 Then str = "0.1"
            work = CDbl(CDbl(CDbl(str) / 100.0)) : yy = TODA(work)
            kyy = kyf - CInt((yy - ymin) * (kyf - kyi) / (ymax - ymin))
            TextSize1 = g.MeasureString(str, f)
            g.DrawLine(LPen, kxi, kyy, kxf, kyy)
            g.DrawString(str, f, Brushes.Black, kxi - TextSize1.Width - 2, kyy - TextSize1.Height / 2)
        Next i
    End Sub


    Private Sub KAIKI(ByVal ndata As Integer, ByRef datax() As Double, ByRef datay() As Double, _
                      ByRef aa As Double, ByRef bb As Double, ByRef rr As Double)

        Dim i As Integer
        Dim x1 As Double : Dim y1 As Double : Dim x2 As Double : Dim xy As Double
        Dim xm As Double : Dim ym As Double
        Dim c1 As Double : Dim c2 As Double : Dim c3 As Double

        '回帰式：y=aa*x+bb
        x1 = 0.0 : y1 = 0.0 : x2 = 0.0 : xy = 0.0
        For i = 0 To ndata
            x1 = x1 + datax(i)
            y1 = y1 + datay(i)
            x2 = x2 + datax(i) * datax(i)
            xy = xy + datax(i) * datay(i)
        Next i
        xm = x1 / CDbl(ndata + 1) : ym = y1 / CDbl(ndata + 1)
        aa = (CDbl(ndata + 1) * xy - x1 * y1) / (CDbl(ndata + 1) * x2 - x1 * x1)
        bb = (x2 * y1 - x1 * xy) / (CDbl(ndata + 1) * x2 - x1 * x1)
        c1 = 0.0 : c2 = 0.0 : c3 = 0.0
        For i = 0 To ndata
            c1 = c1 + (datax(i) - xm) * (datay(i) - ym)
            c2 = c2 + (datax(i) - xm) * (datax(i) - xm)
            c3 = c3 + (datay(i) - ym) * (datay(i) - ym)
        Next i
        rr = c1 / Math.Sqrt(c2 * c3)

    End Sub

    Private Function TODA(ByVal rx As Double) As Double

        Dim ry As Double
        Dim ay As Double
        Dim uay As Double
        Dim b0 As Double
        Dim b1 As Double
        Dim b2 As Double
        Dim b3 As Double
        Dim b4 As Double
        Dim b5 As Double
        Dim b6 As Double
        Dim b7 As Double
        Dim b8 As Double
        Dim b9 As Double
        Dim b10 As Double

        If rx = 0.5 Then
            uay = 0
        Else
            If rx < 0.5 Then
                ry = rx
            Else
                ry = 1 - rx
            End If
            ay = -Math.Log(4.0 * ry * (1.0 - ry))
            b0 = 1.570796288
            b1 = 0.03706987906
            b2 = -0.0008364353589
            b3 = -0.0002250947176
            b4 = 0.000006841218299
            b5 = 0.000005824238515
            b6 = -0.00000104527497
            b7 = 0.00000008360937017
            b8 = -0.000000003231081277
            b9 = 0.00000000003657763036
            b10 = 0.0000000000006936233982
            uay = Math.Sqrt(ay * (b0 + b1 * ay + b2 * ay ^ 2.0 + b3 * ay ^ 3.0 + b4 * ay ^ 4.0 + b5 * ay ^ 5.0 + b6 * ay ^ 6.0 + b7 * ay ^ 7.0 + b8 * ay ^ 8.0 + b9 * ay ^ 9.0 + b10 * ay ^ 10.0))
            If rx < 0.5 Then uay = -uay
        End If
        TODA = uay

    End Function

End Class