]>
Commit | Line | Data |
---|---|---|
7c673cae FG |
1 | // Copyright John Maddock 2008. |
2 | // Use, modification and distribution are subject to the | |
3 | // Boost Software License, Version 1.0. (See accompanying file | |
4 | // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt) | |
5 | ||
6 | #ifndef BOOST_MATH_DISTIBUTIONS_DETAIL_GENERIC_QUANTILE_HPP | |
7 | #define BOOST_MATH_DISTIBUTIONS_DETAIL_GENERIC_QUANTILE_HPP | |
8 | ||
9 | namespace boost{ namespace math{ namespace detail{ | |
10 | ||
11 | template <class Dist> | |
12 | struct generic_quantile_finder | |
13 | { | |
14 | typedef typename Dist::value_type value_type; | |
15 | typedef typename Dist::policy_type policy_type; | |
16 | ||
17 | generic_quantile_finder(const Dist& d, value_type t, bool c) | |
18 | : dist(d), target(t), comp(c) {} | |
19 | ||
20 | value_type operator()(const value_type& x) | |
21 | { | |
22 | return comp ? | |
23 | value_type(target - cdf(complement(dist, x))) | |
24 | : value_type(cdf(dist, x) - target); | |
25 | } | |
26 | ||
27 | private: | |
28 | Dist dist; | |
29 | value_type target; | |
30 | bool comp; | |
31 | }; | |
32 | ||
33 | template <class T, class Policy> | |
34 | inline T check_range_result(const T& x, const Policy& pol, const char* function) | |
35 | { | |
36 | if((x >= 0) && (x < tools::min_value<T>())) | |
37 | return policies::raise_underflow_error<T>(function, 0, pol); | |
38 | if(x <= -tools::max_value<T>()) | |
39 | return -policies::raise_overflow_error<T>(function, 0, pol); | |
40 | if(x >= tools::max_value<T>()) | |
41 | return policies::raise_overflow_error<T>(function, 0, pol); | |
42 | return x; | |
43 | } | |
44 | ||
45 | template <class Dist> | |
46 | typename Dist::value_type generic_quantile(const Dist& dist, const typename Dist::value_type& p, const typename Dist::value_type& guess, bool comp, const char* function) | |
47 | { | |
48 | typedef typename Dist::value_type value_type; | |
49 | typedef typename Dist::policy_type policy_type; | |
50 | typedef typename policies::normalise< | |
51 | policy_type, | |
52 | policies::promote_float<false>, | |
53 | policies::promote_double<false>, | |
54 | policies::discrete_quantile<>, | |
55 | policies::assert_undefined<> >::type forwarding_policy; | |
56 | ||
57 | // | |
58 | // Special cases first: | |
59 | // | |
60 | if(p == 0) | |
61 | { | |
62 | return comp | |
63 | ? check_range_result(range(dist).second, forwarding_policy(), function) | |
64 | : check_range_result(range(dist).first, forwarding_policy(), function); | |
65 | } | |
66 | if(p == 1) | |
67 | { | |
68 | return !comp | |
69 | ? check_range_result(range(dist).second, forwarding_policy(), function) | |
70 | : check_range_result(range(dist).first, forwarding_policy(), function); | |
71 | } | |
72 | ||
73 | generic_quantile_finder<Dist> f(dist, p, comp); | |
74 | tools::eps_tolerance<value_type> tol(policies::digits<value_type, forwarding_policy>() - 3); | |
75 | boost::uintmax_t max_iter = policies::get_max_root_iterations<forwarding_policy>(); | |
76 | std::pair<value_type, value_type> ir = tools::bracket_and_solve_root( | |
77 | f, guess, value_type(2), true, tol, max_iter, forwarding_policy()); | |
78 | value_type result = ir.first + (ir.second - ir.first) / 2; | |
79 | if(max_iter >= policies::get_max_root_iterations<forwarding_policy>()) | |
80 | { | |
81 | return policies::raise_evaluation_error<value_type>(function, "Unable to locate solution in a reasonable time:" | |
82 | " either there is no answer to quantile" | |
83 | " or the answer is infinite. Current best guess is %1%", result, forwarding_policy()); | |
84 | } | |
85 | return result; | |
86 | } | |
87 | ||
88 | }}} // namespaces | |
89 | ||
90 | #endif // BOOST_MATH_DISTIBUTIONS_DETAIL_GENERIC_QUANTILE_HPP | |
91 |