// VTK Common
#include "vtkMath.h"
#include "vtkImageData.h"
#include "vtkTransform.h"
// VTK Graphics
#include "vtkContourFilter.h"
// VTK IO
#include "vtkVolume16Reader.h"
// VTK Rendering
#include "vtkPolyDataMapper.h"
#include "vtkProperty.h"
#include "vtkActor.h"
#include "vtkCamera.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"

int
main (int argc, char* argv[])
{
  const char* dataRoot = getenv("VTK_DATA_ROOT");
  const char* subdir = "/Data/";
  const char* file1 = "headsq/quarter";
  char* path1 = new char [strlen(dataRoot)+strlen(subdir)+strlen(file1)+1];

  strcpy(path1, dataRoot);
  strcat(path1, subdir);
  strcat(path1, file1);

  //
  // Create the RenderWindow, Renderer and Interactor
  //
  vtkRenderer* ren1 = vtkRenderer::New();
  vtkRenderWindow* renWin = vtkRenderWindow::New();
    {
    renWin->AddRenderer(ren1);
    ren1->Delete();
    }
  vtkRenderWindowInteractor* iren = vtkRenderWindowInteractor::New();
    {
    iren->SetRenderWindow(renWin);
    renWin->Delete();
    }

  int     RESOLUTION    = 64;
  int     START_SLICE   = 50;
  int     END_SLICE     = 60;
  double  PIXEL_SIZE    = 3.2;
  int     centerX       = RESOLUTION / 2;
  int     centerY       = RESOLUTION / 2;
  int     centerZ       = (END_SLICE - START_SLICE) / 2;
  int     endX          = RESOLUTION - 1;
  int     endY          = RESOLUTION - 1;
  int     endZ          = END_SLICE - 1;
  double  origin        = (RESOLUTION / 2.0) * PIXEL_SIZE * -1.0;

  //
  // these transformations permute medical image data to maintain proper 
  // orientation regardless of the acqusition order. After applying these 
  // transforms with vtkTransformFilter, a view up of 0,-1,0 will result 
  // in the body part facing the viewer. 
  // NOTE: Some transformations have a -1 scale factor for one of the 
  //       components. To ensure proper polygon orientation and normal 
  //       direction, you must apply the vtkPolyDataNormals filter.
  //
  // Naming:
  //   si - superior to inferior (top,to,bottom);
  //   is - inferior to superior (bottom,to,top);
  //   ap - anterior to posterior (front,to,back);
  //   pa - posterior to anterior (back,to,front);
  //   lr - left to right
  //   rl - right to left
  //
  double matrixSI[16] = {1,0, 0,0,  0, 0, 1,0,  0,-1,0,0,  0,0,0,1};
  double matrixIS[16] = {1,0, 0,0,  0, 0,-1,0,  0,-1,0,0,  0,0,0,1};
  double matrixLR[16] = {0,0,-1,0,  0,-1, 0,0,  1, 0,0,0,  0,0,0,1};
  double matrixRL[16] = {0,0, 1,0,  0,-1, 0,0,  1, 0,0,0,  0,0,0,1};

  vtkTransform* si = vtkTransform::New();
    {
    si->SetMatrix(matrixSI);
    }
  vtkTransform* is = vtkTransform::New();
    {
    is->SetMatrix(matrixIS);
    }
  vtkTransform* ap = vtkTransform::New();
    {
    ap->Scale(1, -1, 1);
    }
  vtkTransform* pa = vtkTransform::New();
    {
    pa->Scale(1, -1, -1);
    }
  vtkTransform* lr = vtkTransform::New();
    {
    lr->SetMatrix(matrixLR);
    }
  vtkTransform* rl = vtkTransform::New();
    {
    rl->SetMatrix(matrixRL);
    }

  vtkTransform* order[6] = { ap, pa, si, is, lr, rl };

  for (int i=0; i<6; i++)
    {
    vtkActor* actor = vtkActor::New();
      {
      vtkPolyDataMapper* mapper = vtkPolyDataMapper::New();
        {
        vtkContourFilter* iso = vtkContourFilter::New();
          {
          vtkVolume16Reader* reader = vtkVolume16Reader::New();
            {
            reader->SetDataDimensions(RESOLUTION, RESOLUTION);
            reader->SetFilePrefix(path1);
            reader->SetDataSpacing(PIXEL_SIZE, PIXEL_SIZE, 1.5);
            reader->SetDataOrigin(origin, origin, 1.5);
            reader->SetImageRange(START_SLICE, END_SLICE);
            reader->SetTransform(order[i]);
            reader->SetHeaderSize(0);
            reader->SetDataMask(0x7fff);
            reader->SetDataByteOrderToLittleEndian();
            reader->GetOutput()->ReleaseDataFlagOn();
            order[i]->Delete();
            }
          iso->SetInput(reader->GetOutput());
          iso->SetValue(0,550.5);
          iso->ComputeScalarsOff();
          iso->ReleaseDataFlagOn();
          reader->Delete();
          }
        mapper->SetInput(iso->GetOutput());
        mapper->ImmediateModeRenderingOn();
        iso->Delete();
        }
      actor->SetMapper(mapper);
      actor->GetProperty()->SetDiffuseColor( vtkMath::Random(.5,1),
                                             vtkMath::Random(.5,1),
                                             vtkMath::Random(.5,1) );
      mapper->Delete();
      }
    ren1->AddActor(actor);
    actor->Delete();
    }

  renWin->SetSize(300,300);
  ren1->GetActiveCamera()->Azimuth(210);
  ren1->GetActiveCamera()->Elevation(30);
  ren1->GetActiveCamera()->Dolly(1.2);
  ren1->ResetCameraClippingRange();

  ren1->SetBackground(.8,.8,.8);

  renWin->Render();

  // Start the event loop.
  iren->Start();

  iren->Delete();

  delete [] path1;

  return 0;
}