Vignetting correction from LCP profile

[StefanGyuzelev]

With current version 26.06.2017 Highlight looks good. With same settings versions 5.2 and current 06.09.2017 Highlight in some images looks very bad. As Аgriggio Alberto assumed the problem is in Vignetting correction from LCP profile.

https://drive.google.com/file/d/0Bwd1yqYI35ModW1pc1A1MzM5aDg/view?usp=sharing https://drive.google.com/file/d/0Bwd1yqYI35MoMlp0X3hYVDZERFE/view?usp=sharing https://drive.google.com/file/d/0Bwd1yqYI35MoRnlNTUhpZ0J5SVU/view?usp=sharing https://drive.google.com/file/d/0Bwd1yqYI35MoaGpLN3ZRWkF3d1k/view?usp=sharing https://drive.google.com/file/d/0Bwd1yqYI35MockJvZ3JrYW9BSFk/view?usp=sharing

[agriggio]

Please test this patch if you can.

diff --git a/rtengine/lcp.cc b/rtengine/lcp.cc
--- a/rtengine/lcp.cc
+++ b/rtengine/lcp.cc
@@ -134,11 +134,20 @@
 }

 // mode: 0=distortion, 1=vignette, 2=CA
-bool LCPPersModel::hasModeData(int mode) const
+bool LCPPersModel::hasModeData(LCPCorrectionMode mode) const
 {
-    return (mode == 0 && !vignette.empty() && !vignette.bad_error) || (mode == 1 && !base.empty() && !base.bad_error)
-           || (mode == 2 && !chromRG.empty() && !chromG.empty() && !chromBG.empty() &&
-               !chromRG.bad_error && !chromG.bad_error && !chromBG.bad_error);
+    switch (mode) {
+    case LCP_MODE_VIGNETTE:
+        return !vignette.empty() && !vignette.bad_error;
+    case LCP_MODE_DISTORTION:
+        return !base.empty() && !base.bad_error;
+    case LCP_MODE_CA:
+        return !chromRG.empty() && !chromG.empty() && !chromBG.empty() &&
+            !chromRG.bad_error && !chromG.bad_error && !chromBG.bad_error;
+    default:
+        assert(false);
+        return false;
+    }
 }

 void LCPPersModel::print() const
@@ -195,11 +204,11 @@
         printf("Vign: %i, fullWidth: %i/%i, focLen %g SwapXY: %i / MirX/Y %i / %i on rot:%i from %i\n",vignette, fullWidth, fullHeight, focalLength, swapXY, mirrorX, mirrorY, rot, rawRotationDeg);
     }

-    pProf->calcParams(vignette ? 0 : 1, focalLength, focusDist, aperture, &mc, nullptr, nullptr);
+    pProf->calcParams(vignette ? LCP_MODE_VIGNETTE : LCP_MODE_DISTORTION, focalLength, focusDist, aperture, &mc, nullptr, nullptr);
     mc.prepareParams(fullWidth, fullHeight, focalLength, focalLength35mm, pProf->sensorFormatFactor, swapXY, mirrorX, mirrorY);

     if (!vignette) {
-        pProf->calcParams(2, focalLength, focusDist, aperture, &chrom[0], &chrom[1], &chrom[2]);
+        pProf->calcParams(LCP_MODE_CA, focalLength, focusDist, aperture, &chrom[0], &chrom[1], &chrom[2]);

         for (int i = 0; i < 3; i++) {
             chrom[i].prepareParams(fullWidth, fullHeight, focalLength, focalLength35mm, pProf->sensorFormatFactor, swapXY, mirrorX, mirrorY);
@@ -410,9 +419,11 @@
     }
     // Two phase filter: first filter out the very rough ones, that distord the average a lot
     // force it, even if there are few frames (community profiles)
-//    filterBadFrames(2.0, 0);
+    filterBadFrames(LCP_MODE_VIGNETTE, 2.0, 0);
+    filterBadFrames(LCP_MODE_CA, 2.0, 0);
     // from the non-distorded, filter again on new average basis, but only if there are enough frames left
-//    filterBadFrames(1.5, 100);
+    filterBadFrames(LCP_MODE_VIGNETTE, 1.5, 50);
+    filterBadFrames(LCP_MODE_CA, 1.5, 50);
 }

@@ -429,67 +440,66 @@
 }

 // from all frames not marked as bad already, take average and filter out frames with higher deviation than this if there are enough values
-int LCPProfile::filterBadFrames(double maxAvgDevFac, int minFramesLeft)
+int LCPProfile::filterBadFrames(LCPCorrectionMode mode, double maxAvgDevFac, int minFramesLeft)
 {
-    // take average error per type, then calculated the maximum deviation allowed
-    double errBase = 0, errChrom = 0, errVignette = 0;
-    int baseCount = 0, chromCount = 0, vignetteCount = 0;
+    // take average error, then calculated the maximum deviation allowed
+    double err = 0;
+    int count = 0;

     for (int pm = 0; pm < MaxPersModelCount && aPersModel[pm]; pm++) {
-        if (aPersModel[pm]->hasModeData(0)) {
-            errVignette += aPersModel[pm]->vignette.mean_error;
-            vignetteCount++;
-        }
-
-        if (aPersModel[pm]->hasModeData(1)) {
-            errBase += aPersModel[pm]->base.mean_error;
-            baseCount++;
-        }
-
-        if (aPersModel[pm]->hasModeData(2)) {
-            errChrom += rtengine::max(aPersModel[pm]->chromRG.mean_error, aPersModel[pm]->chromG.mean_error, aPersModel[pm]->chromBG.mean_error);
-            chromCount++;
+        if (aPersModel[pm]->hasModeData(mode)) {
+            count++;
+            switch (mode) {
+            case LCP_MODE_VIGNETTE:
+                err += aPersModel[pm]->vignette.mean_error;
+                break;
+            case LCP_MODE_DISTORTION:
+                err += aPersModel[pm]->base.mean_error;
+                break;
+            case LCP_MODE_CA:
+                err += rtengine::max(aPersModel[pm]->chromRG.mean_error, aPersModel[pm]->chromG.mean_error, aPersModel[pm]->chromBG.mean_error);
+                break;
+            }
         }
     }

     // Only if we have enough frames, filter out errors
     int filtered = 0;

-    if (baseCount + chromCount + vignetteCount >= minFramesLeft) {
-        if (baseCount > 0) {
-            errBase /= (double)baseCount;
-        }
-
-        if (chromCount > 0) {
-            errChrom /= (double)chromCount;
-        }
-
-        if (vignetteCount > 0) {
-            errVignette /= (double)vignetteCount;
+    if (count >= minFramesLeft) {
+        if (count > 0) {
+            err /= (double)count;
         }

         // Now mark all the bad ones as bad, and hasModeData will return false;
         for (int pm = 0; pm < MaxPersModelCount && aPersModel[pm]; pm++) {
-            if (aPersModel[pm]->hasModeData(0) && aPersModel[pm]->vignette.mean_error > maxAvgDevFac * errVignette) {
-                aPersModel[pm]->vignette.bad_error = true;
-                filtered++;
-            }
-
-            if (aPersModel[pm]->hasModeData(1) && aPersModel[pm]->base.mean_error > maxAvgDevFac * errBase) {
-                aPersModel[pm]->base.bad_error = true;
-                filtered++;
-            }
-
-            if (aPersModel[pm]->hasModeData(2) &&
-                    (aPersModel[pm]->chromRG.mean_error > maxAvgDevFac * errChrom || aPersModel[pm]->chromG.mean_error > maxAvgDevFac * errChrom
-                     || aPersModel[pm]->chromBG.mean_error > maxAvgDevFac * errChrom)) {
-                aPersModel[pm]->chromRG.bad_error = aPersModel[pm]->chromG.bad_error = aPersModel[pm]->chromBG.bad_error = true;
-                filtered++;
+            if (aPersModel[pm]->hasModeData(mode)) {
+                switch (mode) {
+                case LCP_MODE_VIGNETTE:
+                    if (aPersModel[pm]->vignette.mean_error > maxAvgDevFac * err) {
+                        aPersModel[pm]->vignette.bad_error = true;
+                        filtered++;
+                    }
+                    break;
+                case LCP_MODE_DISTORTION:
+                    if (aPersModel[pm]->base.mean_error > maxAvgDevFac * err) {
+                        aPersModel[pm]->base.bad_error = true;
+                        filtered++;
+                    }
+                    break;
+                case LCP_MODE_CA:
+                    if ((aPersModel[pm]->chromRG.mean_error > maxAvgDevFac * err || aPersModel[pm]->chromG.mean_error > maxAvgDevFac * err
+                     || aPersModel[pm]->chromBG.mean_error > maxAvgDevFac * err)) {
+                        aPersModel[pm]->chromRG.bad_error = aPersModel[pm]->chromG.bad_error = aPersModel[pm]->chromBG.bad_error = true;
+                        filtered++;
+                    }
+                    break;
+                }
             }
         }

         if (settings->verbose) {
-            printf("Filtered %.1f%% frames for maxAvgDevFac %g leaving %i\n", filtered*100./(baseCount+chromCount+vignetteCount), maxAvgDevFac, baseCount+chromCount+vignetteCount-filtered);
+            printf("Filtered %.1f%% frames for maxAvgDevFac %g leaving %i\n", filtered*100./count, maxAvgDevFac, count-filtered);
         }
     }

@@ -497,8 +507,7 @@
 }

-// mode: 0=vignette, 1=distortion, 2=CA
-void LCPProfile::calcParams(int mode, float focalLength, float focusDist, float aperture, LCPModelCommon *pCorr1, LCPModelCommon *pCorr2, LCPModelCommon *pCorr3) const
+void LCPProfile::calcParams(LCPCorrectionMode mode, float focalLength, float focusDist, float aperture, LCPModelCommon *pCorr1, LCPModelCommon *pCorr2, LCPModelCommon *pCorr3) const
 {
     float euler = exp(1.0);

@@ -541,24 +550,24 @@
             if (aPersModel[pm]->hasModeData(mode)) {
                 double lowMeanErr, highMeanErr;
                 switch (mode) {
-                case 0:
+                case LCP_MODE_VIGNETTE:
                     meanErr = aPersModel[pm]->vignette.mean_error;
                     lowMeanErr = pLow->vignette.mean_error;
                     highMeanErr = pHigh->vignette.mean_error;
                     break;
-                case 1:
+                case LCP_MODE_DISTORTION:
                     meanErr = aPersModel[pm]->base.mean_error;
                     lowMeanErr = pLow->base.mean_error;
                     highMeanErr = pHigh->base.mean_error;
                     break;
-                default: //case 2:
+                default: // LCP_MODE_CA
                     meanErr = aPersModel[pm]->chromG.mean_error;
                     lowMeanErr = pLow->chromG.mean_error;
                     highMeanErr = pHigh->chromG.mean_error;
                     break;
                 }

-                if (aperture > 0 && mode != 2) {
+                if (aperture > 0 && mode != LCP_MODE_CA) {
                     if (aPersModel[pm]->focLen == bestFocLenLow && (
                             (aper == aperture && lowMeanErr > meanErr)
                             || (aper >= aperture && aper < pLow->aperture && pLow->aperture > aperture)
@@ -572,7 +581,7 @@
                             || (aper >= aperture && (pHigh->aperture < aperture || fabs(aperture - aper) < fabs(aperture - pHigh->aperture))))) {
                         pHigh = aPersModel[pm];
                     }
-                } else if (focusDist > 0 && mode != 0) {
+                } else if (focusDist > 0 && mode != LCP_MODE_VIGNETTE) {
                     // by focus distance
                     if (aPersModel[pm]->focLen == bestFocLenLow && (
                             (focDist == focusDist && lowMeanErr > meanErr)
@@ -615,26 +624,26 @@
         }

         // and average the other factor if available
-        if (mode == 0 && pLow->aperture < aperture && pHigh->aperture > aperture) {
+        if (mode == LCP_MODE_VIGNETTE && pLow->aperture < aperture && pHigh->aperture > aperture) {
             // Mix in aperture
             float facAperLow = (pHigh->aperture - aperture) / (pHigh->aperture - pLow->aperture);
             facLow = focLenOnSpot ? facAperLow : (0.5 * facLow + 0.5 * facAperLow);
-        } else if (mode != 0 && focusDist > 0 && pLow->focDist < focusDist && pHigh->focDist > focusDist) {
+        } else if (mode != LCP_MODE_VIGNETTE && focusDist > 0 && pLow->focDist < focusDist && pHigh->focDist > focusDist) {
             // focus distance for all else (if focus distance is given)
             float facDistLow = (log(pHigh->focDist) + euler - focusDistLog) / (log(pHigh->focDist) - log(pLow->focDist));
             facLow = focLenOnSpot ? facDistLow : (0.8 * facLow + 0.2 * facDistLow);
         }

         switch (mode) {
-        case 0:  // vignette
+        case LCP_MODE_VIGNETTE:
             pCorr1->merge(pLow->vignette, pHigh->vignette, facLow);
             break;

-        case 1:  // distortion
+        case LCP_MODE_DISTORTION:
             pCorr1->merge(pLow->base, pHigh->base, facLow);
             break;

-        case 2:  // CA
+        case LCP_MODE_CA:
             pCorr1->merge(pLow->chromRG, pHigh->chromRG, facLow);
             pCorr2->merge(pLow->chromG,  pHigh->chromG,  facLow);
             pCorr3->merge(pLow->chromBG, pHigh->chromBG, facLow);
@@ -646,7 +655,7 @@
         }
     } else {
         if (settings->verbose) {
-            printf("Error: LCP file contained no %s parameters\n", mode == 0 ? "vignette" : mode == 1 ? "distortion" : "CA" );
+            printf("Error: LCP file contained no %s parameters\n", mode == LCP_MODE_VIGNETTE ? "vignette" : mode == LCP_MODE_DISTORTION ? "distortion" : "CA" );
         }
     }
 }
diff --git a/rtengine/lcp.h b/rtengine/lcp.h
--- a/rtengine/lcp.h
+++ b/rtengine/lcp.h
@@ -33,6 +33,12 @@
 namespace rtengine
 {

+enum LCPCorrectionMode {
+    LCP_MODE_VIGNETTE = 0,
+    LCP_MODE_DISTORTION = 1,
+    LCP_MODE_CA = 2
+};
+
 // Perspective model common data, also used for Vignette and Fisheye
 class LCPModelCommon final
 {
@@ -76,7 +82,7 @@
     LCPModelCommon vignette;  // vignette (may be empty)

     LCPPersModel();
-    bool hasModeData(int mode) const;
+    bool hasModeData(LCPCorrectionMode mode) const;
     void print() const;
 };

@@ -93,7 +99,7 @@
     static void XMLCALL XmlTextHandler (void *pLCPProfile, const XML_Char *s, int len);
     static void XMLCALL XmlEndHandler  (void *pLCPProfile, const char *el);

-    int filterBadFrames(double maxAvgDevFac, int minFramesLeft);
+    int filterBadFrames(LCPCorrectionMode mode, double maxAvgDevFac, int minFramesLeft);

     void handle_text(std::string text);
     std::ostringstream textbuf;
@@ -112,7 +118,7 @@
     explicit LCPProfile(const Glib::ustring &fname);
     ~LCPProfile();

-    void calcParams(int mode, float focalLength, float focusDist, float aperture, LCPModelCommon *pCorr1, LCPModelCommon *pCorr2, LCPModelCommon *pCorr3) const;  // Interpolates between the persModels frames
+    void calcParams(LCPCorrectionMode mode, float focalLength, float focusDist, float aperture, LCPModelCommon *pCorr1, LCPModelCommon *pCorr2, LCPModelCommon *pCorr3) const;  // Interpolates between the persModels frames

     void print() const;
 };

[agriggio]

Also pushed as branch lcp-vignetting-issue4062

[StefanGyuzelev]

Hello! The version downloaded today and yours patch for lcp.cc and lcp.h works well! The problem is eliminated incredibly fast!

[Floessie]

@agriggio Alberto, before you merge: May I do some little cleanups on your branch? Just scopes for the cases and enum class for LCPCorrectionMode...

[agriggio]

@Floessie be my guest! :smile: BTW, I'm still not 100% sure this is the right way to go. I would still like to take a look at the related issue mentioned by @heckflosse on pixls. It might just be that this patch (which sort of restores the old behaviour of LCP vignetting correction) is simply hiding the real problem. I say this because what the patch does essentially is discarding some entries from the LCP file, when the mean error exceeds a certain threshold. I would assume that for the "official" LCP files from Adobe, there shouldn't be a need of discarding anything. So, this whole thing "smells" a little bit...

[heckflosse]

@agriggio Alberto, as vignetting correction by flatfield file, vignetting correction by lcp also multiplies a value to each pixel value. If that multiplicator can not be > 1.0, we do not need auto clip control. If the multiplicator can be > 1.0 (which should be the case, as the corrected example is brighter than the uncorrected example) it can happen that we get clipped values in formerly unclipped regions.

[Beep6581]

I would assume that for the "official" LCP files from Adobe, there shouldn't be a need of discarding anything

@agriggio I remember Oliver Duis, who originally implemented LCP support in RT, writing that LCP files (including Adobe's ones) contain erroneous control points and that one of the things the spec didn't define (well?) was how to filter them out. It also lacked information on interpolation between frames.

--edit e.g. https://github.com/Beep6581/RawTherapee/issues/1359#issuecomment-129963291 https://github.com/Beep6581/RawTherapee/issues/1413#issuecomment-129965209

Maybe that helps you.

[agriggio]

thanks, I'll have a look. I confirm the spec is definitely vague on many points. it doesn't even say how exactly to match the correction parameters. what is apparent is that the old code was too aggressive in discarding parameters for distortion, and the new one was too conservative in keeping those for vignetting :-) the above patch tries to get to a middle ground, but (as in the previous case, when I changed the code for 5.2) I have no idea how general this is... hopefully, with a bit of luck I might be able to integrate lensfun soon, so that the LCP business might become less critical...

[Floessie]

@agriggio Alberto, I took a bit more time to clean lcp.* up. Please test, if I messed up something. 😄

[agriggio]

@Floessie thanks for the cleanup! The changeset is a hard to review because there are a lot of "spurious" changes to indentation/formatting (I presume this is astyle, right?). But I trust you that everything is fine :-) I'll merge with the lensfun branch and test.

[Floessie]

@agriggio Alberto, sorry for the spurious changes, maybe I overdid it a bit. No, it wasn't astyled yet.

But I trust you that everything is fine :-)

Better don't. I was in a kind of hurry after I saw that not even the definition order followed the declaration order, so a lot had to be shuffled (I changed some whitespace while comparing the header to the implementation). I only checked that RT still started and could load RAWs... :)

One functional change I made was removing code that would exchange , with . for atof(). The first thing we do in main() is

setlocale (LC_NUMERIC, "C"); // to set decimal point to "."

so that step seemed unnecessary to me.

But I might have broken other parts, so please check. You have a better insight into that part of RT.

Best, Flössie

BTW: I don't even mind if you revert my changes. It's only a cleanup after all...

[agriggio]

@Floessie from my (quick) tests everything seems to work -- thanks! I'll merge this into lensfun and commit there

[Beep6581]

Can we close this issue and delete branch lcp-vignetting-issue4062?

[agriggio]

@Beep6581 sure, go ahead