[3870] | 1 | /* -*- mode: c++; c-basic-offset: 4; indent-tabs-mode: nil -*- */ |
---|
| 2 | /* |
---|
| 3 | * Copyright (C) 2004-2013 HUBzero Foundation, LLC |
---|
| 4 | * |
---|
| 5 | */ |
---|
| 6 | |
---|
| 7 | #include <iostream> |
---|
| 8 | #include <float.h> |
---|
| 9 | |
---|
| 10 | #include <vtkSmartPointer.h> |
---|
| 11 | #include <vtkDataSet.h> |
---|
| 12 | #include <vtkProbeFilter.h> |
---|
| 13 | #include <vtkCutter.h> |
---|
| 14 | #include <vtkImageData.h> |
---|
| 15 | #include <vtkImageResample.h> |
---|
| 16 | #include <vtkPolyData.h> |
---|
| 17 | #include <vtkUnstructuredGrid.h> |
---|
| 18 | #include <vtkDelaunay2D.h> |
---|
| 19 | #include <vtkDelaunay3D.h> |
---|
| 20 | #include <vtkTransform.h> |
---|
| 21 | #include <vtkGaussianSplatter.h> |
---|
| 22 | #include <vtkExtractVOI.h> |
---|
| 23 | #include <vtkShepardMethod.h> |
---|
| 24 | |
---|
| 25 | #include "DataSetResample.h" |
---|
| 26 | #include "Trace.h" |
---|
| 27 | |
---|
| 28 | enum PrincipalPlane { |
---|
| 29 | PLANE_XY, |
---|
| 30 | PLANE_ZY, |
---|
| 31 | PLANE_XZ |
---|
| 32 | }; |
---|
| 33 | |
---|
| 34 | static bool is2D(vtkDataSet *dataSet, PrincipalPlane *plane, double *offset) |
---|
| 35 | { |
---|
| 36 | double bounds[6]; |
---|
| 37 | dataSet->GetBounds(bounds); |
---|
| 38 | if (bounds[4] == bounds[5]) { |
---|
| 39 | // Zmin = Zmax, XY plane |
---|
| 40 | if (plane != NULL) { |
---|
| 41 | *plane = PLANE_XY; |
---|
| 42 | } |
---|
| 43 | if (offset != NULL) |
---|
| 44 | *offset = bounds[4]; |
---|
| 45 | return true; |
---|
| 46 | } else if (bounds[0] == bounds[1]) { |
---|
| 47 | // Xmin = Xmax, ZY plane |
---|
| 48 | if (plane != NULL) { |
---|
| 49 | *plane = PLANE_ZY; |
---|
| 50 | } |
---|
| 51 | if (offset != NULL) |
---|
| 52 | *offset = bounds[0]; |
---|
| 53 | return true; |
---|
| 54 | } else if (bounds[2] == bounds[3]) { |
---|
| 55 | // Ymin = Ymax, XZ plane |
---|
| 56 | if (plane != NULL) { |
---|
| 57 | *plane = PLANE_XZ; |
---|
| 58 | } |
---|
| 59 | if (offset != NULL) |
---|
| 60 | *offset = bounds[2]; |
---|
| 61 | return true; |
---|
| 62 | } |
---|
| 63 | return false; |
---|
| 64 | } |
---|
| 65 | |
---|
| 66 | /** |
---|
| 67 | * \brief Determines if mesh is a point cloud (has no cells) |
---|
| 68 | */ |
---|
| 69 | static bool isCloud(vtkDataSet *dataSet) |
---|
| 70 | { |
---|
| 71 | vtkPolyData *pd = vtkPolyData::SafeDownCast(dataSet); |
---|
| 72 | if (pd) { |
---|
| 73 | // If PolyData, is a cloud if there are no cells other than vertices |
---|
| 74 | if (pd->GetNumberOfLines() == 0 && |
---|
| 75 | pd->GetNumberOfPolys() == 0 && |
---|
| 76 | pd->GetNumberOfStrips() == 0) { |
---|
| 77 | return true; |
---|
| 78 | } else { |
---|
| 79 | // Has cells |
---|
| 80 | return false; |
---|
| 81 | } |
---|
| 82 | } else { |
---|
| 83 | return (dataSet->GetNumberOfCells() == 0); |
---|
| 84 | } |
---|
| 85 | } |
---|
| 86 | |
---|
| 87 | vtkSmartPointer<vtkImageData> |
---|
| 88 | nv::resampleVTKDataSet(vtkDataSet *dataSetIn, int maxDim, nv::CloudStyle cloudStyle) |
---|
| 89 | { |
---|
| 90 | double bounds[6]; |
---|
| 91 | dataSetIn->GetBounds(bounds); |
---|
| 92 | |
---|
| 93 | int xdim, ydim, zdim; |
---|
| 94 | double xSpacing, ySpacing, zSpacing; |
---|
| 95 | |
---|
| 96 | vtkSmartPointer<vtkImageData> dataSetOut; |
---|
| 97 | if (vtkImageData::SafeDownCast(dataSetIn) != NULL) { |
---|
| 98 | // We have a uniform grid, use vtkImageResample |
---|
| 99 | vtkImageData *imgIn = vtkImageData::SafeDownCast(dataSetIn); |
---|
| 100 | |
---|
| 101 | TRACE("Image input Dims: %d,%d,%d Spacing: %g,%g,%g", |
---|
| 102 | imgIn->GetDimensions()[0], |
---|
| 103 | imgIn->GetDimensions()[1], |
---|
| 104 | imgIn->GetDimensions()[2], |
---|
| 105 | imgIn->GetSpacing()[0], |
---|
| 106 | imgIn->GetSpacing()[1], |
---|
| 107 | imgIn->GetSpacing()[2]); |
---|
| 108 | |
---|
| 109 | // Input DataSet is vtkImageData |
---|
| 110 | xdim = (imgIn->GetDimensions()[0] > maxDim) ? maxDim : imgIn->GetDimensions()[0]; |
---|
| 111 | ydim = (imgIn->GetDimensions()[1] > maxDim) ? maxDim : imgIn->GetDimensions()[1]; |
---|
| 112 | zdim = (imgIn->GetDimensions()[2] > maxDim) ? maxDim : imgIn->GetDimensions()[2]; |
---|
| 113 | |
---|
| 114 | xSpacing = (bounds[1]-bounds[0])/((double)(xdim-1)); |
---|
| 115 | ySpacing = (bounds[3]-bounds[2])/((double)(ydim-1)); |
---|
| 116 | zSpacing = (bounds[5]-bounds[4])/((double)(zdim-1)); |
---|
| 117 | |
---|
| 118 | vtkSmartPointer<vtkImageResample> filter = vtkSmartPointer<vtkImageResample>::New(); |
---|
| 119 | filter->SetInputData(dataSetIn); |
---|
| 120 | filter->SetAxisOutputSpacing(0, xSpacing); |
---|
| 121 | filter->SetAxisOutputSpacing(1, ySpacing); |
---|
| 122 | filter->SetAxisOutputSpacing(2, zSpacing); |
---|
| 123 | filter->Update(); |
---|
| 124 | dataSetOut = filter->GetOutput(); |
---|
| 125 | } else { |
---|
| 126 | double xLen = bounds[1] - bounds[0]; |
---|
| 127 | double yLen = bounds[3] - bounds[2]; |
---|
| 128 | double zLen = bounds[5] - bounds[4]; |
---|
| 129 | |
---|
| 130 | xdim = (xLen == 0.0 ? 1 : maxDim); |
---|
| 131 | ydim = (yLen == 0.0 ? 1 : maxDim); |
---|
| 132 | zdim = (zLen == 0.0 ? 1 : maxDim); |
---|
| 133 | |
---|
| 134 | xSpacing = (xLen == 0.0 ? 0.0 : xLen/((double)(xdim-1))); |
---|
| 135 | ySpacing = (yLen == 0.0 ? 0.0 : yLen/((double)(ydim-1))); |
---|
| 136 | zSpacing = (zLen == 0.0 ? 0.0 : zLen/((double)(zdim-1))); |
---|
| 137 | |
---|
| 138 | vtkSmartPointer<vtkImageData> imageData = vtkSmartPointer<vtkImageData>::New(); |
---|
| 139 | imageData->SetDimensions(xdim, ydim, zdim); |
---|
| 140 | imageData->SetOrigin(bounds[0], bounds[2], bounds[4]); |
---|
| 141 | imageData->SetSpacing(xSpacing, ySpacing, zSpacing); |
---|
| 142 | |
---|
| 143 | vtkSmartPointer<vtkProbeFilter> filter = vtkSmartPointer<vtkProbeFilter>::New(); |
---|
| 144 | filter->SetInputData(imageData); |
---|
| 145 | |
---|
| 146 | if (isCloud(dataSetIn)) { |
---|
| 147 | // Need to mesh or splat cloud |
---|
| 148 | if (cloudStyle == CLOUD_MESH) { |
---|
| 149 | PrincipalPlane plane; |
---|
| 150 | double offset; |
---|
| 151 | if (is2D(dataSetIn, &plane, &offset)) { |
---|
| 152 | vtkSmartPointer<vtkDelaunay2D> mesher = vtkSmartPointer<vtkDelaunay2D>::New(); |
---|
| 153 | if (plane == PLANE_ZY) { |
---|
| 154 | vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New(); |
---|
| 155 | trans->RotateWXYZ(90, 0, 1, 0); |
---|
| 156 | if (offset != 0.0) { |
---|
| 157 | trans->Translate(-offset, 0, 0); |
---|
| 158 | } |
---|
| 159 | mesher->SetTransform(trans); |
---|
| 160 | } else if (plane == PLANE_XZ) { |
---|
| 161 | vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New(); |
---|
| 162 | trans->RotateWXYZ(-90, 1, 0, 0); |
---|
| 163 | if (offset != 0.0) { |
---|
| 164 | trans->Translate(0, -offset, 0); |
---|
| 165 | } |
---|
| 166 | mesher->SetTransform(trans); |
---|
| 167 | } else if (offset != 0.0) { |
---|
| 168 | // XY with Z offset |
---|
| 169 | vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New(); |
---|
| 170 | trans->Translate(0, 0, -offset); |
---|
| 171 | mesher->SetTransform(trans); |
---|
| 172 | } |
---|
| 173 | mesher->SetInputData(dataSetIn); |
---|
| 174 | filter->SetSourceConnection(mesher->GetOutputPort()); |
---|
| 175 | //filter->SpatialMatchOn(); |
---|
| 176 | filter->Update(); |
---|
| 177 | dataSetOut = filter->GetImageDataOutput(); |
---|
| 178 | } else { |
---|
| 179 | vtkSmartPointer<vtkDelaunay3D> mesher = vtkSmartPointer<vtkDelaunay3D>::New(); |
---|
| 180 | mesher->SetInputData(dataSetIn); |
---|
| 181 | filter->SetSourceConnection(mesher->GetOutputPort()); |
---|
| 182 | //filter->SpatialMatchOn(); |
---|
| 183 | filter->Update(); |
---|
| 184 | dataSetOut = filter->GetImageDataOutput(); |
---|
| 185 | } |
---|
| 186 | } else if (cloudStyle == CLOUD_SPLAT) { |
---|
| 187 | // CLOUD_SPLAT |
---|
| 188 | vtkSmartPointer<vtkGaussianSplatter> splatter = vtkSmartPointer<vtkGaussianSplatter>::New(); |
---|
| 189 | splatter->SetInputData(dataSetIn); |
---|
| 190 | PrincipalPlane plane; |
---|
| 191 | double offset; |
---|
| 192 | int dims[3]; |
---|
| 193 | dims[0] = xdim; |
---|
| 194 | dims[1] = ydim; |
---|
| 195 | dims[2] = zdim; |
---|
| 196 | if (is2D(dataSetIn, &plane, &offset)) { |
---|
| 197 | vtkSmartPointer<vtkExtractVOI> volumeSlicer = vtkSmartPointer<vtkExtractVOI>::New(); |
---|
| 198 | if (plane == PLANE_ZY) { |
---|
| 199 | dims[0] = 3; |
---|
| 200 | volumeSlicer->SetVOI(1, 1, 0, dims[1]-1, 0, dims[1]-1); |
---|
| 201 | } else if (plane == PLANE_XZ) { |
---|
| 202 | dims[1] = 3; |
---|
| 203 | volumeSlicer->SetVOI(0, dims[0]-1, 1, 1, 0, dims[2]-1); |
---|
| 204 | } else { |
---|
| 205 | dims[2] = 3; |
---|
| 206 | volumeSlicer->SetVOI(0, dims[0]-1, 0, dims[1]-1, 1, 1); |
---|
| 207 | } |
---|
| 208 | splatter->SetSampleDimensions(dims); |
---|
| 209 | volumeSlicer->SetInputConnection(splatter->GetOutputPort()); |
---|
| 210 | volumeSlicer->SetSampleRate(1, 1, 1); |
---|
| 211 | volumeSlicer->Update(); |
---|
| 212 | dataSetOut = volumeSlicer->GetOutput(); |
---|
| 213 | } else { |
---|
| 214 | splatter->SetSampleDimensions(dims); |
---|
| 215 | splatter->Update(); |
---|
| 216 | dataSetOut = splatter->GetOutput(); |
---|
| 217 | } |
---|
| 218 | } else { |
---|
| 219 | // CLOUD_SHEAPRD_METHOD |
---|
| 220 | vtkSmartPointer<vtkShepardMethod> splatter = vtkSmartPointer<vtkShepardMethod>::New(); |
---|
| 221 | splatter->SetInputData(dataSetIn); |
---|
| 222 | PrincipalPlane plane; |
---|
| 223 | double offset; |
---|
| 224 | int dims[3]; |
---|
| 225 | dims[0] = xdim; |
---|
| 226 | dims[1] = ydim; |
---|
| 227 | dims[2] = zdim; |
---|
| 228 | if (is2D(dataSetIn, &plane, &offset)) { |
---|
| 229 | vtkSmartPointer<vtkExtractVOI> volumeSlicer = vtkSmartPointer<vtkExtractVOI>::New(); |
---|
| 230 | if (plane == PLANE_ZY) { |
---|
| 231 | dims[0] = 3; |
---|
| 232 | volumeSlicer->SetVOI(1, 1, 0, dims[1]-1, 0, dims[1]-1); |
---|
| 233 | } else if (plane == PLANE_XZ) { |
---|
| 234 | dims[1] = 3; |
---|
| 235 | volumeSlicer->SetVOI(0, dims[0]-1, 1, 1, 0, dims[2]-1); |
---|
| 236 | } else { |
---|
| 237 | dims[2] = 3; |
---|
| 238 | volumeSlicer->SetVOI(0, dims[0]-1, 0, dims[1]-1, 1, 1); |
---|
| 239 | } |
---|
| 240 | splatter->SetSampleDimensions(dims); |
---|
| 241 | volumeSlicer->SetInputConnection(splatter->GetOutputPort()); |
---|
| 242 | volumeSlicer->SetSampleRate(1, 1, 1); |
---|
| 243 | volumeSlicer->Update(); |
---|
| 244 | dataSetOut = volumeSlicer->GetOutput(); |
---|
| 245 | } else { |
---|
| 246 | splatter->SetSampleDimensions(dims); |
---|
| 247 | splatter->Update(); |
---|
| 248 | dataSetOut = splatter->GetOutput(); |
---|
| 249 | } |
---|
| 250 | } |
---|
| 251 | } else { |
---|
| 252 | filter->SetSourceData(dataSetIn); |
---|
| 253 | //filter->SpatialMatchOn(); |
---|
| 254 | filter->Update(); |
---|
| 255 | dataSetOut = filter->GetImageDataOutput(); |
---|
| 256 | } |
---|
| 257 | } |
---|
| 258 | |
---|
| 259 | TRACE("Image output Dims: %d,%d,%d Spacing: %g,%g,%g", |
---|
| 260 | dataSetOut->GetDimensions()[0], |
---|
| 261 | dataSetOut->GetDimensions()[1], |
---|
| 262 | dataSetOut->GetDimensions()[2], |
---|
| 263 | dataSetOut->GetSpacing()[0], |
---|
| 264 | dataSetOut->GetSpacing()[1], |
---|
| 265 | dataSetOut->GetSpacing()[2]); |
---|
| 266 | |
---|
| 267 | return dataSetOut; |
---|
| 268 | } |
---|