matRad_exportDicomRTStruct

 matRad function to export dicom RT structure set. 
 Class method of matRad_DicomExporter
 
 call
   matRad_DicomExporter.matRad_exportDicomRTStruct()
          

 References
   -

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 Copyright 2015 the matRad development team. 
 
 This file is part of the matRad project. It is subject to the license 
 terms in the LICENSE file found in the top-level directory of this 
 distribution and at https://github.com/e0404/matRad/LICENSES.txt. No part 
 of the matRad project, including this file, may be copied, modified, 
 propagated, or distributed except according to the terms contained in the 
 LICENSE file.

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

0001 function obj = matRad_exportDicomRTStruct(obj)
0002 % matRad function to export dicom RT structure set.
0003 % Class method of matRad_DicomExporter
0004 %
0005 % call
0006 %   matRad_DicomExporter.matRad_exportDicomRTStruct()
0007 %
0008 %
0009 % References
0010 %   -
0011 %
0012 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0013 %
0014 % Copyright 2015 the matRad development team.
0015 %
0016 % This file is part of the matRad project. It is subject to the license
0017 % terms in the LICENSE file found in the top-level directory of this
0018 % distribution and at https://github.com/e0404/matRad/LICENSES.txt. No part
0019 % of the matRad project, including this file, may be copied, modified,
0020 % propagated, or distributed except according to the terms contained in the
0021 % LICENSE file.
0022 %
0023 % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0024 
0025  matRad_cfg = MatRad_Config.instance();
0026 matRad_cfg.dispInfo('Exporting DICOM RTStruct...\n');
0027 
0028 env = matRad_getEnvironment();
0029 isOctave = strcmp(env,'OCTAVE');
0030 
0031 if isOctave
0032     matRad_cfg.dispWarning('RTStruct export currently not supported by matRad running in Octave due to crashing dicomwrite! Skipping...');
0033     return;
0034 end
0035 
0036 %% Metadata
0037 %Class UID
0038 ClassUID = '1.2.840.10008.5.1.4.1.1.481.3'; %RT Structure Set
0039 meta.MediaStorageSOPClassUID = ClassUID;
0040 meta.SOPClassUID             = ClassUID;
0041 %TransferSyntaxUID = '1.2.840.10008.1.2.1'; %Explicit VR Little Endian - correct?
0042 %meta.TransferSyntaxUID = TransferSyntaxUID;
0043 
0044 %Identifiers
0045 meta.SOPInstanceUID             = dicomuid;
0046 meta.MediaStorageSOPInstanceUID = meta.SOPInstanceUID;
0047 meta.SeriesInstanceUID          = dicomuid;
0048 meta.SeriesNumber               = 1;
0049 meta.InstanceNumber             = 1;
0050 
0051 %Remaining Meta Data
0052 meta.Modality = 'RTSTRUCT';
0053 meta.Manufacturer = '';
0054 meta.ReferringPhysicianName = obj.dicomName();
0055 meta.OperatorsName = obj.OperatorsName;
0056 meta.StationName = '';
0057 meta = obj.assignDefaultMetaValue(meta,'ManufacturerModelName','matRad DicomExport');
0058 
0059 meta.PatientName = obj.PatientName;
0060 meta.PatientID = obj.PatientID;
0061 meta.PatientBirthDate = obj.PatientBirthDate;
0062 meta.PatientSex = obj.PatientSex;
0063 
0064 %Name
0065 meta = matRad_DicomExporter.assignDefaultMetaValue(meta,'StructureSetLabel','matRad_cst');
0066 meta = matRad_DicomExporter.assignDefaultMetaValue(meta,'StructureSetName','matRad exported cst');
0067 
0068 
0069 
0070 %ID of the Study
0071 meta.StudyInstanceUID = obj.StudyInstanceUID;
0072 meta.StudyID          = obj.StudyID; 
0073 
0074 %Dates & Times
0075 currDate = now;
0076 currDateStr = datestr(currDate,'yyyymmdd');
0077 currTimeStr = datestr(currDate,'HHMMSS');
0078 meta.InstanceCreationDate = currDateStr;
0079 meta.InstanceCreationTime = currTimeStr;
0080 meta.StudyDate = obj.StudyDate;
0081 meta.StudyTime = obj.StudyTime;
0082 meta = matRad_DicomExporter.assignDefaultMetaValue(meta,'StructureSetDate',currDateStr);
0083 meta = matRad_DicomExporter.assignDefaultMetaValue(meta,'StructureSetTime',currTimeStr);
0084 
0085 
0086 %Remaining stuff
0087 meta.AccessionNumber = '';
0088 
0089 %meta.PatientBirthDate = '';
0090 %meta.PatientSex = 'O';
0091 %meta.SoftwareVersion = '';
0092 
0093 ct = obj.ct;
0094 
0095 %Create X Y Z vectors if not present
0096 if ~any(isfield(ct,{'x','y','z'}))
0097     %positionOffset = transpose(ct.cubeDim ./ 2);
0098     positionOffset = ct.cubeDim ./ 2;
0099     ct.x = ct.resolution.x*[0:ct.cubeDim(2)-1] - positionOffset(2);
0100     ct.y = ct.resolution.y*[0:ct.cubeDim(1)-1] - positionOffset(1);
0101     ct.z = ct.resolution.z*[0:ct.cubeDim(3)-1] - positionOffset(3);
0102 end
0103 
0104 
0105 %Since we are exporting HU directly --> no rescaling in any case
0106 %meta.SliceThickness = ct.resolution.z;
0107 %meta.PixelSpacing = [ct.resolution.y; ct.resolution.x];
0108 %meta.ImageOrientationPatient = [1;0;0;0;1;0]; %lps
0109 
0110 %meta.RescaleSlope = 1;
0111 %meta.RescaleIntercept = 0;
0112 
0113 obj.cst = matRad_computeVoiContoursWrapper(obj.cst,ct);
0114 
0115 for i = 1:size(obj.cst,1)
0116    
0117     fprintf('Processinging VOI ''%s''...',obj.cst{i,2});
0118    %Select contours in axial slices
0119     contours        = obj.cst{i,7}(:,3);
0120     contourSliceIx  = find(~cellfun(@isempty,contours));    
0121     contourSlicePos = ct.z(contourSliceIx);
0122     contours        = contours(contourSliceIx);
0123         
0124     % consider multiple contours of the same object
0125     for slice = 1:numel(contours)
0126        
0127        tmpContour = [];
0128        lower         = 1; % lower marks the beginning of a section
0129        cnt           = 0;
0130        
0131        while lower-1 ~= size(contours{slice,1},2)
0132          
0133           steps = contours{slice,1}(2,lower); % number of elements of current line section
0134           
0135           tmpContour(:,lower-cnt:lower-cnt+steps-1) = contours{slice,1}(:,lower+1:lower+steps);
0136           
0137           lower = lower+steps+1;
0138           cnt   = cnt + 1;
0139        end
0140        
0141        contours{slice,1} = tmpContour;
0142        
0143     end
0144     
0145     contours = cellfun(@(c) [ct.resolution.x; ct.resolution.y].* c + [ct.x(1) - ct.resolution.x;...
0146                              ct.y(1) - ct.resolution.y],contours,'UniformOutput',false);
0147       
0148     
0149     contours = cellfun(@(c,pos) addSlicePos(c,pos),contours,num2cell(contourSlicePos)','UniformOutput',false);
0150     %contours = cellfun(@transpose,contours,'UniformOutput',false);
0151     
0152     %Structure Definition
0153     ROISequenceItem.ROINumber = i;    
0154     ROISequenceItem.ReferencedFrameOfReferenceUID = obj.FrameOfReferenceUID;
0155     ROISequenceItem.ROIName = obj.cst{i,2};
0156     ROISequenceItem.ROIGenerationAlgorithm = '';
0157     
0158     meta.StructureSetROISequence.(['Item_' num2str(i)]) = ROISequenceItem;
0159     
0160     %Contour Sequence
0161     if ~isOctave
0162         ROIContourSequenceItem.ROIDisplayColor = int32(round(255 * obj.cst{i,5}.visibleColor));
0163     end
0164     
0165     %Now create Contour subitems
0166     ContourSequence = struct;
0167     for c = 1:numel(contours)
0168         ContourItem = struct;
0169         
0170         %First store Image Sequence
0171         currCtSliceMeta = obj.ctSliceMetas(contourSliceIx(c));
0172         %currCtSliceMeta = matRad_DicomExporter.assignDefaultMetaValue(currCtSliceMeta,'SOPClassUID','1.2.840.10008.5.1.4.1.1.2');
0173         
0174         %Not sure about this
0175         ContourItem.ContourImageSequence.Item_1.ReferencedSOPClassUID = currCtSliceMeta.SOPClassUID; %We are referencing the CT
0176         ContourItem.ContourImageSequence.Item_1.ReferencedSOPInstanceUID = currCtSliceMeta.SOPInstanceUID; %TODO: ID of the respective ctSlice
0177         
0178         %Now remaing data
0179         ContourItem.ContourGeometricType = 'CLOSED_PLANAR';
0180         ContourItem.NumberOfContourPoints = size(contours{c},2);
0181         ContourItem.ContourData = contours{c}(:);
0182         
0183         
0184         %Now store Item
0185         ContourSequence.(['Item_' num2str(c)]) = ContourItem;
0186     end
0187     ROIContourSequenceItem.ContourSequence = ContourSequence;
0188     
0189     ROIContourSequenceItem.ReferencedROINumber = i;
0190     
0191     %Store to meta data
0192     meta.ROIContourSequence.(['Item_' num2str(i)]) = ROIContourSequenceItem;
0193     
0194     %RTROI Observation Sequence
0195     RTROIObservationsSequenceItem.ObservationNumber = i;
0196     RTROIObservationsSequenceItem.ReferencedROINumber = i;
0197     RTROIObservationsSequenceItem.ROIObservationLabel = obj.cst{i,2};
0198         
0199     
0200     if strcmp(obj.cst{i,3},'TARGET')
0201         if ~isempty(regexpi(obj.cst{i,2},['(' strjoin(obj.targetPtvDict) ')']))
0202             RTROIObservationsSequenceItem.RTROIInterpretedType = 'PTV';
0203             fprintf('identified target type as PTV...');
0204         elseif ~isempty(regexpi(obj.cst{i,2},['(' strjoin(obj.targetGtvDict) ')']))
0205             RTROIObservationsSequenceItem.RTROIInterpretedType = 'GTV';
0206             fprintf('identified target type as GTV...');
0207         elseif ~isempty(regexpi(obj.cst{i,2},['(' strjoin(obj.targetGtvDict) ')']))
0208             RTROIObservationsSequenceItem.RTROIInterpretedType = 'CTV';
0209             fprintf('identified target type as CTV...');
0210         else
0211             RTROIObservationsSequenceItem.RTROIInterpretedType = 'CTV';
0212             fprintf('Defaulting target type to CTV...');
0213         end
0214     else
0215         if ~isempty(regexpi(obj.cst{i,2},['(' strjoin(obj.externalContourDict) ')']))
0216             fprintf('automatically identified as External Contour...');
0217             RTROIObservationsSequenceItem.RTROIInterpretedType = 'EXTERNAL';
0218         else
0219             RTROIObservationsSequenceItem.RTROIInterpretedType = 'AVOIDANCE';
0220         end
0221     end
0222     
0223     RTROIObservationsSequenceItem.ROIInterpreter = obj.dicomName();
0224 
0225     meta.RTROIObservationsSequence.(['Item_' num2str(i)]) = RTROIObservationsSequenceItem;
0226     
0227     fprintf('Done!\n');
0228     %matRad_progress(i,size(obj.cst,1));
0229 end
0230 
0231 
0232 
0233 for i = 1:numel(obj.ctSliceMetas)
0234     ImageSequenceItem.ReferencedSOPClassUID = obj.ctSliceMetas(i).SOPClassUID;
0235     ImageSequenceItem.ReferencedSOPInstanceUID = obj.ctSliceMetas(i).SOPInstanceUID;
0236     ContourImageSequence.(['Item_' num2str(i)]) = ImageSequenceItem;
0237 end
0238 
0239 RTReferencedSeriesSequenceItem.SeriesInstanceUID = obj.ctSliceMetas(1).SeriesInstanceUID;
0240 RTReferencedSeriesSequenceItem.ContourImageSequence = ContourImageSequence;
0241 
0242 RTReferencedStudySequenceItem.ReferencedSOPClassUID = '1.2.840.10008.3.1.2.3.2'; %Apparently this class UID is deprecated in DICOM standard - what to use instead?
0243 RTReferencedStudySequenceItem.ReferencedSOPInstanceUID = obj.StudyInstanceUID;
0244 RTReferencedStudySequenceItem.RTReferencedSeriesSequence.Item_1 = RTReferencedSeriesSequenceItem;
0245 
0246 meta.ReferencedFrameOfReferenceSequence.Item_1.FrameOfReferenceUID = obj.FrameOfReferenceUID;
0247 meta.ReferencedFrameOfReferenceSequence.Item_1.RTReferencedStudySequence.Item_1 = RTReferencedStudySequenceItem;
0248 
0249 filename = 'RTstruct.dcm';
0250 filepath = obj.dicomDir;
0251 filename = fullfile(filepath,filename);
0252 
0253 
0254 if isOctave
0255     dicomwrite(int16(zeros(2)),filename,meta);
0256 else
0257     obj.rtssExportStatus = dicomwrite([],filename,meta,'CreateMode','copy');%,'TransferSyntax',TransferSyntaxUID);
0258 end
0259 
0260 obj.rtssMeta = meta;
0261 end 
0262 
0263 function c = addSlicePos(c,slicePos)
0264     nPoints  = size(c,2);
0265     slicePos = slicePos * ones(1,nPoints);
0266     c(3,:)   = slicePos;
0267 end