VTK  9.3.0
vtkBlockSortHelper.h
Go to the documentation of this file.
1// SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2// SPDX-License-Identifier: BSD-3-Clause
8#ifndef vtkBlockSortHelper_h
9#define vtkBlockSortHelper_h
10
11#include "vtkCamera.h"
12#include "vtkDataSet.h"
13#include "vtkMatrix4x4.h"
14#include "vtkNew.h"
15#include "vtkRenderer.h"
16#include "vtkVector.h"
17#include "vtkVectorOperators.h"
18
19#include <vector>
20
22{
23VTK_ABI_NAMESPACE_BEGIN
24
25template <typename T>
26inline void GetBounds(T a, double bds[6])
27{
28 a.GetBounds(bds);
29}
30
31template <>
32inline void GetBounds(vtkDataSet* first, double bds[6])
33{
34 first->GetBounds(bds);
35}
36
43template <typename T>
45{
49
50 //----------------------------------------------------------------------------
52 {
53 vtkCamera* cam = ren->GetActiveCamera();
54 this->CameraIsParallel = (cam->GetParallelProjection() != 0);
55
56 double camWorldPos[4];
57 cam->GetPosition(camWorldPos);
58 camWorldPos[3] = 1.0;
59
60 double camWorldFocalPoint[4];
61 cam->GetFocalPoint(camWorldFocalPoint);
62 camWorldFocalPoint[3] = 1.0;
63
64 // Transform the camera position to the volume (dataset) coordinate system.
65 vtkNew<vtkMatrix4x4> InverseVolumeMatrix;
66 InverseVolumeMatrix->DeepCopy(volMatrix);
67 InverseVolumeMatrix->Invert();
68 InverseVolumeMatrix->MultiplyPoint(camWorldPos, camWorldPos);
69 InverseVolumeMatrix->MultiplyPoint(camWorldFocalPoint, camWorldFocalPoint);
70
71 this->CameraPosition = vtkVector3d(camWorldPos[0], camWorldPos[1], camWorldPos[2]);
72 this->CameraPosition = this->CameraPosition / vtkVector3d(camWorldPos[3]);
73
74 vtkVector3d camFP(camWorldFocalPoint[0], camWorldFocalPoint[1], camWorldFocalPoint[2]);
75 camFP = camFP / vtkVector3d(camWorldFocalPoint[3]);
76
77 this->CameraViewDirection = camFP - this->CameraPosition;
78 }
79
80 //----------------------------------------------------------------------------
81 // camPos is used when is_parallel is false, else viewdirection is used.
82 // thus a valid camPos is only needed if is_parallel is false, and a valid viewdirection
83 // is only needed if is_parallel is true.
84 BackToFront(const vtkVector3d& camPos, const vtkVector3d& viewdirection, bool is_parallel)
85 : CameraPosition(camPos)
86 , CameraViewDirection(viewdirection)
87 , CameraIsParallel(is_parallel)
88 {
89 }
90
91 // -1 if first is closer than second
92 // 0 if unknown
93 // 1 if second is farther than first
94 template <typename TT>
95 inline int CompareOrderWithUncertainty(TT& first, TT& second)
96 {
97 double abounds[6], bbounds[6];
98 vtkBlockSortHelper::GetBounds<TT>(first, abounds);
99 vtkBlockSortHelper::GetBounds<TT>(second, bbounds);
100 return CompareBoundsOrderWithUncertainty(abounds, bbounds);
101 }
102
103 // -1 if first is closer than second
104 // 0 if unknown
105 // 1 if second is farther than first
106 inline int CompareBoundsOrderWithUncertainty(const double abounds[6], const double bbounds[6])
107 {
108 double bboundsP[6];
109 double aboundsP[6];
110 for (int i = 0; i < 6; ++i)
111 {
112 int low = 2 * (i / 2);
113 bboundsP[i] = bbounds[i];
114 if (bboundsP[i] < abounds[low])
115 {
116 bboundsP[i] = abounds[low];
117 }
118 if (bboundsP[i] > abounds[low + 1])
119 {
120 bboundsP[i] = abounds[low + 1];
121 }
122 aboundsP[i] = abounds[i];
123 if (aboundsP[i] < bbounds[low])
124 {
125 aboundsP[i] = bbounds[low];
126 }
127 if (aboundsP[i] > bbounds[low + 1])
128 {
129 aboundsP[i] = bbounds[low + 1];
130 }
131 }
132
133 int dims = 0;
134 int degenDims = 0;
135 int degenAxes[3];
136 int dimSize[3];
137 for (int i = 0; i < 6; i += 2)
138 {
139 if (aboundsP[i] != aboundsP[i + 1])
140 {
141 dimSize[dims] = aboundsP[i + 1] - aboundsP[i];
142 dims++;
143 }
144 else
145 {
146 degenAxes[degenDims] = i / 2;
147 degenDims++;
148 }
149 }
150
151 // overlapping volumes?
152 // in that case find the two largest dimensions
153 // geneally this should not happen
154 if (dims == 3)
155 {
156 if (dimSize[0] < dimSize[1])
157 {
158 if (dimSize[0] < dimSize[2])
159 {
160 degenAxes[0] = 0;
161 }
162 else
163 {
164 degenAxes[0] = 2;
165 }
166 }
167 else
168 {
169 if (dimSize[1] < dimSize[2])
170 {
171 degenAxes[0] = 1;
172 }
173 else
174 {
175 degenAxes[0] = 2;
176 }
177 }
178 dims = 2;
179 degenDims = 1;
180 }
181
182 // compute the direction from a to b
183 double atobdir[3] = { bbounds[0] + bbounds[1] - abounds[0] - abounds[1],
184 bbounds[2] + bbounds[3] - abounds[2] - abounds[3],
185 bbounds[4] + bbounds[5] - abounds[4] - abounds[5] };
186 double atoblength = vtkMath::Normalize(atobdir);
187
188 // no comment on blocks that do not touch
189 if (fabs(aboundsP[degenAxes[0] * 2] - bboundsP[degenAxes[0] * 2]) > 0.01 * atoblength)
190 {
191 return 0;
192 }
193
194 // compute the camera direction
195 vtkVector3d dir;
196 if (this->CameraIsParallel)
197 {
198 dir = this->CameraViewDirection;
199 }
200 else
201 {
202 // compute a point for the half space plane
203 vtkVector3d planePoint(0.25 * (aboundsP[0] + aboundsP[1] + bboundsP[0] + bboundsP[1]),
204 0.25 * (aboundsP[2] + aboundsP[3] + bboundsP[2] + bboundsP[3]),
205 0.25 * (aboundsP[4] + aboundsP[5] + bboundsP[4] + bboundsP[5]));
206 dir = planePoint - this->CameraPosition;
207 }
208 dir.Normalize();
209
210 // planar interface
211 if (dims == 2)
212 {
213 double plane[3] = { 0, 0, 0 };
214 plane[degenAxes[0]] = 1.0;
215 // dot product with camera direction
216 double dot = dir[0] * plane[0] + dir[1] * plane[1] + dir[2] * plane[2];
217 if (dot == 0)
218 {
219 return 0;
220 }
221 // figure out what side we are on
222 double side = plane[0] * atobdir[0] + plane[1] * atobdir[1] + plane[2] * atobdir[2];
223 return (dot * side < 0 ? 1 : -1);
224 }
225
226 return 0;
227 }
228};
229
230#ifdef MB_DEBUG
231template <class RandomIt>
232class gnode
233{
234public:
235 RandomIt Value;
236 bool Visited = false;
237 std::vector<gnode<RandomIt>*> Closer;
238};
239
240template <class RandomIt>
241bool operator==(gnode<RandomIt> const& lhs, gnode<RandomIt> const& rhs)
242{
243 return lhs.Value == rhs.Value && lhs.Closer == rhs.Closer;
244}
245
246template <class RandomIt>
247bool findCycle(gnode<RandomIt>& start, std::vector<gnode<RandomIt>>& graph,
248 std::vector<gnode<RandomIt>>& active, std::vector<gnode<RandomIt>>& loop)
249{
250 if (start.Visited)
251 {
252 return false;
253 }
254
255 // add the current node to active list
256 active.push_back(start);
257
258 // traverse the closer nodes one by one depth first
259 for (auto& close : start.Closer)
260 {
261 if (close->Visited)
262 {
263 continue;
264 }
265
266 // is the node already in the active list? if so we have a loop
267 for (auto ait = active.begin(); ait != active.end(); ++ait)
268 {
269 if (ait->Value == close->Value)
270 {
271 loop.push_back(*ait);
272 return true;
273 }
274 }
275 // otherwise recurse
276 if (findCycle(*close, graph, active, loop))
277 {
278 // a loop was detected, build the loop output
279 loop.push_back(*close);
280 return true;
281 }
282 }
283
284 active.erase(std::find(active.begin(), active.end(), start));
285 start.Visited = true;
286 return false;
287}
288#endif
289
290template <class RandomIt, typename T>
291inline void Sort(RandomIt bitr, RandomIt eitr, BackToFront<T>& me)
292{
293 auto start = bitr;
294
295 // brute force for testing
296
297 std::vector<typename RandomIt::value_type> working;
298 std::vector<typename RandomIt::value_type> result;
299 working.assign(bitr, eitr);
300 size_t numNodes = working.size();
301
302#ifdef MB_DEBUG
303 // check for any short loops and debug
304 for (auto it = working.begin(); it != working.end(); ++it)
305 {
306 auto it2 = it;
307 it2++;
308 for (; it2 != working.end(); ++it2)
309 {
310 int comp1 = me.CompareOrderWithUncertainty(*it, *it2);
311 int comp2 = me.CompareOrderWithUncertainty(*it2, *it);
312 if (comp1 * comp2 > 0)
313 {
314 me.CompareOrderWithUncertainty(*it, *it2);
315 me.CompareOrderWithUncertainty(*it2, *it);
316 }
317 }
318 }
319
320 // build the graph
321 std::vector<gnode<RandomIt>> graph;
322 for (auto it = working.begin(); it != working.end(); ++it)
323 {
324 gnode<RandomIt> anode;
325 anode.Value = it;
326 graph.push_back(anode);
327 }
328 for (auto& git : graph)
329 {
330 for (auto& next : graph)
331 {
332 if (git.Value != next.Value && me.CompareOrderWithUncertainty(*git.Value, *next.Value) > 0)
333 {
334 git.Closer.push_back(&next);
335 }
336 }
337 }
338
339 // graph constructed, now look for a loop
340 std::vector<gnode<RandomIt>> active;
341 std::vector<gnode<RandomIt>> loop;
342 for (auto& gval : graph)
343 {
344 loop.clear();
345 if (findCycle(gval, graph, active, loop))
346 {
348 dir.Normalize();
349 vtkGenericWarningMacro("found a loop cam dir: " << dir[0] << " " << dir[1] << " " << dir[2]);
350 for (auto& lval : loop)
351 {
352 double bnds[6];
353 vtkBlockSortHelper::GetBounds((*lval.Value), bnds);
354 vtkGenericWarningMacro(<< bnds[0] << " " << bnds[1] << " " << bnds[2] << " " << bnds[3]
355 << " " << bnds[4] << " " << bnds[5]);
356 }
357 }
358 }
359#endif
360
361 // loop over items and find the first that is not in front of any others
362 for (auto it = working.begin(); it != working.end();)
363 {
364 auto it2 = working.begin();
365 for (; it2 != working.end(); ++it2)
366 {
367 // if another block is in front of this block then this is not the
368 // closest block
369 if (it != it2 && me.CompareOrderWithUncertainty(*it, *it2) > 0)
370 {
371 // not a winner
372 break;
373 }
374 }
375 if (it2 == working.end())
376 {
377 // found a winner, add it to the results, remove from the working set and then restart
378 result.push_back(*it);
379 working.erase(it);
380 it = working.begin();
381 }
382 else
383 {
384 ++it;
385 }
386 }
387
388 if (result.size() != numNodes)
389 {
390 vtkGenericWarningMacro("sorting failed");
391 }
392
393 // copy results to original container
394 std::reverse_copy(result.begin(), result.end(), start);
395};
396VTK_ABI_NAMESPACE_END
397}
398
399#endif // vtkBlockSortHelper_h
400// VTK-HeaderTest-Exclude: vtkBlockSortHelper.h
a virtual camera for 3D rendering
Definition vtkCamera.h:41
virtual double * GetFocalPoint()
Set/Get the focal of the camera in world coordinates.
virtual double * GetPosition()
Set/Get the position of the camera in world coordinates.
virtual vtkTypeBool GetParallelProjection()
Set/Get the value of the ParallelProjection instance variable.
abstract class to specify dataset behavior
Definition vtkDataSet.h:53
double * GetBounds()
Return a pointer to the geometry bounding box in the form (xmin,xmax, ymin,ymax, zmin,...
static float Normalize(float v[3])
Normalize (in place) a 3-vector.
Definition vtkMath.h:1731
represent and manipulate 4x4 transformation matrices
Allocate and hold a VTK object.
Definition vtkNew.h:51
abstract specification for renderers
Definition vtkRenderer.h:59
vtkCamera * GetActiveCamera()
Get the current camera.
double Normalize()
Normalize the vector in place.
Definition vtkVector.h:77
Collection of comparison functions for std::sort.
void Sort(RandomIt bitr, RandomIt eitr, BackToFront< T > &me)
void GetBounds(T a, double bds[6])
@ loop
Definition vtkX3D.h:407
operator() for back-to-front sorting.
int CompareOrderWithUncertainty(TT &first, TT &second)
int CompareBoundsOrderWithUncertainty(const double abounds[6], const double bbounds[6])
BackToFront(vtkRenderer *ren, vtkMatrix4x4 *volMatrix)
BackToFront(const vtkVector3d &camPos, const vtkVector3d &viewdirection, bool is_parallel)
bool VTKCOMMONCORE_EXPORT operator==(const std::string &a, const vtkStringToken &b)