Skip to content

Commit

Permalink
Refs #7363. Change to vtkStructuredGrid.
Browse files Browse the repository at this point in the history
  • Loading branch information
quantumsteve committed Dec 9, 2014
1 parent 2149a09 commit 72faff7
Show file tree
Hide file tree
Showing 11 changed files with 275 additions and 127 deletions.
Expand Up @@ -117,20 +117,10 @@ int vtkMDHWSource::RequestData(vtkInformation *, vtkInformationVector **, vtkInf

vtkDataSet* product = m_presenter->execute(factory, loadingProgressUpdate, drawingProgressUpdate);

//-------------------------------------------------------- Corrects problem whereby boundaries not set propertly in PV.
vtkBox* box = vtkBox::New();
box->SetBounds(product->GetBounds());
vtkPVClipDataSet* clipper = vtkPVClipDataSet::New();
clipper->SetInputData(product);
clipper->SetClipFunction(box);
clipper->SetInsideOut(true);
clipper->Update();
vtkDataSet* clipperOutput = clipper->GetOutput();
//--------------------------------------------------------

vtkUnstructuredGrid *output = vtkUnstructuredGrid::SafeDownCast(
outInfo->Get(vtkDataObject::DATA_OBJECT()));
output->ShallowCopy(clipperOutput);
vtkDataSet* output = vtkDataSet::GetData(outInfo);
output->ShallowCopy(product);
product->Delete();

try
{
m_presenter->makeNonOrthogonal(output);
Expand All @@ -142,8 +132,7 @@ int vtkMDHWSource::RequestData(vtkInformation *, vtkInformationVector **, vtkInf
<< "plot non-orthogonal axes. " << error);
}
m_presenter->setAxisLabels(output);

clipper->Delete();

}
return 1;
}
Expand All @@ -153,16 +142,23 @@ int vtkMDHWSource::RequestInformation(vtkInformation *vtkNotUsed(request), vtkIn
if(m_presenter == NULL && !m_wsName.empty())
{
m_presenter = new MDHWInMemoryLoadingPresenter(new MDLoadingViewAdapter<vtkMDHWSource>(this), new ADSWorkspaceProvider<Mantid::API::IMDHistoWorkspace>, m_wsName);
if(!m_presenter->canReadFile())
{
vtkErrorMacro(<<"Cannot fetch the specified workspace from Mantid ADS.");
}
else
{
m_presenter->executeLoadMetadata();
setTimeRange(outputVector);
}
}
if (m_presenter == NULL)
{
// updater information has been called prematurely. We will reexecute once all attributes are setup.
return 1;
}
if(!m_presenter->canReadFile())
{
vtkErrorMacro(<<"Cannot fetch the specified workspace from Mantid ADS.");
return 0;
}

m_presenter->executeLoadMetadata();
setTimeRange(outputVector);
std::vector<int> extents = dynamic_cast<MDHWInMemoryLoadingPresenter*>(m_presenter)->getExtents();
outputVector->GetInformationObject(0)->Set(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT(),&extents[0],extents.size());

return 1;
}

Expand Down
@@ -1,7 +1,7 @@
#ifndef _vtkMDHWSource_h
#define _vtkMDHWSource_h

#include "vtkUnstructuredGridAlgorithm.h"
#include "vtkStructuredGridAlgorithm.h"
#include <string>

namespace Mantid
Expand Down Expand Up @@ -40,11 +40,11 @@ namespace Mantid
*/

// cppcheck-suppress class_X_Y
class VTK_EXPORT vtkMDHWSource : public vtkUnstructuredGridAlgorithm
class VTK_EXPORT vtkMDHWSource : public vtkStructuredGridAlgorithm
{
public:
static vtkMDHWSource *New();
vtkTypeMacro(vtkMDHWSource, vtkUnstructuredGridAlgorithm);
vtkTypeMacro(vtkMDHWSource, vtkStructuredGridAlgorithm);
void PrintSelf(ostream& os, vtkIndent indent);

void SetWsName(std::string wsName);
Expand Down
Expand Up @@ -3,6 +3,7 @@

#include "MantidVatesAPI/MDHWLoadingPresenter.h"
#include <boost/scoped_ptr.hpp>
#include <vector>

class vtkDataSet;
namespace Mantid
Expand Down Expand Up @@ -48,6 +49,7 @@ namespace Mantid
virtual bool canReadFile() const;
virtual std::string getWorkspaceTypeName();
virtual int getSpecialCoordinates();
std::vector<int> getExtents();
private:
/// Repository for accessing workspaces. At this level, does not specify how or where from.
boost::scoped_ptr<WorkspaceProvider> m_repository;
Expand Down
Expand Up @@ -12,7 +12,6 @@
#include <vector>

class vtkDataSet;
class vtkUnstructuredGrid;

namespace Mantid
{
Expand Down Expand Up @@ -75,7 +74,7 @@ namespace VATES
/// Reduce the dimensionality of matrix by 1
void stripMatrix(Kernel::DblMatrix &mat);
/// Add the skew basis to metadata
void updateMetaData(vtkUnstructuredGrid *ugrid);
void updateMetaData(vtkDataSet *ugrid);
vtkDataSet *m_dataSet; ///< Pointer to VTK dataset to modify
std::string m_wsName; ///< The name of the workspace to fetch
//FIXME: Temp var for getting hardcoded stuff back
Expand Down
Expand Up @@ -31,7 +31,7 @@
#include "MantidVatesAPI/vtkDataSetFactory.h"
#include "MantidAPI/IMDWorkspace.h"
#include "MantidVatesAPI/ThresholdRange.h"
#include <vtkUnstructuredGrid.h>
//#include <vtkUnstructuredGrid.h>
#include <vtkFloatArray.h>
#include <vtkCellData.h>
#include <vtkHexahedron.h>
Expand Down
14 changes: 14 additions & 0 deletions Code/Mantid/Vates/VatesAPI/src/MDHWInMemoryLoadingPresenter.cpp
Expand Up @@ -125,5 +125,19 @@ namespace Mantid
{
return m_specialCoords;
}

std::vector<int> MDHWInMemoryLoadingPresenter::getExtents()
{
// Hack which only works in 3D. Needs to be updated for 4 dimensions!
using namespace Mantid::API;
Workspace_sptr ws = m_repository->fetchWorkspace(m_wsName);
IMDHistoWorkspace_sptr histoWs = boost::dynamic_pointer_cast<Mantid::API::IMDHistoWorkspace>(ws);
std::vector<int> extents(6, 0);
extents[1] = histoWs->getDimension(0)->getNBins();
extents[3] = histoWs->getDimension(1)->getNBins();
extents[5] = histoWs->getDimension(2)->getNBins();

return extents;
}
}
}
Expand Up @@ -73,7 +73,7 @@ vtkDataSetToNonOrthogonalDataSet::~vtkDataSetToNonOrthogonalDataSet()
void vtkDataSetToNonOrthogonalDataSet::execute()
{
// Downcast to a vtkUnstructuredGrid
vtkUnstructuredGrid *data = vtkUnstructuredGrid::SafeDownCast(m_dataSet);
vtkPointSet *data = vtkPointSet::SafeDownCast(m_dataSet);
if (NULL == data)
{
throw std::runtime_error("VTK dataset does not inherit from vtkPointSet");
Expand Down Expand Up @@ -373,7 +373,7 @@ void vtkDataSetToNonOrthogonalDataSet::copyToRaw(double *arr, MantidVec vec)
* VTK dataset.
* @param ugrid : The VTK dataset to add the metadata to
*/
void vtkDataSetToNonOrthogonalDataSet::updateMetaData(vtkUnstructuredGrid *ugrid)
void vtkDataSetToNonOrthogonalDataSet::updateMetaData(vtkDataSet *ugrid)
{
double baseX[3];
double baseY[3];
Expand Down
1 change: 1 addition & 0 deletions Code/Mantid/Vates/VatesAPI/src/vtkMDHexFactory.cpp
Expand Up @@ -4,6 +4,7 @@
#include "MantidVatesAPI/vtkMDHexFactory.h"
#include "MantidVatesAPI/Common.h"
#include "MantidVatesAPI/ProgressAction.h"
#include "MantidVatesAPI/NoThresholdRange.h"
#include <vtkCellData.h>
#include <vtkFloatArray.h>
#include <vtkHexahedron.h>
Expand Down
130 changes: 36 additions & 94 deletions Code/Mantid/Vates/VatesAPI/src/vtkMDHistoHexFactory.cpp
Expand Up @@ -7,6 +7,10 @@
#include "MantidAPI/NullCoordTransform.h"
#include "MantidKernel/ReadLock.h"

#include "vtkNew.h"
#include "vtkStructuredGrid.h"
#include "vtkDoubleArray.h"

using Mantid::API::IMDWorkspace;
using Mantid::API::IMDHistoWorkspace;
using Mantid::Kernel::CPUTimer;
Expand Down Expand Up @@ -106,11 +110,22 @@ namespace VATES
vtkPoints *points = vtkPoints::New();
points->Allocate(static_cast<int>(imageSize));

vtkFloatArray * signal = vtkFloatArray::New();
signal->Allocate(imageSize);
//vtkFloatArray * signal = vtkFloatArray::New();
//signal->Allocate(imageSize);
//signal->SetName(m_scalarName.c_str());
//signal->SetNumberOfComponents(1);

vtkNew<vtkDoubleArray> signal;
signal->SetName(m_scalarName.c_str());
signal->SetNumberOfComponents(1);
// (argument #3) tell VTK to not delete this array
signal->SetArray(m_workspace->getSignalArray(),imageSize,1);
assert(imageSize == signal->GetNumberOfTuples());

vtkNew<vtkStructuredGrid> visualDataset;
visualDataset->SetDimensions(nBinsX+1,nBinsY+1,nBinsZ+1);


double signalScalar;
const int nPointsX = nBinsX+1;
const int nPointsY = nBinsY+1;
Expand All @@ -122,15 +137,13 @@ namespace VATES
create the points that will be needed; so an array of pointNeeded
is set so that all required vertices are marked, and created in a second step. */

// Array of the points that should be created, set to false
bool * pointNeeded = new bool[nPointsX*nPointsY*nPointsZ];
memset(pointNeeded, 0, nPointsX*nPointsY*nPointsZ*sizeof(bool));
// Array with true where the voxel should be shown
bool * voxelShown = new bool[nBinsX*nBinsY*nBinsZ];
double progressFactor = 0.5/double(nBinsZ);
double progressOffset = 0.5;

size_t index = 0;


vtkIdType index = 0;
for (int z = 0; z < nBinsZ; z++)
{
//Report progress updates for the first 50%
Expand All @@ -150,27 +163,10 @@ namespace VATES
else
signalScalar = m_workspace->getSignalNormalizedAt(x,y,z);

if (isSpecial( signalScalar ) || !m_thresholdRange->inRange(signalScalar))
{
// out of range
voxelShown[index] = false;
}
else
bool maskValue = (isSpecial( signalScalar ) || !m_thresholdRange->inRange(signalScalar));
if (maskValue)
{
// Valid data
voxelShown[index] = true;
signal->InsertNextValue(static_cast<float>(signalScalar));

// Make sure all 8 neighboring points are set to true
size_t pointIndex = x + (nPointsX * y) + (nPointsX*nPointsY*z); //(Note this index is different then the other one)
pointNeeded[pointIndex] = true; pointIndex++;
pointNeeded[pointIndex] = true; pointIndex += nPointsX-1;
pointNeeded[pointIndex] = true; pointIndex++;
pointNeeded[pointIndex] = true; pointIndex += nPointsX*nPointsY - nPointsX - 1;
pointNeeded[pointIndex] = true; pointIndex++;
pointNeeded[pointIndex] = true; pointIndex += nPointsX-1;
pointNeeded[pointIndex] = true; pointIndex++;
pointNeeded[pointIndex] = true;
visualDataset->BlankCell(index);
}
index++;
}
Expand All @@ -188,7 +184,6 @@ namespace VATES
Mantid::coord_t out[3];

// Array with the point IDs (only set where needed)
vtkIdType * pointIDs = new vtkIdType[nPointsX*nPointsY*nPointsZ];
index = 0;
progressFactor = 0.5/static_cast<double>(nPointsZ);

Expand All @@ -203,84 +198,31 @@ namespace VATES
for (int x = 0; x < nPointsX; x++)
{
// Create the point only when needed
if (pointNeeded[index])
in[0] = (minX + (static_cast<coord_t>(x) * incrementX)); //Calculate increment in x;
if (transform)
{
in[0] = (minX + (static_cast<coord_t>(x) * incrementX)); //Calculate increment in x;
if (transform)
{
transform->apply(in, out);
pointIDs[index] = points->InsertNextPoint(out);
}
else
{
pointIDs[index] = points->InsertNextPoint(in);
}
transform->apply(in, out);
points->InsertNextPoint(out);
}
index++;
}
}
}

std::cout << tim << " to create the needed points." << std::endl;

vtkUnstructuredGrid *visualDataSet = vtkUnstructuredGrid::New();
visualDataSet->Allocate(imageSize);
visualDataSet->SetPoints(points);
visualDataSet->GetCellData()->SetScalars(signal);

// ------ Hexahedron creation ----------------
// It is approx. 40 x faster to create the hexadron only once, and reuse it for each voxel.
vtkHexahedron *theHex = vtkHexahedron::New();
index = 0;

for (int z = 0; z < nBinsZ; z++)
{
for (int y = 0; y < nBinsY; y++)
{
for (int x = 0; x < nBinsX; x++)
{
if (voxelShown[index])
else
{
//Only create topologies for those cells which are not sparse.
// create a hexahedron topology
vtkIdType id_xyz = pointIDs[(x) + (y)*nPointsX + z*nPointsX*nPointsY];
vtkIdType id_dxyz = pointIDs[(x+1) + (y)*nPointsX + z*nPointsX*nPointsY];
vtkIdType id_dxdyz = pointIDs[(x+1) + (y+1)*nPointsX + z*nPointsX*nPointsY];
vtkIdType id_xdyz = pointIDs[(x) + (y+1)*nPointsX + z*nPointsX*nPointsY];

vtkIdType id_xydz = pointIDs[(x) + (y)*nPointsX + (z+1)*nPointsX*nPointsY];
vtkIdType id_dxydz = pointIDs[(x+1) + (y)*nPointsX + (z+1)*nPointsX*nPointsY];
vtkIdType id_dxdydz = pointIDs[(x+1) + (y+1)*nPointsX + (z+1)*nPointsX*nPointsY];
vtkIdType id_xdydz = pointIDs[(x) + (y+1)*nPointsX + (z+1)*nPointsX*nPointsY];

//create the hexahedron
theHex->GetPointIds()->SetId(0, id_xyz);
theHex->GetPointIds()->SetId(1, id_dxyz);
theHex->GetPointIds()->SetId(2, id_dxdyz);
theHex->GetPointIds()->SetId(3, id_xdyz);
theHex->GetPointIds()->SetId(4, id_xydz);
theHex->GetPointIds()->SetId(5, id_dxydz);
theHex->GetPointIds()->SetId(6, id_dxdydz);
theHex->GetPointIds()->SetId(7, id_xdydz);

visualDataSet->InsertNextCell(VTK_HEXAHEDRON, theHex->GetPointIds());
points->InsertNextPoint(in);
}
index++;
}
}
}
theHex->Delete();

std::cout << tim << " to create the needed points." << std::endl;
visualDataset->SetPoints(points);
visualDataset->GetCellData()->SetScalars(signal.GetPointer());
visualDataset->Register(NULL);

std::cout << tim << " to create and add the hexadrons." << std::endl;


points->Delete();
signal->Delete();
visualDataSet->Squeeze();
delete [] pointIDs;
delete [] voxelShown;
delete [] pointNeeded;
return visualDataSet;
visualDataset->Squeeze();
return visualDataset.GetPointer();

}

Expand Down

0 comments on commit 72faff7

Please sign in to comment.