1 /*
2 * Mathlib : A C Library of Special Functions
3 * Copyright (C) 1998 Ross Ihaka
4 * Copyright (C) 2000 The R Core Team
5 * Copyright (C) 2003 The R Foundation
6 *
7 * This program is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, a copy is available at
19 * http://www.r-project.org/Licenses/
20 *
21 * SYNOPSIS
22 *
23 * double dnorm4(double x, double mu, double sigma, int give_log)
24 * {dnorm (..) is synonymous and preferred inside R}
25 *
26 * DESCRIPTION
27 *
28 * Compute the density of the normal distribution.
29 */
30
31 #include "nmath.h"
32 #include "dpq.h"
33
dnorm4(double x,double mu,double sigma,int give_log)34 double dnorm4(double x, double mu, double sigma, int give_log)
35 {
36 #ifdef IEEE_754
37 if (ISNAN(x) || ISNAN(mu) || ISNAN(sigma))
38 return x + mu + sigma;
39 #endif
40 if(!R_FINITE(sigma)) return R_D__0;
41 if(!R_FINITE(x) && mu == x) return ML_NAN;/* x-mu is NaN */
42 if (sigma <= 0) {
43 if (sigma < 0) ML_ERR_return_NAN;
44 /* sigma == 0 */
45 return (x == mu) ? ML_POSINF : R_D__0;
46 }
47 x = (x - mu) / sigma;
48
49 if(!R_FINITE(x)) return R_D__0;
50 return (give_log ?
51 -(M_LN_SQRT_2PI + 0.5 * x * x + log(sigma)) :
52 M_1_SQRT_2PI * exp(-0.5 * x * x) / sigma);
53 /* M_1_SQRT_2PI = 1 / sqrt(2 * pi) */
54 }
55