Commit 27eb09b1 authored by Haocheng LIU's avatar Haocheng LIU
Browse files

VTK Bug #16736: Add large data support&test for vtkMPIImageReader

The "ReadSlice" method in 'vtkMPIImageReader' can only read raw data if
the dataset is within 2^31 elements. Using the same technique in Merge
Request !1510, we can use a while loop to extend its ability to read the large
data slice by slice. A test file is also added for vtkMPIImageReader
class.
parent 8bc5f126
......@@ -3,5 +3,6 @@ include(vtkMPI)
vtk_add_test_mpi(${vtk-module}CxxTests-MPI tests
TESTING_DATA
ParallelIso.cxx
ParallelIso2.cxx
)
vtk_test_mpi_executable(${vtk-module}CxxTests-MPI tests)
......@@ -44,6 +44,8 @@
#include "vtkDebugLeaks.h"
namespace {
static const float ISO_START=4250.0;
static const float ISO_STEP=-1250.0;
static const int ISO_NUM=3;
......@@ -231,6 +233,8 @@ void MyMain( vtkMultiProcessController *controller, void *arg )
}
} // end anon namespace
int ParallelIso( int argc, char* argv[] )
{
// This is here to avoid false leak messages from vtkDebugLeaks when
......
/*=========================================================================
Program: Visualization Toolkit
Module: ParallelIso2.cxx
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
This software is distributed WITHOUT ANY WARRANTY; without even
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notice for more information.
=========================================================================*/
// This example demonstrates the use of data parallelism in VTK. The
// pipeline ( vtkMPIImageReader -> vtkContourFilter -> vtkElevationFilter )
// is created in parallel and each process is assigned 1 piece to process.
// All satellite processes send the result to the first process which
// collects and renders them.
#include <mpi.h>
#include "vtkActor.h"
#include "vtkAppendPolyData.h"
#include "vtkCamera.h"
#include "vtkConeSource.h"
#include "vtkContourFilter.h"
#include "vtkDataSet.h"
#include "vtkElevationFilter.h"
#include "vtkMath.h"
#include "vtkMPIController.h"
#include "vtkMPIImageReader.h"
#include "vtkPolyData.h"
#include "vtkPolyDataMapper.h"
#include "vtkTestUtilities.h"
#include "vtkRegressionTestImage.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkRenderer.h"
#include "vtkWindowToImageFilter.h"
#include "vtkImageData.h"
#include "vtkStreamingDemandDrivenPipeline.h"
#include "vtkInformation.h"
#include "vtkDebugLeaks.h"
namespace {
static const float ISO_START=4250.0;
static const float ISO_STEP=-1250.0;
static const int ISO_NUM=3;
// Just pick a tag which is available
static const int ISO_VALUE_RMI_TAG=300;
static const int ISO_OUTPUT_TAG=301;
struct ParallelIsoArgs_tmp
{
int* retVal;
int argc;
char** argv;
};
struct ParallelIsoRMIArgs_tmp
{
vtkContourFilter* ContourFilter;
vtkMultiProcessController* Controller;
vtkElevationFilter* Elevation;
};
// call back to set the iso surface value.
void SetIsoValueRMI(void *localArg, void* vtkNotUsed(remoteArg),
int vtkNotUsed(remoteArgLen), int vtkNotUsed(id))
{
ParallelIsoRMIArgs_tmp* args = (ParallelIsoRMIArgs_tmp*)localArg;
float val;
vtkMultiProcessController* contrl = args->Controller;
int myid = contrl->GetLocalProcessId();
int numProcs = contrl->GetNumberOfProcesses();
vtkContourFilter *iso = args->ContourFilter;
val = iso->GetValue(0);
iso->SetValue(0, val + ISO_STEP);
args->Elevation->UpdatePiece(myid, numProcs, 0);
contrl->Send(args->Elevation->GetOutput(), 0, ISO_OUTPUT_TAG);
}
// This will be called by all processes
void MyMain( vtkMultiProcessController *controller, void *arg )
{
vtkMPIImageReader *reader;
vtkContourFilter *iso;
vtkElevationFilter *elev;
int myid, numProcs;
float val;
ParallelIsoArgs_tmp* args = reinterpret_cast<ParallelIsoArgs_tmp*>(arg);
// Obtain the id of the running process and the total
// number of processes
myid = controller->GetLocalProcessId();
numProcs = controller->GetNumberOfProcesses();
// Create the reader, the data file name might have
// to be changed depending on where the data files are.
char* fname = vtkTestUtilities::ExpandDataFileName(args->argc, args->argv,
"Data/headsq/quarter");
reader = vtkMPIImageReader::New();
reader->SetDataByteOrderToLittleEndian();
reader->SetDataExtent(0, 63, 0, 63, 1, 93);
reader->SetFilePrefix(fname);
reader->SetDataSpacing(3.2, 3.2, 1.5);
delete[] fname;
// Iso-surface.
iso = vtkContourFilter::New();
iso->SetInputConnection(reader->GetOutputPort());
iso->SetValue(0, ISO_START);
iso->ComputeScalarsOff();
iso->ComputeGradientsOff();
// Compute a different color for each process.
elev = vtkElevationFilter::New();
elev->SetInputConnection(iso->GetOutputPort());
val = (myid+1) / static_cast<float>(numProcs);
elev->SetScalarRange(val, val+0.001);
// Make sure all processes update at the same time.
elev->UpdatePiece(myid, numProcs, 0);
if (myid != 0)
{
// If I am not the root process
ParallelIsoRMIArgs_tmp args2;
args2.ContourFilter = iso;
args2.Controller = controller;
args2.Elevation = elev;
// Last, set up a RMI call back to change the iso surface value.
// This is done so that the root process can let this process
// know that it wants the contour value to change.
controller->AddRMI(SetIsoValueRMI, (void *)&args2, ISO_VALUE_RMI_TAG);
controller->ProcessRMIs();
}
else
{
// Create the rendering part of the pipeline
vtkAppendPolyData *app = vtkAppendPolyData::New();
vtkRenderer *ren = vtkRenderer::New();
vtkRenderWindow *renWindow = vtkRenderWindow::New();
vtkRenderWindowInteractor *iren = vtkRenderWindowInteractor::New();
vtkPolyDataMapper *mapper = vtkPolyDataMapper::New();
vtkActor *actor = vtkActor::New();
vtkCamera *cam = vtkCamera::New();
renWindow->AddRenderer(ren);
iren->SetRenderWindow(renWindow);
ren->SetBackground(0.9, 0.9, 0.9);
renWindow->SetSize( 400, 400);
mapper->SetInputConnection(app->GetOutputPort());
actor->SetMapper(mapper);
ren->AddActor(actor);
cam->SetFocalPoint(100, 100, 65);
cam->SetPosition(100, 450, 65);
cam->SetViewUp(0, 0, -1);
cam->SetViewAngle(30);
cam->SetClippingRange(177.0, 536.0);
ren->SetActiveCamera(cam);
// loop through some iso surface values.
for (int j = 0; j < ISO_NUM; ++j)
{
for (int i = 1; i < numProcs; ++i)
{
// trigger the RMI to change the iso surface value.
controller->TriggerRMI(i, ISO_VALUE_RMI_TAG);
}
// set the local value
iso->SetValue(0, iso->GetValue(0) + ISO_STEP);
elev->UpdatePiece(myid, numProcs, 0);
for (int i = 1; i < numProcs; ++i)
{
vtkPolyData* pd = vtkPolyData::New();
controller->Receive(pd, i, ISO_OUTPUT_TAG);
if (j == ISO_NUM - 1)
{
app->AddInputData(pd);
}
pd->Delete();
}
}
// Tell the other processors to stop processing RMIs.
for (int i = 1; i < numProcs; ++i)
{
controller->TriggerRMI(i, vtkMultiProcessController::BREAK_RMI_TAG);
}
vtkPolyData* outputCopy = vtkPolyData::New();
outputCopy->ShallowCopy(elev->GetOutput());
app->AddInputData(outputCopy);
outputCopy->Delete();
app->Update();
outputCopy->RemoveGhostCells();
renWindow->Render();
*(args->retVal) =
vtkRegressionTester::Test(args->argc, args->argv, renWindow, 10);
if ( *(args->retVal) == vtkRegressionTester::DO_INTERACTOR)
{
iren->Start();
}
// Clean up
app->Delete();
ren->Delete();
renWindow->Delete();
iren->Delete();
mapper->Delete();
actor->Delete();
cam->Delete();
}
// clean up objects in all processes.
reader->Delete();
iso->Delete();
elev->Delete();
}
} // end anon namespace
int ParallelIso2( int argc, char* argv[] )
{
// This is here to avoid false leak messages from vtkDebugLeaks when
// using mpich. It appears that the root process which spawns all the
// main processes waits in MPI_Init() and calls exit() when
// the others are done, causing apparent memory leaks for any objects
// created before MPI_Init().
MPI_Init(&argc, &argv);
// Note that this will create a vtkMPIController if MPI
// is configured, vtkThreadedController otherwise.
vtkMPIController* controller = vtkMPIController::New();
controller->Initialize(&argc, &argv, 1);
// Added for regression test.
// ----------------------------------------------
int retVal = 1;
ParallelIsoArgs_tmp args;
args.retVal = &retVal;
args.argc = argc;
args.argv = argv;
// ----------------------------------------------
controller->SetSingleMethod(MyMain, &args);
controller->SingleMethodExecute();
controller->Finalize();
controller->Delete();
return !retVal;
}
......@@ -52,6 +52,7 @@
#include "vtkPointData.h"
#include "vtkTransform.h"
#include "vtkType.h"
#include <algorithm>
// This macro can be wrapped around MPI function calls to easily report errors.
// Reporting errors is more important with file I/O because, unlike network I/O,
......@@ -295,15 +296,24 @@ void vtkMPIImageReader::ReadSlice(int slice, const int extent[6], void *buffer)
length *= extent[3]-extent[2]+1;
if (this->GetFileDimensionality() == 3) length *= extent[5]-extent[4]+1;
if (length > VTK_INT_MAX)
vtkIdType pos = 0;
while (length > pos)
{
vtkErrorMacro(<< "Cannot read more than " << VTK_INT_MAX
<< " bytes at a time.");
MPI_Status stat;
// we know this will fit in an int because it can't exceed VTK_INT_MAX.
const int remaining = static_cast<int>(std::min(length - pos,
static_cast<vtkIdType>(VTK_INT_MAX)));
MPICall(MPI_File_read(file.Handle, (static_cast<char*>(buffer)) + pos, remaining,
MPI_BYTE, &stat));
int rd = 0;
MPICall(MPI_Get_elements(&stat, MPI_BYTE, &rd));
if (MPI_UNDEFINED == rd)
{
vtkErrorMacro("Error obtaining number of values read in " << remaining <<
"-byte read.");
}
pos += static_cast<vtkIdType>(rd);
}
// Do the read. This is a coordinated parallel operation for efficiency.
MPICall(MPI_File_read_all(file.Handle, buffer, static_cast<int>(length),
MPI_BYTE, MPI_STATUS_IGNORE));
MPICall(MPI_File_close(&file.Handle));
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment