OpenCV 3.0 CLAHE with 16bit images

edit

Hello, I'm trying to use the openCV CLAHE implementation on 16 bits medical images and the results are not what I expected, it looks like the algorithm overlflows, resulting in a very dark image. If I convert the original image to CV_8U using cv::Mat:: convertTo() with an alpha value = 1.0 / (1 << (BitsStored - 8)), this is because images are really 10 or 12 bit unsigned int images, the result is correct and what I expect, see the examples below.

Original image

8 bits CLAHE

16 bits CLAHE

The code is very simpe

int main(int argc, char *argv[])  
{
   QApplication a(argc, argv);

   // Load DICOM Image
   DcmObject raw;
   if( !raw.load("rawimage.dcm") )
       qCritical(qPrintable(QString("Error loading rawimage.dcm: %1").arg(raw.errorString())));

   // Convert to cv::Mat, internally can covert to 8 bit if wanted
   cv::Mat mat = dcmImageToCvMat(DcmImage(raw));
   cv::imshow("GRAYSCALE", mat);

   cv::Ptr<cv::CLAHE> clahe = cv::createCLAHE(3.0, cv::Size(8,8));
   cv::Mat dst;
   clahe->apply(mat, dst);
   cv::imshow("CLAHE", dst);

   cv::waitKey();
}

Below is the code I use to convert from DICOM to cv::Mat:, it uses a dicom toolkit written by me:

cv::Mat dcmImageToCvMat(const DcmImage &img, bool convertTo8Bits)
{
   cv::Mat mat;
   switch( img.pixelRepresentation() )
   {
   case DcmImage::PixelUInt8:
       mat = cv::Mat(img.rows(), img.columns(), CV_8U);
       break;
   case DcmImage::PixelInt8:
       mat = cv::Mat(img.rows(), img.columns(), CV_8S);
       break;
   case DcmImage::PixelUInt16:
       mat = cv::Mat(img.rows(), img.columns(), CV_16U);
       break;
   case DcmImage::PixelInt16:
       mat = cv::Mat(img.rows(), img.columns(), CV_16S);
       break;
   case DcmImage::PixelInt32:
   case DcmImage::PixelUInt32:
       mat =cv::Mat(img.rows(), img.columns(),  CV_32S);
       break;
   }

   // Get RAW Pixel data
   img.raw(mat.data, 0, 0, img.columns(), img.rows(), mat.step);

   // Invert if Photometric interpretation is MONOCHROME1 (white is 0)
   if( img.photometricInterpretation()==DcmImage::PmiMonochrome1 )
   {
       cv::Mat max(mat.size(), mat.type());
       max.setTo(img.maxPossibleValue());
       mat = max - mat;
   }

   // Convert to CV_16U with correct scaling of values
   switch( img.pixelRepresentation() )
   {
   case DcmImage::PixelInt8:
       mat.convertTo(mat, CV_8U, 1.0, -img.minPossibleValue());
       break;
   case DcmImage::PixelInt16:
       mat.convertTo(mat, CV_16U, 1.0, -img.minPossibleValue());
       break;
   case DcmImage::PixelInt32:
   case DcmImage::PixelUInt32:
       mat.convertTo(mat, CV_16U, 1.0 / double(1 << (img.bitsStored() - 16)), -img.minPossibleValue());
       break;
   }

   // Converts to 8 bits if wanted
   if( convertTo8Bits )
   {
       cv::Mat mat8;
       float scale = 1 << (img.bitsStored() - 8);
       float shift = img.isSigned() ? 127.0 : 0.0;
       cv::convertScaleAbs(mat, mat8, 1.0 / scale, shift);
       return mat8;
   }
   else
       return mat;
}

Thank you very much,

Gianluca Ghelli