EVE/Reve/NLTProjections.cxx

   1 #include "NLTProjections.h"
   2 #include "Reve.h"
   3
   4 using namespace Reve;
   5
   6
   7 //______________________________________________________________________________
   8 // NLTProjection
   9 //
  10 // Base-class for non-linear projection of 3D point.
  11 // Enables to define an external center of distortion and a scale to
  12 // fixate a bounding box of a projected point.
  13
  14 ClassImp(Reve::NLTProjection)
  15
  16 Float_t NLTProjection::fgEps = 0.005f;
  17 //______________________________________________________________________________
  18 NLTProjection::NLTProjection(Vector& center) :
  19   fType(PT_Unknown),
  20   fGeoMode(GM_Unknown),
  21   fName(0),
  22   fCenter(center.x, center.y, center.z),
  23   fDistortion(0.0f),
  24   fFixedRadius(300),
  25   fScale(1.0f)
  26 {
  27   // Constructor.
  28 }
  29
  30 //______________________________________________________________________________
  31 void NLTProjection::ProjectVector(Vector& v)
  32 {
  33   // Project Reve::Vector.
  34
  35   ProjectPoint(v.x, v.y, v.z);
  36 }
  37
  38 //______________________________________________________________________________
  39 void NLTProjection::UpdateLimit()
  40 {
  41   // Update convergence in +inf and -inf.
  42
  43   if ( fDistortion == 0.0f )
  44     return;
  45
  46   Float_t lim =  1.0f/fDistortion + fFixedRadius;
  47   Float_t* c = GetProjectedCenter();
  48   fUpLimit.Set(lim + c[0], lim + c[1], c[2]);
  49   fLowLimit.Set(-lim + c[0], -lim + c[1], c[2]);
  50 }
  51
  52 //______________________________________________________________________________
  53 void NLTProjection::SetDistortion(Float_t d)
  54 {
  55   // Set distortion.
  56
  57   fDistortion=d;
  58   fScale = 1+fFixedRadius*fDistortion;
  59   UpdateLimit();
  60 }
  61
  62 //______________________________________________________________________________
  63 void NLTProjection::SetFixedRadius(Float_t r)
  64 {
  65   // Set fixed radius.
  66
  67   fFixedRadius=r;
  68   fScale = 1 + fFixedRadius*fDistortion;
  69   UpdateLimit();
  70 }
  71
  72 //______________________________________________________________________________
  73 void NLTProjection::SetDirectionalVector(Int_t screenAxis, Vector& vec)
  74 {
  75   // Get vector for axis in a projected space.
  76
  77   for (Int_t i=0; i<3; i++)
  78   {
  79     vec[i] = (i==screenAxis) ? 1. : 0.;
  80   }
  81 }
  82
  83 //______________________________________________________________________________
  84 Float_t NLTProjection::GetValForScreenPos(Int_t i, Float_t sv)
  85 {
  86   // Inverse projection.
  87
  88   static const Exc_t eH("NLTProjection::GetValForScreenPos ");
  89
  90   Float_t xL, xM, xR;
  91   Vector V, DirVec;
  92   SetDirectionalVector(i, DirVec);
  93   if (fDistortion > 0.0f && ((sv > 0 && sv > fUpLimit[i]) || (sv < 0 && sv < fLowLimit[i])))
  94     throw(eH + Form("screen value '%f' out of limit '%f'.", sv, sv > 0 ? fUpLimit[i] : fLowLimit[i]));
  95
  96   Vector zero; ProjectVector(zero);
  97   // search from -/+ infinity according to sign of screen value
  98   if (sv > zero[i])
  99   {
 100     xL = 0; xR = 1000;
 101     while (1)
 102     {
 103       V.Mult(DirVec, xR); ProjectVector(V);
 104       // printf("positive projected %f, value %f,xL, xR ( %f, %f)\n", V[i], sv, xL, xR);
 105       if (V[i] > sv || V[i] == sv) break;
 106       xL = xR; xR *= 2;
 107     }
 108   }
 109   else if (sv < zero[i])
 110   {
 111     xR = 0; xL = -1000;
 112     while (1)
 113     {
 114       V.Mult(DirVec, xL); ProjectVector(V);
 115       // printf("negative projected %f, value %f,xL, xR ( %f, %f)\n", V[i], sv, xL, xR);
 116       if (V[i] < sv || V[i] == sv) break;
 117       xR = xL; xL *= 2;
 118     }
 119   }
 120   else
 121   {
 122     return 0.0f;
 123   }
 124
 125   do
 126   {
 127     xM = 0.5f * (xL + xR);
 128     V.Mult(DirVec, xM);
 129     ProjectVector(V);
 130     if (V[i] > sv)
 131       xR = xM;
 132     else
 133       xL = xM;
 134   } while(TMath::Abs(V[i] - sv) >= fgEps);
 135
 136   return xM;
 137 }
 138
 139 //______________________________________________________________________________
 140 Float_t NLTProjection::GetScreenVal(Int_t i, Float_t x)
 141 {
 142   // Project point on given axis and return projected value.
 143
 144   Vector dv;
 145   SetDirectionalVector(i, dv); dv = dv*x;
 146   ProjectVector(dv);
 147   return dv[i];
 148 }
 149
 150
 151 //______________________________________________________________________________
 152 // NLTRhoZ
 153 //
 154 // Transformation from 3D to 2D. X axis represent Z coordinate. Y axis have value of
 155 // radius with a sign of Y coordinate.
 156
 157 ClassImp(Reve::NLTRhoZ)
 158
 159 //______________________________________________________________________________
 160 void NLTRhoZ::SetCenter(Vector& v)
 161 {
 162   // Set center of distortion (virtual method).
 163
 164   fCenter = v;
 165
 166   Float_t R = TMath::Sqrt(v.x*v.x+v.y*v.y);
 167   fProjectedCenter.x = fCenter.z;
 168   fProjectedCenter.y = TMath::Sign(R, fCenter.y);
 169   fProjectedCenter.z = 0;
 170   UpdateLimit();
 171 }
 172
 173 //______________________________________________________________________________
 174 void NLTRhoZ::ProjectPoint(Float_t& x, Float_t& y, Float_t& z,  PProc_e proc )
 175 {
 176   // Project point.
 177
 178   using namespace TMath;
 179
 180   if(proc == PP_Plane || proc == PP_Full)
 181   {
 182     // project
 183     y = Sign((Float_t)Sqrt(x*x+y*y), y);
 184     x = z;
 185   }
 186   if(proc == PP_Distort || proc == PP_Full)
 187   {
 188     // move to center
 189     x -= fProjectedCenter.x;
 190     y -= fProjectedCenter.y;
 191     // distort
 192     y = (y*fScale) / (1.0f + Abs(y)*fDistortion);
 193     x = (x*fScale) / (1.0f + Abs(x)*fDistortion);
 194     // move back from center
 195     x += fProjectedCenter.x;
 196     y += fProjectedCenter.y;
 197   }
 198   z = 0.0f;
 199 }
 200
 201 //______________________________________________________________________________
 202 void NLTRhoZ::SetDirectionalVector(Int_t screenAxis, Vector& vec)
 203 {
 204   // Get direction in the unprojected space for axis index in the projected space.
 205   // This is virtual method from base-class NLTProjection.
 206
 207   if(screenAxis == 0)
 208     vec.Set(0., 0., 1);
 209   else if (screenAxis == 1)
 210     vec.Set(0., 1., 0);
 211
 212 }
 213 //______________________________________________________________________________
 214 Bool_t NLTRhoZ::AcceptSegment(Vector& v1, Vector& v2, Float_t tolerance)
 215 {
 216   // Check if segment of two projected points is valid.
 217
 218   Float_t a = fProjectedCenter.y;
 219   Bool_t val = kTRUE;
 220   if((v1.y <  a && v2.y > a) || (v1.y > a && v2.y < a))
 221   {
 222     val = kFALSE;
 223     if (tolerance > 0)
 224     {
 225       Float_t a1 = TMath::Abs(v1.y - a), a2 = TMath::Abs(v2.y - a);
 226       if (a1 < a2)
 227       {
 228         if (a1 < tolerance) { v1.y = a; val = kTRUE; }
 229       }
 230       else
 231       {
 232         if (a2 < tolerance) { v2.y = a; val = kTRUE; }
 233       }
 234     }
 235   }
 236   return val;
 237 }
 238
 239
 240 //______________________________________________________________________________
 241 // NLTCircularFishEye
 242 //
 243 // XY projection with distortion around given center.
 244
 245 ClassImp(Reve::NLTCircularFishEye)
 246
 247 //______________________________________________________________________________
 248 void NLTCircularFishEye::ProjectPoint(Float_t& x, Float_t& y, Float_t& z,
 249                                       PProc_e proc)
 250 {
 251   // Project point.
 252
 253   using namespace TMath;
 254
 255   if (proc != PP_Plane)
 256   {
 257     x -= fCenter.x;
 258     y -= fCenter.y;
 259     Float_t phi = x == 0.0 && y == 0.0 ? 0.0 : ATan2(y,x);
 260     Float_t R = Sqrt(x*x+y*y);
 261     // distort
 262     Float_t NR = (R*fScale) / (1.0f + R*fDistortion);
 263     x = NR*Cos(phi) + fCenter.x;
 264     y = NR*Sin(phi) + fCenter.y;
 265   }
 266   z = 0.0f;
 267 }
 268