Additional hunks of the patch from Brent Townshend to improve HDR support.
This needs some additional work to make the annealer deal with various MismatchImagePixelTypes.
- Provide a generic definition of "infinite" cost.
- Make the annealer work efficiently with non-integral cost function values.
1.1 --- a/src/mask.h Thu Jun 11 15:20:52 2009 +0000
1.2 +++ b/src/mask.h Sun Jun 28 02:31:28 2009 -0700
1.3 @@ -107,7 +107,7 @@
1.4 {
1.5 typedef typename EnblendNumericTraits<PixelType>::ImagePixelComponentType PixelComponentType;
1.6 typedef typename EnblendNumericTraits<ResultType>::ImagePixelComponentType ResultPixelComponentType;
1.7 - typedef LinearIntensityTransform<ResultType> RangeMapper;
1.8 + typedef LinearIntensityTransform<ResultType,typename NumericTraits<ResultType>::RealPromote> RangeMapper;
1.9
1.10 public:
1.11
1.12 @@ -118,11 +118,12 @@
1.13
1.14 inline ResultType operator()(const PixelType& a, const PixelType& b) const {
1.15 typedef typename NumericTraits<PixelType>::isScalar src_is_scalar;
1.16 - return diff(a, b, src_is_scalar());
1.17 + typedef typename NumericTraits<PixelType>::isIntegral src_is_integral;
1.18 + return diff(a, b, src_is_scalar(), src_is_integral());
1.19 }
1.20
1.21 protected:
1.22 - inline ResultType diff(const PixelType& a, const PixelType& b, VigraFalseType) const {
1.23 + inline ResultType diff(const PixelType& a, const PixelType& b, VigraFalseType, VigraTrueType) const {
1.24 PixelComponentType aLum = a.luminance();
1.25 PixelComponentType bLum = b.luminance();
1.26 PixelComponentType aHue = a.hue();
1.27 @@ -134,20 +135,77 @@
1.28 return rm(std::max(hueDiff, lumDiff));
1.29 }
1.30
1.31 - inline ResultType diff(const PixelType& a, const PixelType& b, VigraTrueType) const {
1.32 - typedef typename NumericTraits<PixelType>::isSigned src_is_signed;
1.33 - return scalar_diff(a, b, src_is_signed());
1.34 + // Convert an RGB triple to XYZ and return in the same data structure
1.35 + // Ref: Reinhard et al, "High Dynamic Range Imaging", pg 80
1.36 + inline PixelType rgb2xyz(const PixelType &a) const {
1.37 + PixelType result;
1.38 + // Use AdobeRGB conversion matrix (assume HDR files have AdobeRGB whitepoint; if not, error in delta E slight anyway)
1.39 + result.setRed(.5767*a.red()+0.1856*a.green()+0.1882*a.blue());
1.40 + result.setGreen(.2974*a.red()+0.6273*a.green()+0.0753*a.blue());
1.41 + result.setBlue(.0270*a.red()+0.0707*a.green()+0.9911*a.blue());
1.42 + return result;
1.43 }
1.44
1.45 - inline ResultType scalar_diff(const PixelType& a, const PixelType& b, VigraTrueType) const {
1.46 + // Convert an XYZ triple to a Lab-like space and return in the same data structure
1.47 + // The full Lab transform has a linear region when X,Y,Z are less than 0.008856, but here we ignore this
1.48 + // since the XYZ are not properly scaled by the illuminant
1.49 + // The objective is only to measure delta-E, so this shouldn't have a significant
1.50 + // Ref: Reinhard et al, "High Dynamic Range Imaging", pg 61
1.51 + inline PixelType xyz2lab(const PixelType &a) const {
1.52 + PixelType result;
1.53 + double fx=pow(a.red(),0.33);
1.54 + double fy=pow(a.green(),0.33);
1.55 + double fz=pow(a.blue(),0.33);
1.56 + // Use only power-law region of conversion since HDR file could be arbitrarily scaled
1.57 + result.setRed(116*fy);
1.58 + result.setGreen(500*(fx-fy));
1.59 + result.setBlue(200*(fy-fz));
1.60 + return result;
1.61 + }
1.62 + // Version for non-integral pixels -- assume to be linear RGB with arbitrary scaling
1.63 + inline ResultType diff(const PixelType & a, const PixelType & b, VigraFalseType, VigraFalseType) const {
1.64 + PixelType axyz=rgb2xyz(a);
1.65 + PixelType bxyz=rgb2xyz(b);
1.66 + PixelType alab=xyz2lab(axyz);
1.67 + PixelType blab=xyz2lab(bxyz);
1.68 + double diff=(alab-blab).magnitude();
1.69 + if (Verbose > VERBOSE_MASK_MESSAGES) {
1.70 + static int firsttime=10;
1.71 + if (firsttime && diff!=0) {
1.72 + cout << "diff(" << a << "," << b << "); alab=" << alab << ',' << "; blab=" << blab << "; diff=" << diff << endl;
1.73 + *(int *)&firsttime = firsttime-1;
1.74 + }
1.75 + }
1.76 + return diff;
1.77 + }
1.78 +
1.79 + inline ResultType diff(const PixelType & a, const PixelType & b, VigraTrueType, VigraTrueType) const {
1.80 + typedef typename NumericTraits<PixelType>::isSigned src_is_signed;
1.81 + typedef typename NumericTraits<PixelType>::isIntegral src_is_integral;
1.82 + return scalar_diff(a, b, src_is_signed(), src_is_integral());
1.83 + }
1.84 +
1.85 + inline ResultType diff(const PixelType & a, const PixelType & b, VigraTrueType, VigraFalseType) const {
1.86 + typedef typename NumericTraits<PixelType>::isSigned src_is_signed;
1.87 + typedef typename NumericTraits<PixelType>::isIntegral src_is_integral;
1.88 + return scalar_diff(a, b, src_is_signed(),src_is_integral());
1.89 + }
1.90 +
1.91 + inline ResultType scalar_diff(const PixelType & a, const PixelType & b, VigraTrueType, VigraTrueType) const {
1.92 return rm(std::abs(a - b));
1.93 }
1.94
1.95 // This appears necessary because NumericTraits<unsigned int>::Promote is an unsigned int instead of an int.
1.96 - inline ResultType scalar_diff(const PixelType& a, const PixelType& b, VigraFalseType) const {
1.97 + inline ResultType scalar_diff(const PixelType & a, const PixelType & b, VigraFalseType, VigraTrueType) const {
1.98 return rm(std::abs(static_cast<int>(a) - static_cast<int>(b)));
1.99 }
1.100
1.101 + // The following ones take non-integral scalar pixels; they don't need to rescale
1.102 + // I don't think this will ever be called
1.103 + inline ResultType scalar_diff(const PixelType & a, const PixelType & b, VigraTrueType, VigraFalseType) const {
1.104 + return std::abs(a - b);
1.105 + }
1.106 +
1.107 RangeMapper rm;
1.108
1.109 };
1.110 @@ -672,14 +730,16 @@
1.111 mismatchImageSize = vBB.size();
1.112 }
1.113
1.114 - typedef UInt8 MismatchImagePixelType;
1.115 - EnblendNumericTraits<MismatchImagePixelType>::ImageType mismatchImage(mismatchImageSize,
1.116 - NumericTraits<MismatchImagePixelType>::max());
1.117 + // To handle HDR data, we need a floating-point mismatch type to handle regions of the image
1.118 + // which may differ by more than an 8-bit pixel can accomodate.
1.119 + typedef typename EnblendNumericTraits<ImagePixelType>::ImagePixelComponentType MismatchImagePixelType;
1.120 +
1.121 + typename EnblendNumericTraits<MismatchImagePixelType>::ImageType mismatchImage(mismatchImageSize, NumericTraits<MismatchImagePixelType>::max());
1.122
1.123 // Visualization of optimization output
1.124 - EnblendNumericTraits<RGBValue<MismatchImagePixelType> >::ImageType* visualizeImage = NULL;
1.125 + typename EnblendNumericTraits<RGBValue<MismatchImagePixelType> >::ImageType *visualizeImage = NULL;
1.126 if (!VisualizeMaskFileName.empty()) {
1.127 - visualizeImage = new EnblendNumericTraits<RGBValue<MismatchImagePixelType> >::ImageType(mismatchImageSize);
1.128 + visualizeImage = new typename EnblendNumericTraits<RGBValue<MismatchImagePixelType> >::ImageType(mismatchImageSize);
1.129 }
1.130
1.131 // mem usage after: Visualize && CoarseMask: iBB * UInt8
1.132 @@ -911,9 +971,34 @@
1.133 }
1.134
1.135 if (visualizeImage) {
1.136 + static int visualizeCount = 0;
1.137 + static char *origname = 0;
1.138 + static char *origsuffix = 0;
1.139 + cout << "Saving visualize mask to " << VisualizeMaskFileName << endl;
1.140 ImageExportInfo visualizeInfo(VisualizeMaskFileName.c_str());
1.141 exportImage(srcImageRange(*visualizeImage), visualizeInfo);
1.142 delete visualizeImage;
1.143 +
1.144 +
1.145 + // In case we are doing multiple blends, augment the filename with an index
1.146 + if (visualizeCount==0) {
1.147 + int flen = strlen(VisualizeMaskFileName)+5;
1.148 + origname=VisualizeMaskFileName;
1.149 + char *lastdot=rindex(VisualizeMaskFileName,'.');
1.150 + if (lastdot != NULL) {
1.151 + *lastdot=0;
1.152 + origsuffix=lastdot+1;
1.153 + } else {
1.154 + origsuffix=0;
1.155 + }
1.156 + VisualizeMaskFileName = new char[flen];
1.157 + }
1.158 + visualizeCount++;
1.159 + sprintf(VisualizeMaskFileName,"%s_%d",origname,visualizeCount);
1.160 + if (origsuffix) {
1.161 + strcat(VisualizeMaskFileName,".");
1.162 + strcat(VisualizeMaskFileName,origsuffix);
1.163 + }
1.164 }
1.165
1.166 // Fill contours to get final optimized mask.