*! version 1.0.0  14may1995
program define df
        version 4.0

        preserve
        clear

        quietly {
                set obs 201
                gen x = 3*((_n-1)/200)
                gen y = .}

        global df1=10
        global df2=10

        wdctl button "Reset"   15  5 25 12 HJN_B0
        wdctl button "df1 +"   15 20 25 12 HJN_B1
        wdctl button "df1 -"   15 35 25 12 HJN_B2
        wdctl button "df2 +"   15 50 25 12 HJN_B3
        wdctl button "df2 -"   15 65 25 12 HJN_B4
        wdctl button "Help"    15 80 25 12 HJN_B5
        wdctl button "Close"   15 95 25 12 HJN_B6

        global D_sm1 "Inc"
        global dfinc=5
        wdctl static D_sm1  3 110 10 10 
        wdctl edit   dfinc 15 110 25 10         
        
        global D_sm2 "df_1:"
        global D_sm3 "df_2:"
        wdctl static D_sm2    15 125 15 12
        wdctl static df1      35 125 15 12
        wdctl static D_sm3    15 140 15 12
        wdctl static df2      35 140 15 12



        global HJN_B0 "dfdr 0"
        global HJN_B1 "dfdr 1"
        global HJN_B2 "dfdr -1"
        global HJN_B3 "dfdr 2"
        global HJN_B4 "dfdr -2"
        global HJN_B5 "whelp df"
        global HJN_B6 "exit 1234"

        gph open

        dfdr 0

        noi wdlg "F Curves"  5 5 60 165

        gph close

        restore

end

program define dfdr
        version 4.0

        local arg=`1'

        if `arg'!=0 {
        graph y x, xlabel ylabel xscale(0,3) yscale(0,3) c(l) s(i) l1(" ") l2(" ") b1(" ") b2(" ") pen(3)
        }
                
        if `arg'==0 { 
                global df1=10
                global df2=10
                gph close 
                gph open}


        if `arg'==1 {
                global df1 = min(200,$df1+$dfinc)}

        if `arg'==-1 {
                global df1 = max(10,$df1-$dfinc)}

        if `arg'==2 {
                global df2 = min(200,$df2+$dfinc)}

        if `arg'==-2 {
                global df2 = max(10,$df2-$dfinc)}



        scalar v1=$df1-1
        scalar v2=$df2-1
        scalar v1half=v1/2
        scalar v2half=v2/2
        scalar vsum=v1half+v2half
        scalar vrat=v1half/v2half
        scalar c1=v1half*log(vrat)+lngamma(vsum)-lngamma(v1half)-lngamma(v2half)
        replace y=exp(c1+(v1half-1)*log(x)-vsum*log(1+vrat*x))

        graph y x, xlabel ylabel xscale(0,3) yscale(0,3) c(l) s(i) l1(" ") l2(" ") b1(" ") b2(" ") pen(2)

*        gph text  4000 4500  0 1 F($df1,$df2)
        
end

exit