source: trunk/gui/src/RpDicomToVtk.cc @ 4979

Last change on this file since 4979 was 4979, checked in by ldelgass, 6 years ago

Add guards around debugging output in DICOM to VTK converter

File size: 13.0 KB
Line 
1/*
2 * ----------------------------------------------------------------------
3 *  RpDicomToVtk -
4 *
5 * ======================================================================
6 *  AUTHOR:  Leif Delgass, Purdue University
7 *  Copyright (c) 2014  HUBzero Foundation, LLC
8 *
9 *  See the file "license.terms" for information on usage and
10 *  redistribution of this file, and for a DISCLAIMER OF ALL WARRANTIES.
11 * ======================================================================
12 */
13
14#include <vtkSmartPointer.h>
15#include <vtkStringArray.h>
16#include <vtkIntArray.h>
17#include <vtkDataSetWriter.h>
18#include <vtkDirectory.h>
19#ifdef USE_VTK_DICOM_PACKAGE
20#include <vtkDICOMSorter.h>
21#include <vtkDICOMReader.h>
22#include <vtkDICOMMetaData.h>
23#else
24#include <vtkDICOMImageReader.h>
25#endif
26
27#include <stdio.h>
28#include "tcl.h"
29
30// #define RP_DICOM_TRACE
31
32static int
33DicomToVtkCmd(ClientData clientData, Tcl_Interp *interp, int objc,
34              Tcl_Obj *const *objv)
35{
36    if (objc != 3) {
37        Tcl_AppendResult(interp, "wrong # args: should be \"",
38                         Tcl_GetString(objv[0]), " dir|files fileList\"", (char*)NULL);
39        return TCL_ERROR;
40    }
41    bool isDir = false;
42    char *arg = Tcl_GetString(objv[1]);
43    if (arg[0] == 'd' && strcmp(arg, "dir") == 0) {
44        isDir = true;
45    }
46
47    char *dirName = NULL;
48    vtkSmartPointer<vtkStringArray> pathArray = vtkSmartPointer<vtkStringArray>::New();
49    if (isDir) {
50        dirName = Tcl_GetString(objv[2]);
51    } else {
52        // Get path array from Tcl list
53        int pathListc;
54        Tcl_Obj **pathListv = NULL;
55        if (Tcl_ListObjGetElements(interp, objv[2], &pathListc, &pathListv) != TCL_OK) {
56            return TCL_ERROR;
57        }
58        pathArray->SetNumberOfValues(pathListc);
59        for (int i = 0; i < pathListc; i++) {
60            pathArray->SetValue(i, Tcl_GetString(pathListv[i]));
61        }
62    }
63
64#ifdef USE_VTK_DICOM_PACKAGE
65    vtkSmartPointer<vtkDICOMReader> reader = vtkSmartPointer<vtkDICOMReader>::New();
66    vtkSmartPointer<vtkDICOMSorter> sorter = vtkSmartPointer<vtkDICOMSorter>::New();
67    if (isDir) {
68        vtkSmartPointer<vtkDirectory> dir = vtkSmartPointer<vtkDirectory>::New();
69        if (!dir->Open(dirName)) {
70            return TCL_ERROR;
71        }
72
73        int numFiles = dir->GetNumberOfFiles();
74        for (int i = 0; i < numFiles; i++) {
75            if (dir->GetFile(i)[0] == '.') {
76                continue;
77            }
78
79            std::string path(dirName);
80            path += "/";
81            path += dir->GetFile(i);
82
83            if (!dir->FileIsDirectory(path.c_str()) &&
84                reader->CanReadFile(path.c_str())) {
85                pathArray->InsertNextValue(path);
86            } else {
87                fprintf(stderr, "Skipping file %s\n", path.c_str());
88            }
89        }
90    }
91
92    sorter->SetInputFileNames(pathArray);
93    sorter->Update();
94
95    int numStudies = sorter->GetNumberOfStudies();
96    int study = 0;
97    int series = 0;
98
99#ifdef RP_DICOM_TRACE
100    fprintf(stderr, "Num Studies: %d\n", numStudies);
101#endif
102    vtkStringArray *files;
103#if 0
104    for (int i = 0; i < numStudies; i++) {
105        int numSeries = sorter->GetNumberOfSeriesInStudy(i);
106#ifdef RP_DICOM_TRACE
107        fprintf(stderr, "Study %d: %d series\n", i, numSeries);
108#endif
109        int k = sorter->GetFirstSeriesInStudy(i);
110        for (int j = 0; j < numSeries; j++) {
111            files = sorter->GetFileNamesForSeries(k+j);
112        }
113    }
114#else
115    if (study < 0 || study > sorter->GetNumberOfStudies()-1) {
116        return TCL_ERROR;
117    }
118    if (series < 0 || series > sorter->GetNumberOfSeriesInStudy(study)-1) {
119        return TCL_ERROR;
120    }
121    files = sorter->GetFileNamesForSeries(sorter->GetFirstSeriesInStudy(study)+series);
122#endif
123
124    Tcl_Obj *objPtr = Tcl_NewListObj(0, NULL);
125    Tcl_Obj *metaDataObj = Tcl_NewListObj(0, NULL);
126    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("num_studies", -1));
127    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(numStudies));
128    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("study", -1));
129    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(study));
130    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("num_series", -1));
131    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(sorter->GetNumberOfSeriesInStudy(study)));
132    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("series", -1));
133    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(series));
134
135    reader->SetFileNames(files);
136#else
137    Tcl_Obj *objPtr = Tcl_NewListObj(0, NULL);
138    Tcl_Obj *metaDataObj = Tcl_NewListObj(0, NULL);
139
140    vtkSmartPointer<vtkDICOMImageReader> reader = vtkSmartPointer<vtkDICOMImageReader>::New();
141    if (isDir) {
142        reader->SetDirectoryName(dirName);
143    } else {
144        reader->SetFileName(pathArray->GetValue(0));
145    }
146#endif
147    reader->Update();
148
149    Tcl_Obj *eltObj;
150#ifdef USE_VTK_DICOM_PACKAGE
151    vtkStringArray *ids = reader->GetStackIDs();
152#ifdef RP_DICOM_TRACE
153    for (int i = 0; i < ids->GetNumberOfValues(); i++) {
154        fprintf(stderr, "Stack: %s\n", ids->GetValue(i).c_str());
155    }
156#endif
157    vtkIntArray *fidxArray = reader->GetFileIndexArray();
158    vtkDICOMMetaData *md = reader->GetMetaData();
159    int numComp = fidxArray->GetNumberOfComponents();
160#if 0
161    for (int i = 0; i < fidxArray->GetNumberOfTuples(); i++) {
162        fprintf(stderr, "%d:", i);
163        for (int j = 0; j < fidxArray->GetNumberOfComponents(); j++) {
164            int idx = (int)fidxArray->GetComponent(i, j);
165            fprintf(stderr, " %d", idx);
166        }
167        fprintf(stderr, "\n");
168    }
169#endif
170#ifdef RP_DICOM_TRACE
171    fprintf(stderr, "Number of data elements: %d\n", md->GetNumberOfDataElements());
172#endif
173    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("num_files", -1));
174    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(md->GetNumberOfInstances()));
175    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("num_components", -1));
176    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(numComp));
177    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("modality", -1));
178    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(md->GetAttributeValue(DC::Modality).AsString().c_str(), -1));
179    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("patient_name", -1));
180    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(md->GetAttributeValue(DC::PatientName).AsString().c_str(), -1));
181    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("transfer_syntax_uid", -1));
182    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(md->GetAttributeValue(DC::TransferSyntaxUID).AsString().c_str(), -1));
183    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("study_uid", -1));
184    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(md->GetAttributeValue(DC::StudyInstanceUID).AsString().c_str(), -1));
185    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("study_id", -1));
186    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(md->GetAttributeValue(DC::StudyID).AsString().c_str(), -1));
187    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("series_number", -1));
188    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(md->GetAttributeValue(DC::SeriesNumber).AsInt()));
189    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("series_in_study", -1));
190    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(md->GetAttributeValue(DC::SeriesInStudy).AsInt()));
191    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("series_uid", -1));
192    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(md->GetAttributeValue(DC::SeriesInstanceUID).AsString().c_str(), -1));
193    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("bits_allocated", -1));
194    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(md->GetAttributeValue(DC::BitsAllocated).AsInt()));
195    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("bits_stored", -1));
196    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(md->GetAttributeValue(DC::BitsStored).AsInt()));
197    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("pixel_representation", -1));
198    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(md->GetAttributeValue(DC::PixelRepresentation).AsString().c_str(), -1));
199    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("rescale_slope", -1));
200    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewDoubleObj(md->GetAttributeValue(DC::RescaleSlope).AsDouble()));
201    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("rescale_intercept", -1));
202    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewDoubleObj(md->GetAttributeValue(DC::RescaleIntercept).AsDouble()));
203    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("time_dimension", -1));
204    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(reader->GetTimeDimension()));
205    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("time_spacing", -1));
206    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewDoubleObj(reader->GetTimeSpacing()));
207#else
208    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("patient_name", -1));
209    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(reader->GetPatientName(), -1));
210    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("transfer_syntax_uid", -1));
211    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(reader->GetTransferSyntaxUID(), -1));
212    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("study_uid", -1));
213    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(reader->GetStudyUID(), -1));
214    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("study_id", -1));
215    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj(reader->GetStudyID(), -1));
216    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("num_components", -1));
217    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(reader->GetNumberOfComponents()));
218    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("bits_allocated", -1));
219    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewIntObj(reader->GetBitsAllocated()));
220    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("pixel_representation", -1));
221    Tcl_ListObjAppendElement(interp, metaDataObj,
222                             Tcl_NewStringObj((reader->GetPixelRepresentation() ? "signed" : "unsigned"), -1));
223    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("rescale_slope", -1));
224    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewDoubleObj(reader->GetRescaleSlope()));
225    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("rescale_offset", -1));
226    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewDoubleObj(reader->GetRescaleOffset()));
227    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewStringObj("gantry_angle", -1));
228    Tcl_ListObjAppendElement(interp, metaDataObj, Tcl_NewDoubleObj(reader->GetGantryAngle()));
229#if 0
230    fprintf(stderr, "Position: %g %g %g\n",
231            reader->GetImagePositionPatient()[0],
232            reader->GetImagePositionPatient()[1],
233            reader->GetImagePositionPatient()[2]);
234    fprintf(stderr, "Orientation: %g %g %g, %g %g %g\n",
235            reader->GetImageOrientationPatient()[0],
236            reader->GetImageOrientationPatient()[1],
237            reader->GetImageOrientationPatient()[2],
238            reader->GetImageOrientationPatient()[3],
239            reader->GetImageOrientationPatient()[4],
240            reader->GetImageOrientationPatient()[5]);
241#endif
242#endif
243
244    Tcl_ListObjAppendList(interp, objPtr, metaDataObj);
245#ifdef RP_DICOM_TRACE
246    fprintf(stderr, "writing VTK\n");
247#endif
248    vtkSmartPointer<vtkDataSetWriter> writer = vtkSmartPointer<vtkDataSetWriter>::New();
249    writer->SetInputConnection(reader->GetOutputPort());
250    writer->SetFileTypeToBinary();
251    writer->WriteToOutputStringOn();
252    writer->Update();
253#ifdef RP_DICOM_TRACE
254    fprintf(stderr, "writing VTK...done\n");
255#endif
256    Tcl_ListObjAppendElement(interp, objPtr, Tcl_NewStringObj("vtkdata", -1));
257
258    eltObj = Tcl_NewByteArrayObj((unsigned char *)writer->GetBinaryOutputString(),
259                                 writer->GetOutputStringLength());
260    Tcl_ListObjAppendElement(interp, objPtr, eltObj);
261
262    Tcl_SetObjResult(interp, objPtr);
263    return TCL_OK;
264}
265
266extern "C" {
267/*
268 * ------------------------------------------------------------------------
269 *  RpDicomToVtk_Init --
270 *
271 *  Invoked when the Rappture GUI library is being initialized
272 *  to install the "DicomToVtk" command into the interpreter.
273 *
274 *  Returns TCL_OK if successful, or TCL_ERROR (along with an error
275 *  message in the interp) if anything goes wrong.
276 * ------------------------------------------------------------------------
277 */
278int
279RpDicomToVtk_Init(Tcl_Interp *interp)
280{
281    /* install the widget command */
282    Tcl_CreateObjCommand(interp, "Rappture::DicomToVtk", DicomToVtkCmd,
283                         NULL, NULL);
284    return TCL_OK;
285}
286
287}
Note: See TracBrowser for help on using the repository browser.