ANIMA  4.0
animaShapesWriter.cxx
Go to the documentation of this file.
1 #include <animaShapesWriter.h>
2 #include <itkMacro.h>
3 
4 #include <vtkPolyDataWriter.h>
5 #include <vtkXMLPolyDataWriter.h>
6 #include <vtksys/SystemTools.hxx>
7 
8 #include <vtkPointData.h>
9 
10 #include <fstream>
11 #include <algorithm>
12 
13 namespace anima {
14 
16 {
17  std::string extensionName = m_FileName.substr(m_FileName.find_last_of('.') + 1);
18  if (extensionName == "vtk")
19  this->WriteFileAsVTKAscii();
20  else if ((extensionName == "vtp")||(extensionName == ""))
21  {
22  if (extensionName == "")
23  m_FileName += ".vtp";
24  this->WriteFileAsVTKXML();
25  }
26  else if (extensionName == "fds")
28  else if (extensionName == "csv")
29  this->WriteFileAsCSV();
30  else
31  throw itk::ExceptionObject(__FILE__, __LINE__,"Unsupported shapes extension.",ITK_LOCATION);
32 }
33 
35 {
36  vtkSmartPointer <vtkPolyDataWriter> vtkWriter = vtkPolyDataWriter::New();
37  vtkWriter->SetInputData(m_InputData);
38  vtkWriter->SetFileName(m_FileName.c_str());
39  vtkWriter->Update();
40 }
41 
43 {
44  vtkSmartPointer <vtkXMLPolyDataWriter> vtkWriter = vtkXMLPolyDataWriter::New();
45  vtkWriter->SetInputData(m_InputData);
46  vtkWriter->SetFileName(m_FileName.c_str());
47  vtkWriter->SetDataModeToBinary();
48  vtkWriter->EncodeAppendedDataOff();
49  vtkWriter->SetCompressorTypeToZLib();
50  vtkWriter->Update();
51 }
52 
54 {
55  std::replace(m_FileName.begin(),m_FileName.end(),'\\','/');
56 
57  std::string baseName;
58  std::size_t lastDotPos = m_FileName.find_last_of('.');
59  baseName.append(m_FileName.begin(),m_FileName.begin() + lastDotPos);
60 
61  std::string noPathName = baseName;
62  std::size_t lastSlashPos = baseName.find_last_of("/");
63 
64  if (lastSlashPos != std::string::npos)
65  {
66  noPathName.clear();
67  noPathName.append(baseName.begin() + lastSlashPos + 1,baseName.end());
68  }
69 
70  vtksys::SystemTools::MakeDirectory(baseName.c_str());
71 
72  std::string vtkFileName = noPathName + "/";
73  vtkFileName += noPathName;
74  vtkFileName += "_0.vtp";
75 
76  std::string vtkWriteFileName = baseName + "/";
77  vtkWriteFileName += noPathName + "_0.vtp";
78 
79  vtkSmartPointer <vtkXMLPolyDataWriter> vtkWriter = vtkXMLPolyDataWriter::New();
80  vtkWriter->SetInputData(m_InputData);
81  vtkWriter->SetFileName(vtkWriteFileName.c_str());
82  vtkWriter->SetDataModeToBinary();
83  vtkWriter->EncodeAppendedDataOff();
84  vtkWriter->SetCompressorTypeToZLib();
85  vtkWriter->Update();
86 
87  std::ofstream outputHeaderFile(m_FileName.c_str());
88  outputHeaderFile << "<?xml version=\"1.0\"?>" << std::endl;
89  outputHeaderFile << "<VTKFile type=\"vtkFiberDataSet\" version=\"1.0\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">" << std::endl;
90  outputHeaderFile << "<vtkFiberDataSet>" << std::endl;
91  outputHeaderFile << "\t<Fibers index=\"0\" file=\"" << vtkFileName << "\">" << std::endl;
92  outputHeaderFile << "\t</Fibers>" << std::endl;
93  outputHeaderFile << "</vtkFiberDataSet>" << std::endl;
94  outputHeaderFile << "</VTKFile>" << std::endl;
95 
96  outputHeaderFile.close();
97 }
98 
100 {
101  vtkSmartPointer<vtkPointData> inputData = m_InputData->GetPointData();
102 
103  // Initialize output file.
104  std::ofstream outputFile;
105  outputFile.open(m_FileName.c_str(), std::ios_base::out);
106  outputFile.precision(std::numeric_limits<long double>::digits10);
107 
108  if (outputFile.bad())
109  throw itk::ExceptionObject(__FILE__, __LINE__, "The output file could not be opened", ITK_LOCATION);
110 
111  // Export data to outputFile.
112  typedef std::vector<int> IndexVectorType;
113  unsigned int numArrays = inputData->GetNumberOfArrays();
114  IndexVectorType arraySizes(numArrays, 0);
115 
116  outputFile << "X,Y,Z,PointId,StreamlineId";
117 
118  for (unsigned int i = 0;i < numArrays;++i)
119  {
120  int arraySize = inputData->GetArray(i)->GetNumberOfComponents();
121  arraySizes[i] = arraySize;
122 
123  if (arraySize == 1)
124  {
125  outputFile << "," << inputData->GetArrayName(i);
126  continue;
127  }
128 
129  for (unsigned int j = 0;j < arraySize;++j)
130  outputFile << "," << inputData->GetArrayName(i) << "#" << j;
131  }
132 
133  //-------------------------------
134  // Setting up streamline geometry
135  //-------------------------------
136 
137  unsigned int numberOfPoints = m_InputData->GetNumberOfPoints();
138  unsigned int numberOfStreamlines = m_InputData->GetNumberOfLines();
139  std::cout << "Number of data points: " << numberOfPoints << std::endl;
140  std::cout << "Number of streamlines: " << numberOfStreamlines << std::endl;
141 
142  // Extract streamline information by point
143  IndexVectorType pointId(numberOfPoints, -1);
144  IndexVectorType streamlineId(numberOfPoints, -1);
145  m_InputData->GetLines()->InitTraversal();
146  vtkSmartPointer<vtkIdList> idList = vtkSmartPointer<vtkIdList>::New();
147  for (unsigned int i = 0;i < numberOfStreamlines;++i)
148  {
149  m_InputData->GetLines()->GetNextCell(idList);
150 
151  unsigned int streamlineSize = idList->GetNumberOfIds();
152 
153  if (streamlineSize == 1)
154  continue;
155 
156  for (unsigned int j = 0;j < streamlineSize;++j)
157  {
158  unsigned int pid = idList->GetId(j);
159  streamlineId[pid] = i+1;
160  pointId[pid] = j+1;
161  }
162  }
163 
164  // Writing table content
165  for (unsigned int i = 0;i < numberOfPoints;++i)
166  {
167  if (numberOfStreamlines != 0)
168  if (streamlineId[i] == -1)
169  continue;
170 
171  outputFile << std::endl;
172 
173  // 1. Write point 3D coordinates
174  double p[3];
175  m_InputData->GetPoint(i, p);
176 
177  for (unsigned int j = 0;j < 3;++j)
178  outputFile << p[j] << ",";
179 
180  // 2. Write streamline index data
181  outputFile << pointId[i] << "," << streamlineId[i];
182 
183  // 3. Write array values if any
184  for (unsigned int k = 0;k < numArrays;++k)
185  for (unsigned int j = 0;j < arraySizes[k];++j)
186  outputFile << "," << inputData->GetArray(k)->GetComponent(i, j);
187  }
188 
189  outputFile << std::endl;
190  outputFile.close();
191 }
192 
193 } // end namespace anima