/[gsl]/gsl/cdf/tdistinv.c
ViewVC logotype

Diff of /gsl/cdf/tdistinv.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.1.2.1 by bjg, Sun Jun 29 13:33:36 2003 UTC revision 1.1.2.2 by bjg, Mon Jul 14 14:33:22 2003 UTC
# Line 17  Line 17 
17   * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307, USA.   * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307, USA.
18   */   */
19    
20    #include <config.h>
21    #include <math.h>
22    #include <gsl/gsl_cdf.h>
23    #include <gsl/gsl_math.h>
24    #include <gsl/gsl_randist.h>
25    #include <gsl/gsl_sf_gamma.h>
26    
27    #include <stdio.h>
28    
29    double
30    gsl_cdf_tdist_Pinv (const double p, const double nu)
31    {
32      double x;
33      
34      if (p == 1.0)
35        {
36          return GSL_POSINF;
37        }
38      else if (p == 0.0)
39        {
40          return GSL_NEGINF;
41        }
42    
43      if (p < 0.05)
44        {
45          
46        }
47      else if (p > 0.95)
48        {
49          
50        }
51      else
52        {
53          double xg = gsl_cdf_ugaussian_Pinv (p);
54          double x0 = xg + xg*(xg*xg+1.0)/(4.0*nu) ;
55          x = x0;
56        }
57    
58      {
59        double dp, phi;
60    
61      start:
62        dp = p - gsl_cdf_tdist_P (x, nu);
63        phi = gsl_ran_tdist_pdf (x, nu);
64    
65        if (dp == 0.0)
66          goto end;
67    
68        {
69          double step = dp / phi;
70    
71          printf("x=%.18e step=%.18e\n", x, step);
72          x += step;
73          
74          if (fabs (step) > 1e-10 * fabs(x))
75            goto start;
76        }
77      }
78    
79    end:
80    
81      return x;
82    }
83    
84  #if 0  #if 0
85  /*  /*
86   * Invert the T distribution. Uses a method as shown in   * Invert the T distribution. Uses a method as shown in
# Line 79  tdist_pdf (const double x, const double Line 143  tdist_pdf (const double x, const double
143    return p;    return p;
144  }  }
145    
 double  
 gsl_cdf_ut_Pinv (double prob, double nu)  
 {  
146    double initial_result;    double initial_result;
147    double result;    double result;
148    double d;    double d;

Legend:
Removed from v.1.1.2.1  
changed lines
  Added in v.1.1.2.2

savannah-hackers-public@gnu.org
ViewVC Help
Powered by ViewVC 1.1.26