File indexing completed on 2024-05-31 04:19:24
0001 #include "DataFormats/Math/interface/FastMath.h"
0002
0003 #include <iostream>
0004
0005 #include <typeinfo>
0006
0007 #include "FWCore/Utilities/interface/HRRealTime.h"
0008
0009 namespace {
0010
0011 template <typename T>
0012 std::pair<T, T> stdatan2r(T x, T y) {
0013 return std::pair<T, T>(std::atan2(x, y), std::sqrt(x * x + y * y));
0014 }
0015
0016 template <typename T>
0017 struct Stat {
0018 std::string name;
0019 size_t n;
0020 size_t npos;
0021 T bias;
0022 T ave;
0023 T rms;
0024 T amax;
0025 Stat(std::string in) : name(in), n(0), npos(0), bias(0), ave(0), rms(0), amax(0) {}
0026 void operator()(T x, T ref);
0027
0028 ~Stat() {
0029 std::cout << name << " " << n << " " << npos << " " << bias / n << " " << ave / n << " "
0030 << (n * rms - bias * bias) / (n * (n - 1)) << " " << amax << std::endl;
0031 }
0032 };
0033
0034 template <typename T>
0035 void Stat<T>::operator()(T x, T ref) {
0036 n++;
0037 if (x > ref)
0038 npos++;
0039 T d = (x - ref) / std::abs(ref);
0040 bias += d;
0041 ave += std::abs(d);
0042 rms += d * d;
0043 amax = std::max(amax, std::abs(d));
0044 }
0045
0046 double dummy;
0047 template <typename T>
0048 inline T eta(T x, T y, T z) {
0049 x = z / std::sqrt(x * x + y * y);
0050 return std::log(x + std::sqrt(x * x + T(1)));
0051 }
0052
0053 template <typename T>
0054 void sampleSquare() {
0055 edm::HRTimeType tf = 0;
0056 edm::HRTimeType ts = 0;
0057 Stat<T> stata("atan2");
0058 Stat<T> statr("r");
0059 T fac[8] = {-8, -5., -2., -1., 1., 2., 5., 8.};
0060 for (int k = 0; k < 100; k++)
0061 for (T x = 1e-15; x < 1.1e+15; x *= 10)
0062 for (T y = 1e-15; y < 1.1e+15; y *= 10)
0063 for (int i = 0; i != 8; ++i)
0064 for (int j = 0; j != 8; ++j) {
0065 T xx = x * fac[i];
0066 T yy = y * fac[j];
0067 edm::HRTimeType sf = edm::hrRealTime();
0068 std::pair<T, T> res = fastmath::atan2r(xx, yy);
0069 tf += (edm::hrRealTime() - sf);
0070 for (int l = 0; l < i + j; ++l)
0071 dummy += yy;
0072 edm::HRTimeType ss = edm::hrRealTime();
0073 std::pair<T, T> ref = stdatan2r(xx, yy);
0074 ts += (edm::hrRealTime() - ss);
0075 stata(res.first, ref.first);
0076 statr(res.second, ref.second);
0077 }
0078 std::cout << typeid(T).name() << " times " << tf << " " << ts << std::endl;
0079 }
0080
0081 }
0082
0083 int main() {
0084
0085
0086 {
0087 std::pair<double, double> res = fastmath::atan2r(-3., 4.);
0088 std::cout << res.first << " " << res.second << std::endl;
0089 std::cout << atan2(-3., 4.) << std::endl;
0090 }
0091 {
0092 std::pair<double, double> res = fastmath::etaphi(-3., 4., 5.);
0093 std::cout << res.first << " " << res.second << std::endl;
0094 std::cout << eta(-3., 4., 5.) << " " << std::atan2(4., -3.) << std::endl;
0095 }
0096
0097 {
0098 std::pair<float, float> res = fastmath::atan2r(3.f, -4.f);
0099 std::cout << res.first << " " << res.second << std::endl;
0100 std::cout << atan2f(3.f, -4.f) << std::endl;
0101 }
0102 {
0103 std::pair<double, double> res = fastmath::etaphi(3.f, -4.f, -5.f);
0104 std::cout << res.first << " " << res.second << std::endl;
0105 std::cout << eta(3.f, -4.f, -5.f) << " " << std::atan2(-4.f, 3.f) << std::endl;
0106 }
0107
0108 sampleSquare<float>();
0109 sampleSquare<double>();
0110 return 0;
0111 }