Mercurial > octave-libtiff
view liboctave/numeric/lo-mappers.cc @ 23219:3ac9f9ecfae5 stable
maint: Update copyright dates.
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Wed, 22 Feb 2017 12:39:29 -0500 |
parents | e9a0469dedd9 |
children | 092078913d54 |
line wrap: on
line source
/* Copyright (C) 1996-2017 John W. Eaton Copyright (C) 2010 VZLU Prague This file is part of Octave. Octave is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version. Octave is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with Octave; see the file COPYING. If not, see <http://www.gnu.org/licenses/>. */ #if defined (HAVE_CONFIG_H) # include "config.h" #endif #include <cfloat> #include "lo-error.h" #include "lo-ieee.h" #include "lo-mappers.h" #include "lo-math.h" #include "lo-specfun.h" #include "lo-utils.h" #include "math-wrappers.h" #include "oct-cmplx.h" #include "f77-fcn.h" // FIXME: We used to have this situation: // // Functions that forward to gnulib belong here so we can keep // gnulib:: out of lo-mappers.h. // // but now we just use std:: and explicit wrappers in C++ code so maybe // some of the forwarding functions can be defined inline here. namespace octave { namespace math { bool is_NA (double x) { return lo_ieee_is_NA (x); } bool is_NA (const Complex& x) { return (is_NA (std::real (x)) || is_NA (std::imag (x))); } bool is_NA (float x) { return lo_ieee_is_NA (x); } bool is_NA (const FloatComplex& x) { return (is_NA (std::real (x)) || is_NA (std::imag (x))); } bool is_NaN_or_NA (const Complex& x) { return (isnan (std::real (x)) || isnan (std::imag (x))); } bool is_NaN_or_NA (const FloatComplex& x) { return (isnan (std::real (x)) || isnan (std::imag (x))); } Complex acos (const Complex& x) { #if defined (HAVE_COMPLEX_STD_ACOS) Complex y = std::acos (x); if (std::imag (x) == 0.0 && std::real (x) > 1.0) return std::conj (y); else return y; #else static Complex i (0, 1); Complex tmp; if (std::imag (x) == 0.0) { // If the imaginary part of X is 0, then avoid generating an // imaginary part of -0 for the expression 1-x*x. // This chooses the same phase of the branch cut as Matlab. double xr = std::real (x); tmp = Complex (1.0 - xr*xr); } else tmp = 1.0 - x*x; return -i * log (x + i * sqrt (tmp)); #endif } FloatComplex acos (const FloatComplex& x) { #if defined (HAVE_COMPLEX_STD_ACOS) FloatComplex y = std::acos (x); if (std::imag (x) == 0.0f && std::real (x) > 1.0f) return std::conj (y); else return y; #else static FloatComplex i (0, 1); FloatComplex tmp; if (std::imag (x) == 0.0f) { // If the imaginary part of X is 0, then avoid generating an // imaginary part of -0 for the expression 1-x*x. // This chooses the same phase of the branch cut as Matlab. float xr = std::real (x); tmp = FloatComplex (1.0f - xr*xr); } else tmp = 1.0f - x*x; return -i * log (x + i * sqrt (tmp)); #endif } Complex asin (const Complex& x) { #if defined (HAVE_COMPLEX_STD_ASIN) Complex y = std::asin (x); if (std::imag (x) == 0.0 && std::real (x) > 1.0) return std::conj (y); else return y; #else static Complex i (0, 1); Complex tmp; if (std::imag (x) == 0.0) { // If the imaginary part of X is 0, then avoid generating an // imaginary part of -0 for the expression 1-x*x. // This chooses the same phase of the branch cut as Matlab. double xr = std::real (x); tmp = Complex (1.0 - xr*xr); } else tmp = 1.0 - x*x; return -i * log (i*x + sqrt (tmp)); #endif } FloatComplex asin (const FloatComplex& x) { #if defined (HAVE_COMPLEX_STD_ASIN) FloatComplex y = std::asin (x); if (std::imag (x) == 0.0f && std::real (x) > 1.0f) return std::conj (y); else return y; #else static FloatComplex i (0, 1); FloatComplex tmp; if (std::imag (x) == 0.0f) { // If the imaginary part of X is 0, then avoid generating an // imaginary part of -0 for the expression 1-x*x. // This chooses the same phase of the branch cut as Matlab. float xr = std::real (x); tmp = FloatComplex (1.0f - xr*xr); } else tmp = 1.0f - x*x; return -i * log (i*x + sqrt (tmp)); #endif } Complex atan (const Complex& x) { #if defined (HAVE_COMPLEX_STD_ATAN) return std::atan (x); #else static Complex i (0, 1); return i * log ((i + x) / (i - x)) / 2.0; #endif } FloatComplex atan (const FloatComplex& x) { #if defined (HAVE_COMPLEX_STD_ATAN) return std::atan (x); #else static FloatComplex i (0, 1); return i * log ((i + x) / (i - x)) / 2.0f; #endif } double log2 (double x) { return std::log2 (x); } float log2 (float x) { return std::log2 (x); } Complex log2 (const Complex& x) { #if defined (M_LN2) static double ln2 = M_LN2; #else static double ln2 = std::log (2.0); #endif return std::log (x) / ln2; } FloatComplex log2 (const FloatComplex& x) { #if defined (M_LN2) static float ln2 = M_LN2; #else static float ln2 = log (2.0f); #endif return std::log (x) / ln2; } double log2 (double x, int& exp) { return frexp (x, &exp); } float log2 (float x, int& exp) { return frexp (x, &exp); } Complex log2 (const Complex& x, int& exp) { double ax = std::abs (x); double lax = log2 (ax, exp); return (ax != lax) ? (x / ax) * lax : x; } FloatComplex log2 (const FloatComplex& x, int& exp) { float ax = std::abs (x); float lax = log2 (ax, exp); return (ax != lax) ? (x / ax) * lax : x; } double exp2 (double x) { #if defined (HAVE_EXP2) return ::exp2 (x); #else # if defined (M_LN2) static double ln2 = M_LN2; # else static double ln2 = std::log (2.0); # endif return exp (x * ln2); #endif } float exp2 (float x) { #if defined (HAVE_EXP2F) return exp2f (x); #elif defined (HAVE_EXP2) return ::exp2 (x); #else # if defined (M_LN2) static float ln2 = M_LN2; # else static float ln2 = log2 (2.0f); # endif return exp (x * ln2); #endif } double copysign (double x, double y) { return std::copysign (x, y); } float copysign (float x, float y) { return std::copysign (x, y); } double signbit (double x) { return std::signbit (x); } float signbit (float x) { return std::signbit (x); } bool negative_sign (double x) { return __lo_ieee_signbit (x); } bool negative_sign (float x) { return __lo_ieee_float_signbit (x); } double trunc (double x) { return std::trunc (x); } float trunc (float x) { return std::trunc (x); } double floor (double x) { return std::floor (x); } float floor (float x) { return std::floor (x); } double round (double x) { return std::round (x); } float round (float x) { return std::round (x); } double frexp (double x, int *expptr) { return octave_frexp_wrapper (x, expptr); } float frexp (float x, int *expptr) { return octave_frexpf_wrapper (x, expptr); } bool isnan (double x) { #if defined (HAVE_CMATH_ISNAN) return std::isnan (x); #else return lo_ieee_isnan (x); #endif } bool isnan (float x) { #if defined (HAVE_CMATH_ISNANF) return std::isnan (x); #else return lo_ieee_isnan (x); #endif } bool finite (double x) { #if defined (HAVE_CMATH_ISFINITE) return std::isfinite (x); #else return lo_ieee_finite (x); #endif } bool finite (float x) { #if defined (HAVE_CMATH_ISFINITEF) return std::isfinite (x); #else return lo_ieee_finite (x); #endif } bool isinf (double x) { #if defined (HAVE_CMATH_ISINF) return std::isinf (x); #else return lo_ieee_isinf (x); #endif } bool isinf (float x) { #if defined (HAVE_CMATH_ISINFF) return std::isinf (x); #else return lo_ieee_isinf (x); #endif } // Sometimes you need a large integer, but not always. octave_idx_type nint_big (double x) { if (x > std::numeric_limits<octave_idx_type>::max ()) return std::numeric_limits<octave_idx_type>::max (); else if (x < std::numeric_limits<octave_idx_type>::min ()) return std::numeric_limits<octave_idx_type>::min (); else return static_cast<octave_idx_type> ((x > 0.0) ? (x + 0.5) : (x - 0.5)); } octave_idx_type nint_big (float x) { if (x > std::numeric_limits<octave_idx_type>::max ()) return std::numeric_limits<octave_idx_type>::max (); else if (x < std::numeric_limits<octave_idx_type>::min ()) return std::numeric_limits<octave_idx_type>::min (); else return static_cast<octave_idx_type> ((x > 0.0f) ? (x + 0.5f) : (x - 0.5f)); } int nint (double x) { if (x > std::numeric_limits<int>::max ()) return std::numeric_limits<int>::max (); else if (x < std::numeric_limits<int>::min ()) return std::numeric_limits<int>::min (); else return static_cast<int> ((x > 0.0) ? (x + 0.5) : (x - 0.5)); } int nint (float x) { if (x > std::numeric_limits<int>::max ()) return std::numeric_limits<int>::max (); else if (x < std::numeric_limits<int>::min ()) return std::numeric_limits<int>::min (); else return static_cast<int> ((x > 0.0f) ? (x + 0.5f) : (x - 0.5f)); } Complex rc_acos (double x) { return fabs (x) > 1.0 ? acos (Complex (x)) : Complex (::acos (x)); } FloatComplex rc_acos (float x) { return fabsf (x) > 1.0f ? acos (FloatComplex (x)) : FloatComplex (::acosf (x)); } Complex rc_acosh (double x) { return x < 1.0 ? acosh (Complex (x)) : Complex (acosh (x)); } FloatComplex rc_acosh (float x) { return x < 1.0f ? acosh (FloatComplex (x)) : FloatComplex (acosh (x)); } Complex rc_asin (double x) { return fabs (x) > 1.0 ? asin (Complex (x)) : Complex (::asin (x)); } FloatComplex rc_asin (float x) { return fabsf (x) > 1.0f ? asin (FloatComplex (x)) : FloatComplex (::asinf (x)); } Complex rc_atanh (double x) { return fabs (x) > 1.0 ? atanh (Complex (x)) : Complex (atanh (x)); } FloatComplex rc_atanh (float x) { return fabsf (x) > 1.0f ? atanh (FloatComplex (x)) : FloatComplex (atanh (x)); } Complex rc_log (double x) { const double pi = 3.14159265358979323846; return x < 0.0 ? Complex (std::log (-x), pi) : Complex (std::log (x)); } FloatComplex rc_log (float x) { const float pi = 3.14159265358979323846f; return x < 0.0f ? FloatComplex (std::log (-x), pi) : FloatComplex (std::log (x)); } Complex rc_log2 (double x) { const double pil2 = 4.53236014182719380962; // = pi / log(2) return x < 0.0 ? Complex (log2 (-x), pil2) : Complex (log2 (x)); } FloatComplex rc_log2 (float x) { const float pil2 = 4.53236014182719380962f; // = pi / log(2) return x < 0.0f ? FloatComplex (log2 (-x), pil2) : FloatComplex (log2 (x)); } Complex rc_log10 (double x) { const double pil10 = 1.36437635384184134748; // = pi / log(10) return x < 0.0 ? Complex (log10 (-x), pil10) : Complex (log10 (x)); } FloatComplex rc_log10 (float x) { const float pil10 = 1.36437635384184134748f; // = pi / log(10) return x < 0.0f ? FloatComplex (log10 (-x), pil10) : FloatComplex (log10f (x)); } Complex rc_sqrt (double x) { return x < 0.0 ? Complex (0.0, sqrt (-x)) : Complex (sqrt (x)); } FloatComplex rc_sqrt (float x) { return x < 0.0f ? FloatComplex (0.0f, sqrtf (-x)) : FloatComplex (sqrtf (x)); } } }