diff --git a/rtengine/amaze_demosaic_RT.cc b/rtengine/amaze_demosaic_RT.cc
index 24c463458..2041f3130 100644
--- a/rtengine/amaze_demosaic_RT.cc
+++ b/rtengine/amaze_demosaic_RT.cc
@@ -1470,48 +1470,48 @@ void RawImageSource::amaze_demosaic_RT(int winx, int winy, int winw, int winh, c
                         vfloat bluev1 = greenv - (temp00v * vdup(LVFU(Dgrb[1][(indx - v1) >> 1])) + (onev - vdup(LVFU(hvwt[(indx + 1 + offset) >> 1]))) * vdup(LVFU(Dgrb[1][(indx + 1 + offset) >> 1])) + (onev - vdup(LVFU(hvwt[(indx - 1 + offset) >> 1]))) * vdup(LVFU(Dgrb[1][(indx - 1 + offset) >> 1])) + temp01v * vdup(LVFU(Dgrb[1][(indx + v1) >> 1]))) * tempv;
                         vfloat redv2  = greenv - vdup(LVFU(Dgrb[0][indx >> 1]));
                         vfloat bluev2 = greenv - vdup(LVFU(Dgrb[1][indx >> 1]));
-                        STVFU(red[row][col], c65535v * vself(selmask, redv1, redv2));
-                        STVFU(blue[row][col], c65535v * vself(selmask, bluev1, bluev2));
+                        STVFU(red[row][col], vmaxf(c65535v * vself(selmask, redv1, redv2), ZEROV));
+                        STVFU(blue[row][col], vmaxf(c65535v * vself(selmask, bluev1, bluev2), ZEROV));
                     }
 
                     if(offset == 0) {
                         for (; indx < rr * ts + cc1 - 16 - (cc1 & 1); indx++, col++) {
                             float temp =  1.f / (hvwt[(indx - v1) >> 1] + 2.f - hvwt[(indx + 1) >> 1] - hvwt[(indx - 1) >> 1] + hvwt[(indx + v1) >> 1]);
-                            red[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
-                                                        temp);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
-                                                         temp);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
+                                                        temp));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
+                                                         temp));
 
                             indx++;
                             col++;
-                            red[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]));
                         }
 
                         if(cc1 & 1) { // width of tile is odd
                             float temp =  1.f / (hvwt[(indx - v1) >> 1] + 2.f - hvwt[(indx + 1) >> 1] - hvwt[(indx - 1) >> 1] + hvwt[(indx + v1) >> 1]);
-                            red[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
-                                                        temp);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
-                                                         temp);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
+                                                        temp));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
+                                                         temp));
                         }
                     } else {
                         for (; indx < rr * ts + cc1 - 16 - (cc1 & 1); indx++, col++) {
-                            red[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]));
 
                             indx++;
                             col++;
                             float temp =  1.f / (hvwt[(indx - v1) >> 1] + 2.f - hvwt[(indx + 1) >> 1] - hvwt[(indx - 1) >> 1] + hvwt[(indx + v1) >> 1]);
-                            red[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
-                                                        temp);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
-                                                         temp);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
+                                                        temp));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
+                                                         temp));
                         }
 
                         if(cc1 & 1) { // width of tile is odd
-                            red[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]));
                         }
                     }
 
@@ -1520,41 +1520,41 @@ void RawImageSource::amaze_demosaic_RT(int winx, int winy, int winw, int winh, c
                     if((fc(cfarray, rr, 2) & 1) == 1) {
                         for (; indx < rr * ts + cc1 - 16 - (cc1 & 1); indx++, col++) {
                             float temp =  1.f / (hvwt[(indx - v1) >> 1] + 2.f - hvwt[(indx + 1) >> 1] - hvwt[(indx - 1) >> 1] + hvwt[(indx + v1) >> 1]);
-                            red[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
-                                                        temp);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
-                                                         temp);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
+                                                        temp));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
+                                                         temp));
 
                             indx++;
                             col++;
-                            red[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]));
                         }
 
                         if(cc1 & 1) { // width of tile is odd
                             float temp =  1.f / (hvwt[(indx - v1) >> 1] + 2.f - hvwt[(indx + 1) >> 1] - hvwt[(indx - 1) >> 1] + hvwt[(indx + v1) >> 1]);
-                            red[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
-                                                        temp);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
-                                                         temp);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
+                                                        temp));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
+                                                         temp));
                         }
                     } else {
                         for (; indx < rr * ts + cc1 - 16 - (cc1 & 1); indx++, col++) {
-                            red[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]));
 
                             indx++;
                             col++;
                             float temp =  1.f / (hvwt[(indx - v1) >> 1] + 2.f - hvwt[(indx + 1) >> 1] - hvwt[(indx - 1) >> 1] + hvwt[(indx + v1) >> 1]);
-                            red[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
-                                                        temp);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
-                                                         temp);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[0][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[0][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[0][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[0][(indx + v1) >> 1]) *
+                                                        temp));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - ((hvwt[(indx - v1) >> 1]) * Dgrb[1][(indx - v1) >> 1] + (1.f - hvwt[(indx + 1) >> 1]) * Dgrb[1][(indx + 1) >> 1] + (1.f - hvwt[(indx - 1) >> 1]) * Dgrb[1][(indx - 1) >> 1] + (hvwt[(indx + v1) >> 1]) * Dgrb[1][(indx + v1) >> 1]) *
+                                                         temp));
                         }
 
                         if(cc1 & 1) { // width of tile is odd
-                            red[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]);
-                            blue[row][col] = 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]);
+                            red[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[0][indx >> 1]));
+                            blue[row][col] = std::max(0.f, 65535.f * (rgbgreen[indx] - Dgrb[1][indx >> 1]));
                         }
                     }
 
@@ -1568,13 +1568,13 @@ void RawImageSource::amaze_demosaic_RT(int winx, int winy, int winw, int winh, c
 #ifdef __SSE2__
 
                     for (; cc < cc1 - 19; cc += 4) {
-                        STVFU(green[row][cc + left], LVF(rgbgreen[rr * ts + cc]) * c65535v);
+                        STVFU(green[row][cc + left], vmaxf(LVF(rgbgreen[rr * ts + cc]) * c65535v, ZEROV));
                     }
 
 #endif
 
                     for (; cc < cc1 - 16; cc++) {
-                        green[row][cc + left] = 65535.f * rgbgreen[rr * ts + cc];
+                        green[row][cc + left] = std::max(0.f, 65535.f * rgbgreen[rr * ts + cc]);
                     }
                 }
 
diff --git a/rtengine/demosaic_algos.cc b/rtengine/demosaic_algos.cc
index c48a002fb..538858a9a 100644
--- a/rtengine/demosaic_algos.cc
+++ b/rtengine/demosaic_algos.cc
@@ -215,7 +215,6 @@ void RawImageSource::border_interpolate( int winw, int winh, int lborders, const
 // Adapted to RawTherapee by Jacques Desmis 3/2013
 // SSE version by Ingo Weyrich 5/2013
 #ifdef __SSE2__
-#define CLIPV(a) vclampf(a,zerov,c65535v)
 void RawImageSource::igv_interpolate(int winw, int winh)
 {
     static const float eps = 1e-5f, epssq = 1e-5f; //mod epssq -10f =>-5f Jacques 3/2013 to prevent artifact (divide by zero)
@@ -284,9 +283,9 @@ void RawImageSource::igv_interpolate(int winw, int winh)
 
             for (col = 0, indx = row * width + col; col < width - 7; col += 8, indx += 8) {
                 temp1v = LVFU( rawData[row][col] );
-                temp1v = CLIPV( temp1v );
+                temp1v = vmaxf(temp1v, ZEROV);
                 temp2v = LVFU( rawData[row][col + 4] );
-                temp2v = CLIPV( temp2v );
+                temp2v = vmaxf(temp2v, ZEROV);
                 tempv = _mm_shuffle_ps( temp1v, temp2v, _MM_SHUFFLE( 2, 0, 2, 0 ) );
                 _mm_storeu_ps( &dest1[indx >> 1], tempv );
                 tempv = _mm_shuffle_ps( temp1v, temp2v, _MM_SHUFFLE( 3, 1, 3, 1 ) );
@@ -294,10 +293,10 @@ void RawImageSource::igv_interpolate(int winw, int winh)
             }
 
             for (; col < width; col++, indx += 2) {
-                dest1[indx >> 1] = CLIP(rawData[row][col]); //rawData = RT data
+                dest1[indx >> 1] = std::max(0.f, rawData[row][col]); //rawData = RT data
                 col++;
                 if(col < width)
-                    dest2[indx >> 1] = CLIP(rawData[row][col]); //rawData = RT data
+                    dest2[indx >> 1] = std::max(0.f, rawData[row][col]); //rawData = RT data
             }
         }
 
@@ -556,42 +555,42 @@ void RawImageSource::igv_interpolate(int winw, int winh)
                 temp2v = LVFU( src2[(indx + 1) >> 1] );
                 tempv = _mm_shuffle_ps( temp1v, temp2v, _MM_SHUFFLE( 1, 0, 1, 0 ) );
                 tempv = _mm_shuffle_ps( tempv, tempv, _MM_SHUFFLE( 3, 1, 2, 0 ) );
-                _mm_storeu_ps( &green[row][col], CLIPV( tempv ));
+                _mm_storeu_ps( &green[row][col], vmaxf(tempv, ZEROV));
                 temp5v = LVFU(redsrc0[indx >> 1]);
                 temp6v = LVFU(redsrc1[(indx + 1) >> 1]);
                 temp3v = _mm_shuffle_ps( temp5v, temp6v, _MM_SHUFFLE( 1, 0, 1, 0 ) );
                 temp3v = _mm_shuffle_ps( temp3v, temp3v, _MM_SHUFFLE( 3, 1, 2, 0 ) );
-                temp3v = CLIPV( tempv - c65535v * temp3v );
+                temp3v = vmaxf(tempv - c65535v * temp3v, ZEROV);
                 _mm_storeu_ps( &red[row][col], temp3v);
                 temp7v = LVFU(bluesrc0[indx >> 1]);
                 temp8v = LVFU(bluesrc1[(indx + 1) >> 1]);
                 temp4v = _mm_shuffle_ps( temp7v, temp8v, _MM_SHUFFLE( 1, 0, 1, 0 ) );
                 temp4v = _mm_shuffle_ps( temp4v, temp4v, _MM_SHUFFLE( 3, 1, 2, 0 ) );
-                temp4v = CLIPV( tempv - c65535v * temp4v );
+                temp4v = vmaxf(tempv - c65535v * temp4v, ZEROV);
                 _mm_storeu_ps( &blue[row][col], temp4v);
 
                 tempv = _mm_shuffle_ps( temp1v, temp2v, _MM_SHUFFLE( 3, 2, 3, 2 ) );
                 tempv = _mm_shuffle_ps( tempv, tempv, _MM_SHUFFLE( 3, 1, 2, 0 ) );
-                _mm_storeu_ps( &green[row][col + 4], CLIPV( tempv ));
+                _mm_storeu_ps( &green[row][col + 4], vmaxf(tempv, ZEROV));
 
                 temp3v = _mm_shuffle_ps( temp5v, temp6v, _MM_SHUFFLE( 3, 2, 3, 2 ) );
                 temp3v = _mm_shuffle_ps( temp3v, temp3v, _MM_SHUFFLE( 3, 1, 2, 0 ) );
-                temp3v = CLIPV( tempv - c65535v * temp3v );
+                temp3v = vmaxf(tempv - c65535v * temp3v, ZEROV);
                 _mm_storeu_ps( &red[row][col + 4], temp3v);
                 temp4v = _mm_shuffle_ps( temp7v, temp8v, _MM_SHUFFLE( 3, 2, 3, 2 ) );
                 temp4v = _mm_shuffle_ps( temp4v, temp4v, _MM_SHUFFLE( 3, 1, 2, 0 ) );
-                temp4v = CLIPV( tempv - c65535v * temp4v );
+                temp4v = vmaxf(tempv - c65535v * temp4v, ZEROV);
                 _mm_storeu_ps( &blue[row][col + 4], temp4v);
             }
 
             for(; col < width - 7; col++, indx += 2) {
-                red  [row][col] = CLIP(src1[indx >> 1] - 65535.f * redsrc0[indx >> 1]);
-                green[row][col] = CLIP(src1[indx >> 1]);
-                blue [row][col] = CLIP(src1[indx >> 1] - 65535.f * bluesrc0[indx >> 1]);
+                red  [row][col] = std::max(0.f, src1[indx >> 1] - 65535.f * redsrc0[indx >> 1]);
+                green[row][col] = std::max(0.f, src1[indx >> 1]);
+                blue [row][col] = std::max(0.f, src1[indx >> 1] - 65535.f * bluesrc0[indx >> 1]);
                 col++;
-                red  [row][col] = CLIP(src2[(indx + 1) >> 1] - 65535.f * redsrc1[(indx + 1) >> 1]);
-                green[row][col] = CLIP(src2[(indx + 1) >> 1]);
-                blue [row][col] = CLIP(src2[(indx + 1) >> 1] - 65535.f * bluesrc1[(indx + 1) >> 1]);
+                red  [row][col] = std::max(0.f, src2[(indx + 1) >> 1] - 65535.f * redsrc1[(indx + 1) >> 1]);
+                green[row][col] = std::max(0.f, src2[(indx + 1) >> 1]);
+                blue [row][col] = std::max(0.f, src2[(indx + 1) >> 1] - 65535.f * bluesrc1[(indx + 1) >> 1]);
             }
         }
     }// End of parallelization
@@ -606,7 +605,6 @@ void RawImageSource::igv_interpolate(int winw, int winh)
     free(vdif);
     free(hdif);
 }
-#undef CLIPV
 #else
 void RawImageSource::igv_interpolate(int winw, int winh)
 {
@@ -649,7 +647,7 @@ void RawImageSource::igv_interpolate(int winw, int winh)
         for (int row = 0; row < height - 0; row++)
             for (int col = 0, indx = row * width + col; col < width - 0; col++, indx++) {
                 int c = FC(row, col);
-                rgb[c][indx] = CLIP(rawData[row][col]); //rawData = RT data
+                rgb[c][indx] = std::max(0.f, rawData[row][col]); //rawData = RT data
             }
 
 #ifdef _OPENMP
@@ -848,9 +846,9 @@ void RawImageSource::igv_interpolate(int winw, int winh)
 
         for(int row = 7; row < height - 7; row++)
             for(int col = 7, indx = row * width + col; col < width - 7; col++, indx++) {
-                red  [row][col] = CLIP(rgb[1][indx] - 65535.f * chr[0][indx]);
-                green[row][col] = CLIP(rgb[1][indx]);
-                blue [row][col] = CLIP(rgb[1][indx] - 65535.f * chr[1][indx]);
+                red  [row][col] = std::max(0.f, rgb[1][indx] - 65535.f * chr[0][indx]);
+                green[row][col] = std::max(0.f, rgb[1][indx]);
+                blue [row][col] = std::max(0.f, rgb[1][indx] - 65535.f * chr[1][indx]);
             }
     }// End of parallelization
     border_interpolate(winw, winh, 8, rawData, red, green, blue);
@@ -1520,9 +1518,9 @@ BENCHFUN
 */
         for(int y = 0; y < TILESIZE && y0 + y < H; y++) {
             for (int j = 0; j < TILESIZE && x0 + j < W; j++) {
-                red[y0 + y][x0 + j]   = tile[(y + TILEBORDER) * CACHESIZE + TILEBORDER + j][0];
-                green[y0 + y][x0 + j] = tile[(y + TILEBORDER) * CACHESIZE + TILEBORDER + j][1];
-                blue[y0 + y][x0 + j]  = tile[(y + TILEBORDER) * CACHESIZE + TILEBORDER + j][2];
+                red[y0 + y][x0 + j]   = std::max(0.f, tile[(y + TILEBORDER) * CACHESIZE + TILEBORDER + j][0]);
+                green[y0 + y][x0 + j] = std::max(0.f, tile[(y + TILEBORDER) * CACHESIZE + TILEBORDER + j][1]);
+                blue[y0 + y][x0 + j]  = std::max(0.f, tile[(y + TILEBORDER) * CACHESIZE + TILEBORDER + j][2]);
             }
         }
 
diff --git a/rtengine/eahd_demosaic.cc b/rtengine/eahd_demosaic.cc
index c470eb297..ad4bda3cd 100644
--- a/rtengine/eahd_demosaic.cc
+++ b/rtengine/eahd_demosaic.cc
@@ -400,11 +400,11 @@ void RawImageSource::eahd_demosaic ()
                 int vc = homv[imx][j];
 
                 if (hc > vc) {
-                    green[i - 1][j] = gh[(i - 1) % 4][j];
+                    green[i - 1][j] = std::max(0.f, gh[(i - 1) % 4][j]);
                 } else if (hc < vc) {
-                    green[i - 1][j] = gv[(i - 1) % 4][j];
+                    green[i - 1][j] = std::max(0.f, gv[(i - 1) % 4][j]);
                 } else {
-                    green[i - 1][j] = (gh[(i - 1) % 4][j] + gv[(i - 1) % 4][j]) / 2;
+                    green[i - 1][j] = std::max(0.f, (gh[(i - 1) % 4][j] + gv[(i - 1) % 4][j]) / 2);
                 }
             }
         }
@@ -421,11 +421,11 @@ void RawImageSource::eahd_demosaic ()
             int vc = homv[(i - 1) % 3][j];
 
             if (hc > vc) {
-                green[i - 1][j] = gh[(i - 1) % 4][j];
+                green[i - 1][j] = std::max(0.f, gh[(i - 1) % 4][j]);
             } else if (hc < vc) {
-                green[i - 1][j] = gv[(i - 1) % 4][j];
+                green[i - 1][j] = std::max(0.f, gv[(i - 1) % 4][j]);
             } else {
-                green[i - 1][j] = (gh[(i - 1) % 4][j] + gv[(i - 1) % 4][j]) / 2;
+                green[i - 1][j] = std::max(0.f, (gh[(i - 1) % 4][j] + gv[(i - 1) % 4][j]) / 2);
             }
         }
 
diff --git a/rtengine/fast_demo.cc b/rtengine/fast_demo.cc
index aaa5ceb02..772096b87 100644
--- a/rtengine/fast_demo.cc
+++ b/rtengine/fast_demo.cc
@@ -452,23 +452,23 @@ void RawImageSource::fast_demosaic()
 #ifdef __SSE2__
 
                     for (j = left + 2, cc = 2; j < right - 5; j += 4, cc += 4) {
-                        _mm_storeu_ps(&red[i][j], LVFU(redtile[rr * TS + cc]));
-                        _mm_storeu_ps(&green[i][j], LVFU(greentile[rr * TS + cc]));
-                        _mm_storeu_ps(&blue[i][j], LVFU(bluetile[rr * TS + cc]));
+                        _mm_storeu_ps(&red[i][j], vmaxf(LVFU(redtile[rr * TS + cc]), ZEROV));
+                        _mm_storeu_ps(&green[i][j], vmaxf(LVFU(greentile[rr * TS + cc]), ZEROV));
+                        _mm_storeu_ps(&blue[i][j], vmaxf(LVFU(bluetile[rr * TS + cc]), ZEROV));
                     }
 
                     for (; j < right - 2; j++, cc++) {
-                        red[i][j] = redtile[rr * TS + cc];
-                        green[i][j] = greentile[rr * TS + cc];
-                        blue[i][j] = bluetile[rr * TS + cc];
+                        red[i][j] = std::max(0.f, redtile[rr * TS + cc]);
+                        green[i][j] = std::max(0.f, greentile[rr * TS + cc]);
+                        blue[i][j] = std::max(0.f, bluetile[rr * TS + cc]);
                     }
 
 #else
 
                     for (int j = left + 2, cc = 2; j < right - 2; j++, cc++) {
-                        red[i][j] = redtile[rr * TS + cc];
-                        green[i][j] = greentile[rr * TS + cc];
-                        blue[i][j] = bluetile[rr * TS + cc];
+                        red[i][j] = std::max(0.f, redtile[rr * TS + cc]);
+                        green[i][j] = std::max(0.f, greentile[rr * TS + cc]);
+                        blue[i][j] = std::max(0.f, bluetile[rr * TS + cc]);
                     }
 
 #endif
diff --git a/rtengine/hphd_demosaic_RT.cc b/rtengine/hphd_demosaic_RT.cc
index 5e05b128e..936c8fa4f 100644
--- a/rtengine/hphd_demosaic_RT.cc
+++ b/rtengine/hphd_demosaic_RT.cc
@@ -217,7 +217,7 @@ void hphd_green(const RawImage *ri, const array2D<float> &rawData, float** hpmap
 
                     const float e4 = 1.f / (dx + (std::fabs(d1) + std::fabs(d2)) + (std::fabs(d3) + std::fabs(d4)) * 0.5f);
 
-                    green[i][j] = rawData[i][j] * 0.5f + (e2 * g2 + e4 * g4) / (e2 + e4);
+                    green[i][j] = std::max(0.f, rawData[i][j] * 0.5f + (e2 * g2 + e4 * g4) / (e2 + e4));
                 } else if (hpmap[i][j] == 2) {
                     const float g1 = rawData[i - 1][j] - rawData[i - 2][j] * 0.5f;
                     const float g3 = rawData[i + 1][j] - rawData[i + 2][j] * 0.5f;
@@ -237,7 +237,7 @@ void hphd_green(const RawImage *ri, const array2D<float> &rawData, float** hpmap
 
                     const float e3 = 1.f / (dy + (std::fabs(d1) + std::fabs(d2)) + (std::fabs(d3) + std::fabs(d4)) * 0.5f);
 
-                    green[i][j] = rawData[i][j] * 0.5f + (e1 * g1 + e3 * g3) / (e1 + e3);
+                    green[i][j] = std::max(0.f, rawData[i][j] * 0.5f + (e1 * g1 + e3 * g3) / (e1 + e3));
                 } else {
                     const float g1 = rawData[i - 1][j] - rawData[i - 2][j] * 0.5f;
                     const float g2 = rawData[i][j + 1] - rawData[i][j + 2] * 0.5f;
@@ -275,7 +275,7 @@ void hphd_green(const RawImage *ri, const array2D<float> &rawData, float** hpmap
 
                     const float e4 = 1.f / (dx + (std::fabs(d1) + std::fabs(d2)) + (std::fabs(d3) + std::fabs(d4)) * 0.5f);
 
-                    green[i][j] = rawData[i][j] * 0.5f + ((e1 * g1 + e2 * g2) + (e3 * g3 + e4 * g4)) / (e1 + e2 + e3 + e4);
+                    green[i][j] = std::max(0.f, rawData[i][j] * 0.5f + ((e1 * g1 + e2 * g2) + (e3 * g3 + e4 * g4)) / (e1 + e2 + e3 + e4));
                 }
             }
         }
diff --git a/rtengine/lmmse_demosaic.cc b/rtengine/lmmse_demosaic.cc
index 6191ca36c..f4af1b801 100644
--- a/rtengine/lmmse_demosaic.cc
+++ b/rtengine/lmmse_demosaic.cc
@@ -131,7 +131,7 @@ void RawImageSource::lmmse_interpolate_omp(int winw, int winh, const array2D<flo
     if (applyGamma) {
         gamtab = &(Color::gammatab_24_17a);
     } else {
-        gamtab = new LUTf(65536, LUT_CLIP_ABOVE | LUT_CLIP_BELOW);
+        gamtab = new LUTf(65536, LUT_CLIP_BELOW);
         gamtab->makeIdentity(65535.f);
     }
 
@@ -623,7 +623,7 @@ void RawImageSource::lmmse_interpolate_omp(int winw, int winh, const array2D<flo
             for (int ii = 0; ii < 3; ii++)
                 if (ii != c) {
                     float *rix = qix[ii] + rr * cc1 + cc;
-                    (*(rgb[ii]))[row][col] = (*gamtab)[65535.f * rix[0]];
+                    (*(rgb[ii]))[row][col] = std::max(0.f, (*gamtab)[65535.f * rix[0]]);
                 } else {
                     (*(rgb[ii]))[row][col] = CLIP(rawData[row][col]);
                 }
@@ -651,9 +651,6 @@ void RawImageSource::lmmse_interpolate_omp(int winw, int winh, const array2D<flo
 
 }
 
-#ifdef __SSE2__
-#define CLIPV(a) vclampf(a,ZEROV,c65535v)
-#endif
 void RawImageSource::refinement(int PassCount)
 {
     int width = W;
@@ -694,7 +691,6 @@ void RawImageSource::refinement(int PassCount)
                 vfloat dLv, dRv, dUv, dDv, v0v;
                 vfloat onev = F2V(1.f);
                 vfloat zd5v = F2V(0.5f);
-                vfloat c65535v = F2V(65535.f);
 
                 for (; col < width - 8; col += 8) {
                     int indx = row * width + col;
@@ -704,7 +700,7 @@ void RawImageSource::refinement(int PassCount)
                     dRv = onev / (onev + vabsf(LC2VFU(pix[c][  2]) - LC2VFU(pix[c][0])) + vabsf(LC2VFU(pix[1][ 1]) - LC2VFU(pix[1][ -1])));
                     dUv = onev / (onev + vabsf(LC2VFU(pix[c][-w2]) - LC2VFU(pix[c][0])) + vabsf(LC2VFU(pix[1][w1]) - LC2VFU(pix[1][-w1])));
                     dDv = onev / (onev + vabsf(LC2VFU(pix[c][ w2]) - LC2VFU(pix[c][0])) + vabsf(LC2VFU(pix[1][w1]) - LC2VFU(pix[1][-w1])));
-                    v0v = CLIPV(LC2VFU(pix[c][0]) + zd5v + ((LC2VFU(pix[1][-1]) - LC2VFU(pix[c][-1])) * dLv + (LC2VFU(pix[1][1]) - LC2VFU(pix[c][1])) * dRv + (LC2VFU(pix[1][-w1]) - LC2VFU(pix[c][-w1])) * dUv + (LC2VFU(pix[1][w1]) - LC2VFU(pix[c][w1])) * dDv ) / (dLv + dRv + dUv + dDv));
+                    v0v = vmaxf(ZEROV, LC2VFU(pix[c][0]) + zd5v + ((LC2VFU(pix[1][-1]) - LC2VFU(pix[c][-1])) * dLv + (LC2VFU(pix[1][1]) - LC2VFU(pix[c][1])) * dRv + (LC2VFU(pix[1][-w1]) - LC2VFU(pix[c][-w1])) * dUv + (LC2VFU(pix[1][w1]) - LC2VFU(pix[c][w1])) * dDv ) / (dLv + dRv + dUv + dDv));
                     STC2VFU(pix[1][0], v0v);
                 }
 
@@ -719,7 +715,7 @@ void RawImageSource::refinement(int PassCount)
                     float dU = 1.f / (1.f + fabsf(pix[c][-w2] - pix[c][0]) + fabsf(pix[1][w1] - pix[1][-w1]));
                     float dD = 1.f / (1.f + fabsf(pix[c][ w2] - pix[c][0]) + fabsf(pix[1][w1] - pix[1][-w1]));
                     float v0 = (pix[c][0] + 0.5f + ((pix[1][ -1] - pix[c][ -1]) * dL + (pix[1][  1] - pix[c][  1]) * dR + (pix[1][-w1] - pix[c][-w1]) * dU + (pix[1][ w1] - pix[c][ w1]) * dD ) / (dL + dR + dU + dD));
-                    pix[1][0] = CLIP(v0);
+                    pix[1][0] = std::max(0.f, v0);
                 }
             }
 
@@ -735,7 +731,6 @@ void RawImageSource::refinement(int PassCount)
                 vfloat dLv, dRv, dUv, dDv, v0v;
                 vfloat onev = F2V(1.f);
                 vfloat zd5v = F2V(0.5f);
-                vfloat c65535v = F2V(65535.f);
 
                 for (; col < width - 8; col += 8) {
                     int indx = row * width + col;
@@ -747,7 +742,7 @@ void RawImageSource::refinement(int PassCount)
                         dRv = onev / (onev + vabsf(LC2VFU(pix[1][  2]) - LC2VFU(pix[1][0])) + vabsf(LC2VFU(pix[c][ 1]) - LC2VFU(pix[c][ -1])));
                         dUv = onev / (onev + vabsf(LC2VFU(pix[1][-w2]) - LC2VFU(pix[1][0])) + vabsf(LC2VFU(pix[c][w1]) - LC2VFU(pix[c][-w1])));
                         dDv = onev / (onev + vabsf(LC2VFU(pix[1][ w2]) - LC2VFU(pix[1][0])) + vabsf(LC2VFU(pix[c][w1]) - LC2VFU(pix[c][-w1])));
-                        v0v = CLIPV(LC2VFU(pix[1][0]) + zd5v - ((LC2VFU(pix[1][-1]) - LC2VFU(pix[c][-1])) * dLv + (LC2VFU(pix[1][1]) - LC2VFU(pix[c][1])) * dRv + (LC2VFU(pix[1][-w1]) - LC2VFU(pix[c][-w1])) * dUv + (LC2VFU(pix[1][w1]) - LC2VFU(pix[c][w1])) * dDv ) / (dLv + dRv + dUv + dDv));
+                        v0v = vmaxf(ZEROV, LC2VFU(pix[1][0]) + zd5v - ((LC2VFU(pix[1][-1]) - LC2VFU(pix[c][-1])) * dLv + (LC2VFU(pix[1][1]) - LC2VFU(pix[c][1])) * dRv + (LC2VFU(pix[1][-w1]) - LC2VFU(pix[c][-w1])) * dUv + (LC2VFU(pix[1][w1]) - LC2VFU(pix[c][w1])) * dDv ) / (dLv + dRv + dUv + dDv));
                         STC2VFU(pix[c][0], v0v);
                     }
                 }
@@ -765,7 +760,7 @@ void RawImageSource::refinement(int PassCount)
                         float dU = 1.f / (1.f + fabsf(pix[1][-w2] - pix[1][0]) + fabsf(pix[c][w1] - pix[c][-w1]));
                         float dD = 1.f / (1.f + fabsf(pix[1][ w2] - pix[1][0]) + fabsf(pix[c][w1] - pix[c][-w1]));
                         float v0 = (pix[1][0] + 0.5f - ((pix[1][ -1] - pix[c][ -1]) * dL + (pix[1][  1] - pix[c][  1]) * dR + (pix[1][-w1] - pix[c][-w1]) * dU + (pix[1][ w1] - pix[c][ w1]) * dD ) / (dL + dR + dU + dD));
-                        pix[c][0] = CLIP(v0);
+                        pix[c][0] = std::max(0.f, v0);
                     }
                 }
             }
@@ -782,7 +777,6 @@ void RawImageSource::refinement(int PassCount)
                 vfloat dLv, dRv, dUv, dDv, v0v;
                 vfloat onev = F2V(1.f);
                 vfloat zd5v = F2V(0.5f);
-                vfloat c65535v = F2V(65535.f);
 
                 for (; col < width - 8; col += 8) {
                     int indx = row * width + col;
@@ -794,7 +788,7 @@ void RawImageSource::refinement(int PassCount)
                     dRv = onev / (onev + vabsf(LC2VFU(pix[d][  2]) - LC2VFU(pix[d][0])) + vabsf(LC2VFU(pix[1][ 1]) - LC2VFU(pix[1][ -1])));
                     dUv = onev / (onev + vabsf(LC2VFU(pix[d][-w2]) - LC2VFU(pix[d][0])) + vabsf(LC2VFU(pix[1][w1]) - LC2VFU(pix[1][-w1])));
                     dDv = onev / (onev + vabsf(LC2VFU(pix[d][ w2]) - LC2VFU(pix[d][0])) + vabsf(LC2VFU(pix[1][w1]) - LC2VFU(pix[1][-w1])));
-                    v0v = CLIPV(LC2VFU(pix[1][0]) + zd5v - ((LC2VFU(pix[1][-1]) - LC2VFU(pix[c][-1])) * dLv + (LC2VFU(pix[1][1]) - LC2VFU(pix[c][1])) * dRv + (LC2VFU(pix[1][-w1]) - LC2VFU(pix[c][-w1])) * dUv + (LC2VFU(pix[1][w1]) - LC2VFU(pix[c][w1])) * dDv ) / (dLv + dRv + dUv + dDv));
+                    v0v = vmaxf(ZEROV, LC2VFU(pix[1][0]) + zd5v - ((LC2VFU(pix[1][-1]) - LC2VFU(pix[c][-1])) * dLv + (LC2VFU(pix[1][1]) - LC2VFU(pix[c][1])) * dRv + (LC2VFU(pix[1][-w1]) - LC2VFU(pix[c][-w1])) * dUv + (LC2VFU(pix[1][w1]) - LC2VFU(pix[c][w1])) * dDv ) / (dLv + dRv + dUv + dDv));
                     STC2VFU(pix[c][0], v0v);
                 }
 
@@ -811,7 +805,7 @@ void RawImageSource::refinement(int PassCount)
                     float dU = 1.f / (1.f + fabsf(pix[d][-w2] - pix[d][0]) + fabsf(pix[1][w1] - pix[1][-w1]));
                     float dD = 1.f / (1.f + fabsf(pix[d][ w2] - pix[d][0]) + fabsf(pix[1][w1] - pix[1][-w1]));
                     float v0 = (pix[1][0] + 0.5f - ((pix[1][ -1] - pix[c][ -1]) * dL + (pix[1][  1] - pix[c][  1]) * dR + (pix[1][-w1] - pix[c][-w1]) * dU + (pix[1][ w1] - pix[c][ w1]) * dD ) / (dL + dR + dU + dD));
-                    pix[c][0] = CLIP(v0);
+                    pix[c][0] = std::max(0.f, v0);
                 }
             }
         } // end parallel
diff --git a/rtengine/rawimagesource_i.h b/rtengine/rawimagesource_i.h
index 47c6b5bab..9dfc87985 100644
--- a/rtengine/rawimagesource_i.h
+++ b/rtengine/rawimagesource_i.h
@@ -90,7 +90,7 @@ inline void RawImageSource::interpolate_row_rb_mul_pp (const array2D<float> &raw
                 }
 
                 b = g_mul * cg[j] + b / std::max(1, n);
-                ab[jx] = b;
+                ab[jx] = std::max(0.f, b);
             } else {
                 // linear R-G interp. horizontally
                 float r;
@@ -103,7 +103,7 @@ inline void RawImageSource::interpolate_row_rb_mul_pp (const array2D<float> &raw
                     r = g_mul * cg[j] + (r_mul * rawData[i][j - 1] - g_mul * cg[j - 1] + r_mul * rawData[i][j + 1] - g_mul * cg[j + 1]) / 2;
                 }
 
-                ar[jx] = r;
+                ar[jx] = std::max(0.f, r);
                 // linear B-G interp. vertically
                 float b;
 
@@ -115,7 +115,7 @@ inline void RawImageSource::interpolate_row_rb_mul_pp (const array2D<float> &raw
                     b = g_mul * cg[j] + (b_mul * rawData[i - 1][j] - g_mul * pg[j] + b_mul * rawData[i + 1][j] - g_mul * ng[j]) / 2;
                 }
 
-                ab[jx] = b;
+                ab[jx] = std::max(0.f, b);
             }
         }
     } else if(pg && ng) {
@@ -150,7 +150,7 @@ inline void RawImageSource::interpolate_row_rb_mul_pp (const array2D<float> &raw
 
                 r = g_mul * cg[j] + r / std::max(n, 1);
 
-                ar[jx] = r;
+                ar[jx] = std::max(0.f, r);
             } else {
                 // linear B-G interp. horizontally
                 float b;
@@ -163,7 +163,7 @@ inline void RawImageSource::interpolate_row_rb_mul_pp (const array2D<float> &raw
                     b = g_mul * cg[j] + (b_mul * rawData[i][j - 1] - g_mul * cg[j - 1] + b_mul * rawData[i][j + 1] - g_mul * cg[j + 1]) / 2;
                 }
 
-                ab[jx] = b;
+                ab[jx] = std::max(0.f, b);
                 // linear R-G interp. vertically
                 float r;
 
@@ -175,7 +175,7 @@ inline void RawImageSource::interpolate_row_rb_mul_pp (const array2D<float> &raw
                     r = g_mul * cg[j] + (r_mul * rawData[i - 1][j] - g_mul * pg[j] + r_mul * rawData[i + 1][j] - g_mul * ng[j]) / 2;
                 }
 
-                ar[jx] = r;
+                ar[jx] = std::max(0.f, r);
             }
         }
     }
diff --git a/rtengine/rcd_demosaic.cc b/rtengine/rcd_demosaic.cc
index 4ceb92b26..5a86aec40 100644
--- a/rtengine/rcd_demosaic.cc
+++ b/rtengine/rcd_demosaic.cc
@@ -275,9 +275,9 @@ void RawImageSource::rcd_demosaic(size_t chunkSize, bool measure)
             for (int row = rowStart + rcdBorder; row < rowEnd - rcdBorder; ++row) {
                 for (int col = colStart + rcdBorder; col < colEnd - rcdBorder; ++col) {
                     int idx = (row - rowStart) * tileSize + col - colStart ;
-                    red[row][col] = CLIP(rgb[0][idx] * 65535.f);
-                    green[row][col] = CLIP(rgb[1][idx] * 65535.f);
-                    blue[row][col] = CLIP(rgb[2][idx] * 65535.f);
+                    red[row][col] = std::max(0.f, rgb[0][idx] * 65535.f);
+                    green[row][col] = std::max(0.f, rgb[1][idx] * 65535.f);
+                    blue[row][col] = std::max(0.f, rgb[2][idx] * 65535.f);
                 }
             }
 
diff --git a/rtengine/vng4_demosaic_RT.cc b/rtengine/vng4_demosaic_RT.cc
index 47982b6da..ef456af3a 100644
--- a/rtengine/vng4_demosaic_RT.cc
+++ b/rtengine/vng4_demosaic_RT.cc
@@ -45,12 +45,12 @@ inline void vng4interpolate_row_redblue (const RawImage *ri, const array2D<float
             // cross interpolation of red/blue
             float rb = (rawData[i - 1][j - 1] - pg[j - 1] + rawData[i + 1][j - 1] - ng[j - 1]);
             rb += (rawData[i - 1][j + 1] - pg[j + 1] + rawData[i + 1][j + 1] - ng[j + 1]);
-            ab[j] = cg[j] + rb * 0.25f;
+            ab[j] = std::max(0.f, cg[j] + rb * 0.25f);
         } else {
             // linear R/B-G interpolation horizontally
-            ar[j] = cg[j] + (rawData[i][j - 1] - cg[j - 1] + rawData[i][j + 1] - cg[j + 1]) / 2;
+            ar[j] = std::max(0.f, cg[j] + (rawData[i][j - 1] - cg[j - 1] + rawData[i][j + 1] - cg[j + 1]) / 2);
             // linear B/R-G interpolation vertically
-            ab[j] = cg[j] + (rawData[i - 1][j] - pg[j] + rawData[i + 1][j] - ng[j]) / 2;
+            ab[j] = std::max(0.f, cg[j] + (rawData[i - 1][j] - pg[j] + rawData[i + 1][j] - ng[j]) / 2);
         }
     }
 }
@@ -354,7 +354,7 @@ void RawImageSource::vng4_demosaic (const array2D<float> &rawData, array2D<float
                         }
                     }
                 }
-                green[row][col] = greenval + (sum1 - sum0) / (2 * num);
+                green[row][col] = std::max(0.f, greenval + (sum1 - sum0) / (2 * num));
             }
             if (row - 1 > firstRow) {
                 vng4interpolate_row_redblue(ri, rawData, red[row - 1], blue[row - 1], green[row - 2], green[row - 1], green[row], row - 1, W);
diff --git a/rtengine/xtrans_demosaic.cc b/rtengine/xtrans_demosaic.cc
index c1a526535..cb1315ed0 100644
--- a/rtengine/xtrans_demosaic.cc
+++ b/rtengine/xtrans_demosaic.cc
@@ -939,9 +939,9 @@ void RawImageSource::xtrans_interpolate (const int passes, const bool useCieLab,
                                 avg[3]++;
                             }
 
-                        red[row + top][col + left] = avg[0] / avg[3];
-                        green[row + top][col + left] = avg[1] / avg[3];
-                        blue[row + top][col + left] = avg[2] / avg[3];
+                        red[row + top][col + left] = std::max(0.f, avg[0] / avg[3]);
+                        green[row + top][col + left] = std::max(0.f, avg[1] / avg[3]);
+                        blue[row + top][col + left] = std::max(0.f, avg[2] / avg[3]);
                     }
 
                 if(plistenerActive && ((++progressCounter) % 32 == 0)) {