GeometryVector/test/BasicVector_t.cpp

0001 #include "DataFormats/GeometryVector/interface/Basic3DVector.h"
0002 #include "DataFormats/GeometryVector/interface/Basic2DVector.h"
0003
0004 #include <vector>
0005
0006 #include <iostream>
0007
0008 // this is a test,
0009 // using namespace mathSSE;
0010
0011 void addScaleddiff(Basic3DVectorF& res, float s, Basic3DVectorF const& a, Basic3DVectorF const& b) {
0012   res += s * (a - b);
0013 }
0014
0015 void addScaleddiff2(Basic3DVectorF& res, float s, Basic3DVectorF const& a, Basic3DVectorF const& b) {
0016   res = res + s * (a - b);
0017 }
0018
0019 void addScaleddiff(Basic3DVectorD& res, float s, Basic3DVectorD const& a, Basic3DVectorD const& b) {
0020   res += s * (a - b);
0021 }
0022
0023 void multiSum(Basic3DVectorF& res, float s, Basic3DVectorF const& a, Basic3DVectorF const& b) {
0024   res = s * (a - b) + s * (a + b);
0025 }
0026
0027 void multiSum(Basic3DVectorD& res, float s, Basic3DVectorD const& a, Basic3DVectorD const& b) {
0028   res = s * (a - b) + s * (a + b);
0029 }
0030
0031 void multiSum(Basic3DVectorLD& res, float s, Basic3DVectorLD const& a, Basic3DVectorLD const& b) {
0032   res = s * (a - b) + s * (a + b);
0033 }
0034
0035 template <typename T, typename U>
0036 typename PreciseFloatType<T, U>::Type dotV(Basic3DVector<T> const& a, Basic3DVector<U> const& b) {
0037   return a * b;
0038 }
0039
0040 template <typename T>
0041 T norm(Basic3DVector<T> const& a) {
0042   return std::sqrt(a * a);
0043 }
0044
0045 template <typename T>
0046 T normV(Basic3DVector<T> const& a) {
0047   return a.mag();
0048 }
0049
0050 template <typename T, typename U>
0051 typename PreciseFloatType<T, U>::Type dotV(Basic2DVector<T> const& a, Basic2DVector<U> const& b) {
0052   return a * b;
0053 }
0054
0055 template <typename T>
0056 T norm(Basic2DVector<T> const& a) {
0057   return std::sqrt(a * a);
0058 }
0059
0060 template <typename T>
0061 T normV(Basic2DVector<T> const& a) {
0062   return a.mag();
0063 }
0064
0065 long aligned(void* p) { return long(p) & 0xf; }
0066
0067 volatile int* vi = 0;
0068
0069 template <typename T>
0070 void verifyAlign() {
0071   int sum = 0;
0072   int nota = 0;
0073   for (int i = 0; i != 100; ++i) {
0074     auto p = new int(3);
0075     sum += (*p);
0076     auto t = new T;
0077     if (aligned(t) != 0)
0078       ++nota;
0079     else
0080       sum += (*t)[0];
0081     delete p;
0082     delete t;
0083     t = new T;
0084     sum += (*t)[0];
0085     if (aligned(t) != 0)
0086       ++nota;
0087     else
0088       sum += (*t)[0];
0089     delete t;
0090     vi = new int[3];
0091     auto vt = new T[3];
0092     if (aligned(vt) != 0)
0093       ++nota;
0094     else
0095       sum += vt[1][1];
0096     delete[] vi;
0097     delete[] vt;
0098   }
0099
0100   std::cout << "sum " << sum << std::endl;
0101   std::cout << "not aligned " << nota << std::endl;
0102 }
0103
0104 int main() {
0105 #if defined(__GNUC__)
0106   std::cout << "gcc " << __GNUC__ << "." << __GNUC_MINOR__ << std::endl;
0107 #endif
0108 #ifdef USE_SSEVECT
0109   std::cout << "sse vector enabled in cmssw" << std::endl;
0110 #endif
0111 #ifdef USE_EXTVECT
0112   std::cout << "extended vector notation enabled in cmssw" << std::endl;
0113 #endif
0114
0115 #ifdef __clang__
0116   struct MxAling {
0117   } __attribute__((__aligned__));
0118   std::cout << "biggest alignment " << alignof(MxAling) << std::endl;
0119 #else
0120   std::cout << "biggest alignment " << __BIGGEST_ALIGNMENT__ << std::endl;
0121 #endif
0122
0123   std::cout << sizeof(Basic2DVectorF) << ' ' << alignof(Basic2DVectorF) << std::endl;
0124   std::cout << sizeof(Basic2DVectorD) << ' ' << alignof(Basic2DVectorD) << std::endl;
0125   std::cout << sizeof(Basic3DVectorF) << ' ' << alignof(Basic3DVectorF) << std::endl;
0126   std::cout << sizeof(Basic3DVectorD) << ' ' << alignof(Basic3DVectorD) << std::endl;
0127   std::cout << sizeof(Basic3DVectorLD) << ' ' << alignof(Basic3DVectorLD) << std::endl;
0128
0129   verifyAlign<Basic3DVectorF>();
0130   // verifyAlign<Basic3DVectorD>();  // crashes for AVX
0131
0132   Basic3DVectorF x(2.0f, 4.0f, 5.0f);
0133   Basic3DVectorF y(-3.0f, 2.0f, -5.0f);
0134   Basic3DVectorD xd(2.0, 4.0, 5.0);
0135   Basic3DVectorD yd = y;
0136
0137   Basic3DVectorLD xld(2.0, 4.0, 5.0);
0138   Basic3DVectorLD yld = y;
0139
0140   Basic2DVectorF x2(2.0f, 4.0f);
0141   Basic2DVectorF y2 = y.xy();
0142   Basic2DVectorD xd2(2.0, 4.0);
0143   Basic2DVectorD yd2 = yd.xy();
0144
0145   {
0146     std::cout << dotV(x, y) << std::endl;
0147     std::cout << normV(x) << std::endl;
0148     std::cout << norm(x) << std::endl;
0149     // std::cout << std::min(x.mathVector(),y.mathVector()) << std::endl;
0150     // std::cout << std::max(x.mathVector(),y.mathVector()) << std::endl;
0151
0152     std::cout << dotV(x, yd) << std::endl;
0153     std::cout << dotV(xd, y) << std::endl;
0154     std::cout << dotV(xd, yd) << std::endl;
0155     std::cout << normV(xd) << std::endl;
0156     std::cout << norm(xd) << std::endl;
0157     std::cout << dotV(xld, yld) << std::endl;
0158     std::cout << normV(xld) << std::endl;
0159     std::cout << norm(xld) << std::endl;
0160
0161     Basic3DVectorF z = x.cross(y);
0162     std::cout << z << std::endl;
0163     std::cout << -z << std::endl;
0164     Basic3DVectorD zd = x.cross(yd);
0165     std::cout << zd << std::endl;
0166     std::cout << -zd << std::endl;
0167     std::cout << xd.cross(y) << std::endl;
0168     std::cout << xd.cross(yd) << std::endl;
0169
0170     Basic3DVectorLD zld = x.cross(yld);
0171     std::cout << zld << std::endl;
0172     std::cout << -zld << std::endl;
0173     std::cout << xld.cross(y) << std::endl;
0174     std::cout << xld.cross(yld) << std::endl;
0175
0176     std::cout << z.eta() << " " << (-z).eta() << std::endl;
0177     std::cout << zd.eta() << " " << (-zd).eta() << std::endl;
0178     std::cout << zld.eta() << " " << (-zld).eta() << std::endl;
0179
0180     auto s = x + xd - 3.1 * z;
0181     std::cout << s << std::endl;
0182     auto s2 = x + xld - 3.1 * zd;
0183     std::cout << s2 << std::endl;
0184   }
0185
0186   {
0187     std::cout << dotV(x2, y2) << std::endl;
0188     std::cout << normV(x2) << std::endl;
0189     std::cout << norm(x2) << std::endl;
0190     // std::cout << std::min(x2.mathVector(),y2.mathVector()) << std::endl;
0191     // std::cout << std::max(x2.mathVector(),y2.mathVector()) << std::endl;
0192
0193     std::cout << dotV(x2, yd2) << std::endl;
0194     std::cout << dotV(xd2, y2) << std::endl;
0195     std::cout << dotV(xd2, yd2) << std::endl;
0196     std::cout << normV(xd2) << std::endl;
0197     std::cout << norm(xd2) << std::endl;
0198
0199     Basic2DVectorF z2(x2);
0200     z2 -= y2;
0201     std::cout << z2 << std::endl;
0202     std::cout << -z2 << std::endl;
0203     Basic2DVectorD zd2 = x2 - yd2;
0204     std::cout << zd2 << std::endl;
0205     std::cout << -zd2 << std::endl;
0206     std::cout << x2.cross(y2) << std::endl;
0207     std::cout << x2.cross(yd2) << std::endl;
0208     std::cout << xd2.cross(y2) << std::endl;
0209     std::cout << xd2.cross(yd2) << std::endl;
0210
0211     auto s2 = x2 + xd2 - 3.1 * z2;
0212     std::cout << s2 << std::endl;
0213   }
0214
0215   {
0216     std::cout << "f" << std::endl;
0217     Basic3DVectorF vx(2.0f, 4.0f, 5.0f);
0218     Basic3DVectorF vy(-3.0f, 2.0f, -5.0f);
0219     vx += vy;
0220     std::cout << vx << std::endl;
0221
0222     Basic3DVectorF vz(1.f, 1.f, 1.f);
0223     addScaleddiff(vz, 0.1f, vx, vy);
0224     std::cout << vz << std::endl;
0225   }
0226
0227   {
0228     std::cout << "d" << std::endl;
0229     Basic3DVectorD vx(2.0, 4.0, 5.0);
0230     Basic3DVectorD vy(-3.0, 2.0, -5.0);
0231     vx += vy;
0232     std::cout << vx << std::endl;
0233
0234     Basic3DVectorD vz(1., 1., 1);
0235     addScaleddiff(vz, 0.1, vx, vy);
0236     std::cout << vz << std::endl;
0237   }
0238
0239   std::cout << "std::vector" << std::endl;
0240   std::vector<Basic3DVectorF> vec1;
0241   vec1.reserve(50);
0242   std::vector<float> vecf(21);
0243   std::vector<Basic3DVectorF> vec2(51);
0244   std::vector<Basic3DVectorF> vec3;
0245   vec3.reserve(23456);
0246 }