r613: Added DICOM Export
[ctsim.git] / libctsim / ctndicom.cpp
index 8501a0f3a34c6c29f5ade8ce2687026d1a087e82..36dbeb162e7ad295949b45ef566226ff06bf7065 100644 (file)
@@ -9,7 +9,7 @@
 **  This is part of the CTSim program
 **  Copyright (c) 1983-2001 Kevin Rosenberg
 **
-**  $Id: ctndicom.cpp,v 1.6 2001/03/05 21:59:55 kevin Exp $
+**  $Id: ctndicom.cpp,v 1.7 2001/03/07 16:34:47 kevin Exp $
 **
 **  This program is free software; you can redistribute it and/or modify
 **  it under the terms of the GNU General Public License (version 2) as
@@ -93,8 +93,22 @@ DicomImporter::loadImage(unsigned short iNRows, unsigned short iNCols, unsigned
     }
   }
 
-  DCM_ELEMENT elemPixelData = {DCM_PXLPIXELDATA, DCM_OT, "", 1, 0, NULL};
+  char szRescaleSlope[17];
+  char szRescaleIntercept[17];
+  double dRescaleSlope = 1;
+  double dRescaleIntercept = 0;
+  DCM_ELEMENT elemRescaleSlope = {DCM_IMGRESCALESLOPE, DCM_DS, "", 1, strlen(szRescaleSlope), szRescaleSlope};
+  DCM_ELEMENT elemRescaleIntercept = {DCM_IMGRESCALEINTERCEPT, DCM_DS, "", 1, strlen(szRescaleIntercept), szRescaleIntercept};
+  if (DCM_ParseObject (&m_pFile, &elemRescaleSlope, 1, NULL, 0, NULL) == DCM_NORMAL) {
+    if (sscanf (szRescaleSlope, "%lf", &dRescaleSlope) != 1)
+      sys_error (ERR_SEVERE, "Error parsing rescale slope");
+  }
+  if (DCM_ParseObject (&m_pFile, &elemRescaleIntercept, 1, NULL, 0, NULL) == DCM_NORMAL) {
+    if (sscanf (szRescaleIntercept, "%lf", &dRescaleIntercept) != 1)
+      sys_error (ERR_SEVERE, "Error parsing rescale intercept");
+  }
 
+  DCM_ELEMENT elemPixelData = {DCM_PXLPIXELDATA, DCM_OT, "", 1, 0, NULL};
   // Get the actual  pixel data (the only other required element)
   if (DCM_GetElementSize (&m_pFile, elemPixelData.tag, &lRtnLength) != DCM_NORMAL) {
     m_bFail = true;
@@ -150,7 +164,7 @@ DicomImporter::loadImage(unsigned short iNRows, unsigned short iNCols, unsigned
         unsigned char cV1 = pRawPixels[lBase];
         unsigned char cV2 = pRawPixels[lBase+1] & iMask;
         int iV = cV1 + (cV2 << 8);
-        v[ix][iy] = iV / dScale;
+        v[ix][iy] = iV * dRescaleSlope + dRescaleIntercept;
       }
     }
   }
@@ -211,8 +225,191 @@ DicomImporter::loadProjections()
 
 DicomImporter::~DicomImporter()
 {
-  DCM_CloseObject (&m_pFile);
+  if (m_pFile)
+    DCM_CloseObject (&m_pFile);
+}
+
+
+DicomExporter::DicomExporter (ImageFile* pImageFile)
+: m_pImageFile(pImageFile), m_pObject(NULL)
+{
+  DCM_Debug (FALSE);
+  if (! pImageFile) {
+    m_bFail = true;
+    m_strFailMessage = "Initialized DicomExported with NULL imagefile";
+    return;
+  }
+  m_bFail = ! createDicomObject();
+}
+
+DicomExporter::~DicomExporter ()
+{
+  if (m_pObject)
+    DCM_CloseObject (&m_pObject);
+}
+
+bool
+DicomExporter::writeFile (const char* const pszFilename)
+{
+  if (! m_pObject)
+    return false;
+
+  m_strFilename = pszFilename;
+
+  CONDITION cond = DCM_WriteFile (&m_pObject, DCM_ORDERLITTLEENDIAN, pszFilename);
+  if (cond != DCM_NORMAL) {
+    m_bFail = true;
+    m_strFailMessage = "Error writing DICOM file ";
+    m_strFailMessage += pszFilename;
+    return false;
+  }
+
+  return true;
+}
+
+bool
+DicomExporter::createDicomObject()
+{
+  CONDITION cond = DCM_CreateObject (&m_pObject, 0);
+  if (cond != DCM_NORMAL) {
+    m_bFail = true;
+    m_strFailMessage = "Error creating DICOM object";
+    return false;
+  }
+
+  double dMin, dMax;
+  m_pImageFile->getMinMax (dMin, dMax);
+  double dWidth = dMax - dMin;
+  if (dWidth == 0.)
+    dWidth = 1E-7;
+  double dScale = 65535. / dWidth;
+
+  double dRescaleIntercept = -dMin;
+  double dRescaleSlope = 1 / dScale;
+  char szRescaleIntercept[17];
+  char szRescaleSlope[17];
+  snprintf (szRescaleIntercept, sizeof(szRescaleIntercept), "%e", dRescaleIntercept);
+  snprintf (szRescaleSlope, sizeof(szRescaleIntercept), "%e", dRescaleSlope);
+
+  char szModality[] = "CT";
+  char szSOPClassUID[65] = "1.2.840.10008.5.4.1.1.2";
+  char szImgPhotometricInterp[] = "MONOCHROME2";
+  char szPixelSpacing[33] = "0\\0";
+  char szRelImageOrientationPatient[100] = "1\\0\\0\\0\\1\\0";
+  char szRelImagePositionPatient[49] = "0\\0\\0";
+  char szAcqKvp[] = "0";
+  char szRelAcquisitionNumber[] = "1";
+  char szRelImageNumber[] = "1";
+  char szIDSOPInstanceUID[] = "";
+  char szIDManufacturer[] = "CTSim";
+  char szRelPositionRefIndicator[] = "0";
+  char szRelFrameOfReferenceUID[] = "";
+  char szRelSeriesNumber[] = "1";
+  char szIDAccessionNumber[] = "0";
+  char szRelStudyID[] = "1";
+  char szIDReferringPhysician[] = "NONE";
+  char szIDStudyTime[] = "000000.0";
+  char szIDStudyDate[] = "00000000";
+  char szRelStudyInstanceUID[] = "";
+  char szPatSex[] = "O";
+  char szPatBirthdate[] = "0000000";
+  char szPatID[] = "NONE";
+  char szPatName[] = "NONE";
+  char szIDImageType[] = "ORIGINAL";
+  char szIDManufacturerModel[65] = "";
+
+#ifdef VERSION
+  snprintf (szIDManufacturerModel, sizeof(szIDManufacturerModel), "VERSION %s", VERSION);
+#endif
+  snprintf (szPixelSpacing, sizeof(szPixelSpacing), "%e\\%e", m_pImageFile->axisIncrementX(), m_pImageFile->axisIncrementY());
+  double minX, maxX, minY, maxY;
+  if (m_pImageFile->getAxisExtent(minX, maxX, minY, maxY)) {
+    minX += m_pImageFile->axisIncrementX() / 2;
+    minY += m_pImageFile->axisIncrementY() / 2;
+    snprintf(szRelImagePositionPatient, sizeof(szRelImagePositionPatient), "%e\\%e\\0", minX, minY);
+  }
+
+  unsigned short iNRows = m_pImageFile->ny();
+  unsigned short iNCols = m_pImageFile->nx();
+  unsigned short iBitsAllocated = 16;
+  unsigned short iBitsStored = 16;
+  unsigned short iHighBit = 15;
+  unsigned short iPixRep = 0;
+  unsigned short iSamplesPerPixel = 1;
+  DCM_ELEMENT aElemRequired[] = {
+    {DCM_IMGROWS, DCM_US, "", 1, sizeof(iNRows), reinterpret_cast<char*>(&iNRows)},
+    {DCM_IMGCOLUMNS, DCM_US, "", 1, sizeof(iNCols), reinterpret_cast<char*>(&iNCols)},
+    {DCM_IMGBITSALLOCATED, DCM_US, "", 1, sizeof(iBitsAllocated), reinterpret_cast<char*>(&iBitsAllocated)},
+    {DCM_IMGBITSSTORED, DCM_US, "", 1, sizeof(iBitsStored), reinterpret_cast<char*>(&iBitsStored)},
+    {DCM_IMGHIGHBIT, DCM_US, "", 1, sizeof(iHighBit), reinterpret_cast<char*>(&iHighBit)},
+    {DCM_IMGPIXELREPRESENTATION, DCM_US, "", 1, sizeof(iPixRep), reinterpret_cast<char*>(&iPixRep)},
+    {DCM_IMGSAMPLESPERPIXEL, DCM_US, "", 1, sizeof(iSamplesPerPixel), reinterpret_cast<char*>(&iSamplesPerPixel)},
+    {DCM_IMGRESCALESLOPE, DCM_DS, "", 1, strlen(szRescaleSlope), szRescaleSlope},
+    {DCM_IMGRESCALEINTERCEPT, DCM_DS, "", 1, strlen(szRescaleIntercept), szRescaleIntercept},
+    {DCM_IMGPHOTOMETRICINTERP, DCM_CS, "", 1, strlen(szImgPhotometricInterp), szImgPhotometricInterp},
+    {DCM_IMGPIXELSPACING, DCM_DS, "", 1, strlen(szPixelSpacing), szPixelSpacing},
+    {DCM_RELIMAGEORIENTATIONPATIENT, DCM_DS, "", 1, strlen(szRelImageOrientationPatient), szRelImageOrientationPatient},
+    {DCM_RELIMAGEPOSITIONPATIENT, DCM_DS, "", 1, strlen(szRelImagePositionPatient), szRelImagePositionPatient},
+    {DCM_ACQKVP, DCM_DS, "", 1, strlen(szAcqKvp), szAcqKvp},
+    {DCM_RELACQUISITIONNUMBER, DCM_IS, "", 1, strlen(szRelAcquisitionNumber), szRelAcquisitionNumber},
+    {DCM_ACQSLICETHICKNESS, DCM_DS, "", 1, strlen(szRelAcquisitionNumber), szRelAcquisitionNumber},
+    {DCM_RELIMAGENUMBER, DCM_IS, "", 1, strlen(szRelImageNumber), szRelImageNumber},
+    {DCM_IDSOPINSTANCEUID, DCM_UI, "", 1, strlen(szIDSOPInstanceUID), szIDSOPInstanceUID},
+    {DCM_IDMANUFACTURER, DCM_LO, "", 1, strlen(szIDManufacturer), szIDManufacturer},
+    {DCM_RELPOSITIONREFINDICATOR, DCM_LO, "", 1, strlen(szRelPositionRefIndicator), szRelPositionRefIndicator},
+    {DCM_RELFRAMEOFREFERENCEUID, DCM_UI, "", 1, strlen(szRelFrameOfReferenceUID), szRelFrameOfReferenceUID},
+    {DCM_RELSERIESNUMBER, DCM_IS, "", 1, strlen(szRelSeriesNumber), szRelSeriesNumber},
+    {DCM_RELSERIESINSTANCEUID, DCM_UI, "", 1, strlen(szIDAccessionNumber), szIDAccessionNumber},
+    {DCM_IDACCESSIONNUMBER, DCM_SH, "", 1, strlen(szIDAccessionNumber), szIDAccessionNumber},
+    {DCM_RELSTUDYID, DCM_SH, "", 1, strlen(szRelStudyID), szRelStudyID},
+    {DCM_IDREFERRINGPHYSICIAN, DCM_PN, "", 1, strlen(szIDReferringPhysician), szIDReferringPhysician},
+    {DCM_IDSTUDYTIME, DCM_TM, "", 1, strlen(szIDStudyTime), szIDStudyTime},
+    {DCM_IDSTUDYDATE, DCM_DA, "", 1, strlen(szIDStudyDate), szIDStudyDate},
+    {DCM_RELSTUDYINSTANCEUID, DCM_UI, "", 1, strlen(szRelStudyInstanceUID), szRelStudyInstanceUID},
+    {DCM_PATSEX, DCM_CS, "", 1, strlen(szPatSex), szPatSex},
+    {DCM_PATBIRTHDATE, DCM_DA, "", 1, strlen(szPatBirthdate), szPatBirthdate},
+    {DCM_PATID, DCM_LO, "", 1, strlen(szPatID), szPatID},
+    {DCM_PATNAME, DCM_PN, "", 1, strlen(szPatName), szPatName},
+    {DCM_IDIMAGETYPE, DCM_CS, "", 1, strlen(szIDImageType), szIDImageType},
+    {DCM_IDMODALITY, DCM_CS, "", 1, strlen(szModality), szModality},
+    {DCM_IDSOPCLASSUID, DCM_UI, "", 1, strlen(szSOPClassUID), szSOPClassUID},
+    {DCM_IDMANUFACTURERMODEL, DCM_LO, "", 1, strlen(szIDManufacturerModel), szIDManufacturerModel},
+  };
+  int nElemRequired = sizeof (aElemRequired) / sizeof(DCM_ELEMENT);
+
+  int iUpdateCount;
+  cond = DCM_ModifyElements (&m_pObject, aElemRequired, nElemRequired, NULL, 0, &iUpdateCount);
+
+
+  DCM_ELEMENT elemPixelData = {DCM_PXLPIXELDATA, DCM_OT, "", 1, 0, NULL};
+
+  unsigned long lRealLength = 2 * m_pImageFile->nx() * m_pImageFile->ny();
+
+  unsigned char* pRawPixels = new unsigned char [lRealLength];
+  elemPixelData.length = lRealLength;
+  elemPixelData.d.ot = pRawPixels;
+  
+  ImageFileArray v = m_pImageFile->getArray();
+  for (int iy = iNRows - 1; iy >= 0; iy--) {
+    for (int ix = 0; ix < iNCols; ix++) {
+        unsigned long lBase = (iy * iNRows + ix) * 2;
+        unsigned int iValue = nearest<int>(dScale * (v[ix][iy] - dMin));
+        pRawPixels[lBase] = iValue & 0xFF;
+        pRawPixels[lBase+1] = (iValue & 0xFF00) >> 8;
+    }
+  }
+  cond = DCM_ModifyElements (&m_pObject, &elemPixelData, 1, NULL, 0, &iUpdateCount);
+  delete pRawPixels;
+
+  if (cond != DCM_NORMAL || iUpdateCount != 1) {
+    m_bFail = true;
+    m_strFailMessage = "Error modifying pixel data";
+    return false;
+  }
+
+  return true;
 }
 
 #endif // HAVE_CTN_DICOM
 
+