1 /*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkShepardKernel.cxx
5
6 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7 All rights reserved.
8 See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9
10 This software is distributed WITHOUT ANY WARRANTY; without even
11 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12 PURPOSE. See the above copyright notice for more information.
13
14 =========================================================================*/
15 #include "vtkShepardKernel.h"
16 #include "vtkAbstractPointLocator.h"
17 #include "vtkObjectFactory.h"
18 #include "vtkIdList.h"
19 #include "vtkDoubleArray.h"
20 #include "vtkDataSet.h"
21 #include "vtkPointData.h"
22 #include "vtkMath.h"
23 #include "vtkMathUtilities.h"
24
25 vtkStandardNewMacro(vtkShepardKernel);
26
27 //----------------------------------------------------------------------------
vtkShepardKernel()28 vtkShepardKernel::vtkShepardKernel()
29 {
30 this->PowerParameter = 2.0;
31 }
32
33
34 //----------------------------------------------------------------------------
35 vtkShepardKernel::~vtkShepardKernel() = default;
36
37
38 //----------------------------------------------------------------------------
39 vtkIdType vtkShepardKernel::
ComputeWeights(double x[3],vtkIdList * pIds,vtkDoubleArray * prob,vtkDoubleArray * weights)40 ComputeWeights(double x[3], vtkIdList *pIds, vtkDoubleArray *prob,
41 vtkDoubleArray *weights)
42 {
43 vtkIdType numPts = pIds->GetNumberOfIds();
44 double d, y[3], sum=0.0;
45 weights->SetNumberOfTuples(numPts);
46 double *p = (prob ? prob->GetPointer(0) : nullptr);
47 double *w = weights->GetPointer(0);
48
49 for (vtkIdType i=0; i<numPts; ++i)
50 {
51 vtkIdType id = pIds->GetId(i);
52 this->DataSet->GetPoint(id,y);
53 if ( this->PowerParameter == 2.0 )
54 {
55 d = vtkMath::Distance2BetweenPoints(x,y);
56 }
57 else
58 {
59 d = pow(sqrt(vtkMath::Distance2BetweenPoints(x,y)), this->PowerParameter);
60 }
61 if ( vtkMathUtilities::FuzzyCompare(d, 0.0, std::numeric_limits<double>::epsilon()*256.0 )) //precise hit on existing point
62 {
63 pIds->SetNumberOfIds(1);
64 pIds->SetId(0,id);
65 weights->SetNumberOfTuples(1);
66 weights->SetValue(0,1.0);
67 return 1;
68 }
69 else
70 {
71 w[i] = (p ? p[i]/d : 1.0/d); //take into account probability if provided
72 sum += w[i];
73 }
74 }//over all points
75
76 // Normalize
77 if ( this->NormalizeWeights && sum != 0.0 )
78 {
79 for (vtkIdType i=0; i<numPts; ++i)
80 {
81 w[i] /= sum;
82 }
83 }
84
85 return numPts;
86 }
87
88 //----------------------------------------------------------------------------
PrintSelf(ostream & os,vtkIndent indent)89 void vtkShepardKernel::PrintSelf(ostream& os, vtkIndent indent)
90 {
91 this->Superclass::PrintSelf(os,indent);
92
93 os << indent << "Power Parameter: "
94 << this->GetPowerParameter() << "\n";
95 }
96