%ffp

Category   :"RC Filters"
Title      :"Very High Radius Unsharp Mask"
Filename   :"VHRUnsharpMask.8bf"
Description:"High-radius low-amount contrast enhancement for large images."
Copyright  :"© 2015 by Russell Cottrell; released under the Gnu General Public License."
Author     :"Russell Cottrell"
Version    :"1.112"
URL        :"http://www.russellcottrell.com"
About      :"Thank you for using !C.\n!t !V by !A.\n!D\nPlease visit russellcottrell.com."

// 12/5/2015

/*
Acknowledgements:

Recursive gaussian demo code
implemented by Tom Fiddaman
from the article: 
Title: Recursive implementation of the Gaussian filter
Authors: I.T. Young, L.J. van Vliet
in: Signal Processing, vol. 44, no. 2, 1995, 139-151.

Interpolation by Harald Heim
The Plugin Site, http://thepluginsite.com
Copyright 2002 by Harald Heim
*/

Dialog: Size=(590,318), MaxBox, MinBox, Theme=on, Initevent
SupportedModes: RGBMode, RGB48Mode, GrayScaleMode, Gray16Mode

ctl(0): Standard, "Radius", Range=(0,200), Divisor=10, pos=(400,20), Val=100
ctl(1): StaticText, "0px", pos=(525,20), size=(60,10)
ctl(2): Standard, "Amount", Range=(-40,100), pos=(400,32), Val=30
ctl(7): Scrollbar, "", pos=(525,31), size=(20,10), range=(0,400), val=0, line=400, invisible

ctl(3): checkbox, text="Show the Gaussian Blur", Pos=(400,52), val=0, Action=Preview

ctl(6): Standard, "Color", Range=(0,200), pos=(400,72), Val=0
ctl(10): checkbox(Pushlike, 3State), text="Normal", Pos=(525,69), Size=(56,12), val=1, Action=Preview

ctl(4): Standard, "Dark/Light Mixer", Range=(-100,100), pos=(425,92), Val=0

ctl(30): StaticText, "Shadow/Highlight Masks:", pos=(365,112), size=(160,10)

ctl(31): Standard, "Tonal Width", Range=(0,127), pos=(400,128), size=(115,*), Val=0, Color=Black
ctl(32): Standard, "Radius", Range=(0,200), Divisor=100, pos=(400,140), size=(115,*), Val=0, Color=Black
ctl(33): Standard, "Amount", Range=(0,200), pos=(400,152), size=(115,*), Val=100, Color=Black

ctl(34): Standard, "Tonal Width", Range=(0,127), pos=(400,168), size=(115,*), Val=0, Color=White
ctl(35): Standard, "Radius", Range=(0,200), Divisor=100, pos=(400,180), size=(115,*), Val=0, Color=White
ctl(36): Standard, "Amount", Range=(0,200), pos=(400,192), size=(115,*), Val=100, Color=White

ctl(37): checkbox(Pushlike), text="Enabled", Pos=(398,205), Size=(34,12), val=1
ctl(39): checkbox(Pushlike), text="Show", Pos=(434,205), Size=(34,12), val=0

ctl(14): checkbox, text="Show Clipping", Pos=(470,208), val=0, Action=Preview, Invisible // now there’s a histogram


ctl(38): OwnerDraw(drawitem), pos=(370,228), size=(256,100)


ctl(41): Ownerdraw(Modalframe), pos=(2,2), size=(448, 278), color=white, Anchor=ANCHOR_LEFT|ANCHOR_TOP, Invisible

ctl(42): StaticText, pos=(10,10), size=(430, 56), Anchor=ANCHOR_LEFT|ANCHOR_TOP, Invisible,
"High-radius low-amount (HiRaLoAm) unsharp masking is a technique for increasing local contrast, or ‘clarity.’  It involves first applying a Gaussian blur to a copy of the image, then subtracting the blurred image from the original.  The blur radius is large compared to that used for conventional sharpening, often an appreciable fraction of the overall image dimensions.  Unfortunately, the blur algorithms used in most software are necessarily approximations to a true mathematical Gaussian blur, and tend to depart from the ideal at large radii.  This limits the usefulness of the technique when editing large images, such as those with dimensions of 4000 pixels or more.  Even when high radii are available, the final result is often very different from the preview, and the intended effect may be impossible to achieve."

ctl(43): StaticText, pos=(10,66), size=(430, 72), Anchor=ANCHOR_LEFT|ANCHOR_TOP, Invisible,
"The Very High Radius Unsharp Mask filter takes advantage of the fact that since the blurred image contains no fine detail, it can be greatly resized with little effect on the final result.  In use, a copy of the original image is reduced in size, comparable to the image preview itself, and a similarly-scaled blur applied.  The blurred image is then enlarged to the size of the original image, and used as the unsharp mask.  The resulting image is very similar to the preview, even if many times its size.\n\n"

"- The Radius setting is a percent of the image diagonal, so it is independent of the actual image size.  The apparent and scaled radii in pixels are also displayed.  Lower values produce visible halos, which fade at higher values.  Sometimes there is a ‘sweet spot’ at low-to-middle values that enhances a certain scale of detail.\n"

ctl(44): StaticText, pos=(10,130), size=(430, 88), Anchor=ANCHOR_LEFT|ANCHOR_TOP, Invisible,
"- A large Amount setting helps you find a good Radius setting, but the intended effect is usually found at lower amounts.  Negative amounts reduce local contrast, an effect sometimes applied to portraits.\n"
"- You may save the Gaussian blurred image if desired; the only control that applies in this case is Radius.\n\n"

"- In Normal mode, the plugin processes each RGB channel separately.  The processed colors often have an unreal appearance but may still be interesting.  Color controls how much of the processed color components are included in the final image.  In Grayscale Only mode, the plugin processes a grayscale copy of the image, then combines it with the color components of the original image.  This is similar, but not identical, to Normal mode with Color set to zero; sometimes the results of the latter are superior.  However, the former runs faster.  In Color Only mode, only the processed color components are used in the final image.\n"

ctl(45): StaticText, pos=(10,210), size=(430, 64), Anchor=ANCHOR_LEFT|ANCHOR_TOP, Invisible, // frame H = Y + H + 4
"- Dark/Light Mixer works as follows:  the central position (zero) displays the processed image; the far left is like Darken blend mode, and the far right is like Lighten.\n"
"- Use the masks to protect blown shadows and highlights.  Try starting with Tonal Width and Radius 0 and Amount 100.  Raise the Tonal Width and Radius to create an appropriate mask, then reduce Amount to 0 and gradually raise it again to achieve the desired effect.  Amounts greater than 100 only affect blurred areas of the mask.\n"
"- The white histogram (background) is the original image; the black one (foreground) is the processed image.\n\n"
"Please visit russellcottrell.com."

ctl[CTL_PREVIEW]: Modify, Size=(333,295) // dialog H - 23

ctl[CTL_PROGRESS]: Modify, Pos=(4,306)
ctl[CTL_ZOOM]: Modify, Pos=(91,306)
ctl(50): Pushbutton, "100%", Pos=(160,303), Anchor=ANCHOR_LEFT|ANCHOR_BOTTOM, Size=(30,14), Action=Preview
ctl(51): Pushbutton, "Fit", Pos=(195,303), Anchor=ANCHOR_LEFT|ANCHOR_BOTTOM, Size=(30,14), Action=Preview

ctl(52): checkbox(Pushlike), text="Help", Pos=(435,301), Size=(34,14), val=0, Anchor=ANCHOR_RIGHT|ANCHOR_BOTTOM // dialog H - 17
ctl(53): Pushbutton, text="Reset", Pos=(474,301), Size=(34,14), val=0, Anchor=ANCHOR_RIGHT|ANCHOR_BOTTOM, Action=Preview
ctl[CTL_OK]: Modify, Pos=(513,301)
ctl[CTL_CANCEL]: Modify, Pos=(552,301) // dialog W - 38


OnFilterstart:
{

switch(imageMode)
{
  case RGBMode:
    sprintf(str0, "8-bit RGB Color");
    break;
  case RGB48Mode:
    sprintf(str0, "16-bit RGB Color");
    break;
  case GrayScaleMode:
    sprintf(str0, "8-bit Grayscale");
    break;
  case Gray16Mode:
    sprintf(str0, "16-bit Grayscale");
    break;
}
sprintf(str1, formatString("!t (%s, !f)"), str0);
setDialogText(str1);

if (e==FME_INIT) {
if ((imageMode == RGBMode) || (imageMode == RGB48Mode)) {
  if (!ctl(3)) {
    enableCtl(6, 3);
    enableCtl(10, 3); } }
else {
  enableCtl(6, 1);
  enableCtl(10, 1); }
j0 = getCtlRange(6, 1);
}

setCtlPos(38, -1, -1, PixelsToHDBUs(259), PixelsToVDBUs(103));

isTileable=false;

refreshWindow();

return false;
}

OnCtl(n):
{

  int jTemp;

  if (n==2 && ctl(2)==100)
    enableCtl(7,3);

  if (n==7 && e==FME_VALUECHANGED)
    setCtlRange(2, -40, 100 + ctl(7));

  if (n==10) {
    switch(ctl(10))
    {
    case 1:
      setCtlText(10, "Normal");
      if ((imageMode == RGBMode) || (imageMode == RGB48Mode))
        enableCtl(6, 3);
      jTemp = ctl(6);
      setCtlVal(6,j0);
      j0 = jTemp;
      break;
    case 2:
      setCtlText(10, "Grayscale Only");
      enableCtl(6, 1);
      break;
    case 0:
      setCtlText(10, "Color Only");
      enableCtl(6, 3);
      jTemp = ctl(6);
      setCtlVal(6,j0);
      j0 = jTemp;
      break;
    }
  }

  if (n==37)
    setCtlText(37, ctl(37) ? "Enabled" : "Off");
  else if (n==39)
    setCtlText(39, ctl(39) ? "Hide" : "Show");

if (ctl(3)) { // gaussian blur
  enableCtl(2, 1);
  enableCtl(6, 1);
  enableCtl(4, 1);
  enableCtl(31, 1);
  enableCtl(32, 1);
  enableCtl(33, 1);
  enableCtl(34, 1);
  enableCtl(35, 1);
  enableCtl(36, 1);
  enableCtl(37, 1);
}
else {
  enableCtl(2, 3);
  if (((imageMode == RGBMode) || (imageMode == RGB48Mode)) && (ctl(10) < 2))
    enableCtl(6, 3);
  enableCtl(4, 3);
  enableCtl(31, 3);
  enableCtl(32, 3);
  enableCtl(33, 3);
  enableCtl(34, 3);
  enableCtl(35, 3);
  enableCtl(36, 3);
  enableCtl(37, 3);
}


  if (e==FME_CLICKED)
  {
    if (ctl(52)) // Help
    {
      setDialogPos(false, -1, 0, -1, max(getDialogHeight(), 299)); // help frame H + 21
      for (i=41; i<=45; i++)
        enableCtl(i, 3);
    }
    else
    {
      for (i=41; i<=45; i++)
        enableCtl(i, 0);
      //Info("%d", getDialogHeight());
    }
  }

    if (ctl(53)) // Reset
    {
      setCtlVal(0,100);
      setCtlVal(2,30);
      setCtlRange(2, -40, 100);
      setCtlVal(7,0);
      setCtlVal(3,0);
      if ((imageMode == RGBMode) || (imageMode == RGB48Mode))
        setCtlVal(6,0);
      setCtlVal(10,1);
      setCtlVal(4,0);
      setCtlVal(31,0);
      setCtlVal(32,0);
      setCtlVal(33,100);
      setCtlVal(34,0);
      setCtlVal(35,0);
      setCtlVal(36,100);
      setCtlVal(37,1);
      setCtlVal(39,0);
      setCtlVal(14,0);
      setCtlVal(52,0);
      setCtlVal(53,0);
      triggerEvent(53, FME_CLICKED, 0);
      j0 = ctl(6);
      triggerEvent(10, FME_CLICKED, 0); // ***AFTER*** 53
      j0 = getCtlRange(6, 1);
      triggerEvent(37, FME_CLICKED, 0);
      triggerEvent(39, FME_CLICKED, 0);
    }

  if (e==FME_INIT)
    setZoom(-888);

  if (e==FME_SIZE)
    setZoom(-888);

  if (n==50)
    setZoom(1);
  else if (n==51)
    setZoom(-888);

return false;

}

ForEveryTile:
{

int C=imageMode<Gray16Mode?255:32768;
int C_2=imageMode<Gray16Mode?128:16384;
int C1=imageMode<Gray16Mode?1:65;
int C254=imageMode<Gray16Mode?254:32703;
double to8bit=imageMode<Gray16Mode?1.0:255.0/32768.0;

int zmax;
bool isRGB, RGBOn, RGBOff, clip;

int diag = sqr(X*X + Y*Y); // image diagonal
double radius = (double)(diag * ctl(0))/1000; // radius as percent of image size
double radiusS = (double)(diag * ctl(32))/10000;
double radiusH = (double)(diag * ctl(35))/10000;
int rc, rcS, rcH, minD;
int x_end_new, y_end_new;

double bg = 0.0;
double fg = 0.0;
double bg2 = 0.0;
double fg2 = 0.0;
double alpha = 1.0;
double alpha1 = (double)ctl(2)/100;
double alpha2 = fabs((double)ctl(4)/100);
double alpha3 = (ctl(10) < 2 ? (double)ctl(6)/100 : 0.0);
double alpha4 = 1.0;
double Yv, Cb, Cr, Cb0, Cb1, Cr0, Cr1, Rv, Gv, Bv;

double twS = (double)ctl(31)/to8bit;
double twH = (double)C - (double)ctl(34)/to8bit;

double amtS = (double)ctl(33)/100;
double amtH = (double)ctl(36)/100;

double p1,q1,fracx,fracy,CalcD;
int Xnew,Ynew,Xold,Yold,p,q,pp,qp,Calc;

double qr,b0,b1,b2,b3;

double b1_b0,b2_b0,b3_b0,BB;

double wn,wn1,wn2,wn3;

int rr, gg, bb, xx, n, xx0, xx1, hh0, hh1, hh2;
int nMax = 0;
int cWhite, cLGray, cGray, cDGray, cBlack;

minD = min(X/2, Y/2);

if (radius > (double)minD)
  radius = (double)minD;
rc = (int)ceil(radius);

if (radiusS > (double)minD)
  radiusS = (double)minD;
rcS = (int)ceil(radiusS);

if (radiusH > (double)minD)
  radiusH = (double)minD;
rcH = (int)ceil(radiusH);

if (doingProxy) { // save the preview dimensions
  k0 = X;
  k1 = Y;
  x0 = radius; }

else if (radius == 0) // can’t blur
  abort();

if (radius*(double)scaleFactor < 10)
  snprintf(str2, 3, "%f", radius*(double)scaleFactor);
else
  sprintf(str2, "%d", (int)(radius*(double)scaleFactor));
if (radius < 10)
  snprintf(str3, 3, "%f", radius);
else
  sprintf(str3, "%d", (int)radius);
strcat(str4, str2);
strcat(str4, " / ");
strcat(str4, str3);
strcat(str4, " px");
setCtlText(1, str4);
strcpy(str4, "");

// Array 99 is a virtual Y channel for the source, used in RGB mode, to minimize runtime calculations.
// In Grayscale Only mode, all changes target the luminosity component of the image to speed things up.
// The RGB Grayscale Only vs. color modes are indicated by the flags RGBOff and RGBOn.
// RGBOff means Grayscale Only mode, the destination 0 (red) channel is used to hold luminosity information only, and zmax = 1.
// RGBOn means one of the other modes, luminosity Y is calculated from all three destination color channels, and zmax = 3.
// The color components of the source image (calculated Cb and Cr) are added back at the very end.

allocArray(99, X, Y, 1, 2); // virtual Y channel

zmax = ((ctl(10)==2) || (imageMode == GrayScaleMode) || (imageMode == Gray16Mode)) ? 1 : 3; // grayscale only

isRGB = ( (imageMode == RGBMode) || (imageMode == RGB48Mode) ) ? true : false;
RGBOn = ( ((imageMode == RGBMode) || (imageMode == RGB48Mode)) && (zmax > 1) ) ? true : false; // RGB with color
RGBOff = ( ((imageMode == RGBMode) || (imageMode == RGB48Mode)) && (zmax == 1) ) ? true : false; // RGB without color

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{

if (isRGB) {
  Yv = (double)src(x,y,0)*0.299 + (double)src(x,y,1)*0.587 + (double)src(x,y,2)*0.114;
  putArray(99,x,y,0, (int)Yv); // saved src virtual Y channel
}

if (RGBOff)
  pset(x,y,0, getArray(99,x,y,0)); // turn the R channel into a Y channel

else {
  for (z=0; z<zmax; z++)
    pset(x,y,z, src(x,y,z));
}

}
}


////////// reduce //////////

if (!doingProxy) {

  Xnew = k0; // scale to the preview dimensions
  Ynew = k1;
  
  for (y=y_start; y<y_end; y++) { // rows

    //Calculate the coordinates in the old image for interpolation
    q1=(double)y*Y/Ynew;
    q=(int)(q1-0.5);

    //Calculate the fractional values
    fracy=q1-q;

  for (x=x_start; x<x_end; x++) { // columns

    //Calculate the coordinates in the old image for interpolation
    p1=(double)x*X/Xnew;
    p=(int)(p1-0.5);

    //Calculate the fractional values
    fracx=p1-p;

  for (z=0; z<zmax; z++) { // channels

    if (y<Ynew && x<Xnew) {

      if (fracx==0 && fracy==0) {
        //Nearest Neighbour
        Calc=pget(p,q,z);          // src

      } else if (fracx==0) {
        //Interpolation Formula
        CalcD= (1.0-fracy)*pget(p,q,z) + fracy*pget(p,q+1,z);  // src
        Calc= (int)CalcD;        

      } else if (fracy==0) {
        //Interpolation Formula
        CalcD=(1.0-fracx)*pget(p,q,z)  + fracx*pget(p+1,q,z);  // src
        Calc=(int)CalcD;
      } else {
        //Interpolation Formula
        CalcD=(1.0-fracx)*((1.0-fracy)*pget(p,q,z) +
          fracy*pget(p,q+1,z)) + fracx*((1.0-fracy)*pget(p+1,q,z) + fracy*pget(p+1,q+1,z)); // src
        Calc=(int)CalcD;
      }

    }

    //else {
      //set the rest of the image to black
    //  Calc=0;
    //}

    pset(x, y, z, Calc);

  }
  }
  if (updateProgress(z*Y+y, Y*zmax*8)) break;
  }

}

updateProgress(1, 8);



////////// blur //////////

if (!doingProxy) {
  radius = x0;
  x_end_new = x_start + k0; // just blur the small image
  y_end_new = y_start + k1;
}
else {
  Xnew = x_end;
  Ynew = y_end;
  x_end_new = x_end;
  y_end_new = y_end;
}

qr = radius;
b0 = 1.57825 + 2.44413*qr + 1.4281*qr*qr + 0.422205*qr*qr*qr;
b1 = 2.44413*qr + 2.85619*qr*qr + 1.26661*qr*qr*qr;
b2 = -(1.4281*qr*qr + 1.26661*qr*qr*qr );
b3 = 0.422205*qr*qr*qr;

b1_b0 = b1/b0;
b2_b0 = b2/b0;
b3_b0 = b3/b0;
BB = 1.0 - (b1_b0+b2_b0+b3_b0);

for (z=0; z<zmax; z++) // loop through all channels
{
for (y=y_start; y<y_end_new; y++) // loop through all rows
{

wn1 = wn2 = wn3 = (double)(pget(x_start+2*rc,y,z)); // src

// warmup forward - no write
for (x=x_start+2*rc; x>x_start; x--) 
{
wn = BB*(double)(pget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3); // src
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_start; x<x_end_new; x++) // real forward pass
{
wn = BB*(double)(pget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3); // src
tset(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_end_new-2; x>x_end_new-2-2*rc; x--) // chill down forward
{
wn = BB*(double)(pget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3); // src
t2set(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

// warmup backward - no write
for (x=x_end_new-1-2*rc; x<x_end_new-1; x++) 
{
wn = BB*(double)(t2get(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_end_new-1; x>=x_start; x--) // backward pass
{
wn = BB*(double)(tget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
pset(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

if (updateProgress(2*((z*Y+y)+(z*X)), (X+Y)*zmax*7)) break;

} // for y

// NOW DO COLUMNS

for (x=x_start; x<x_end_new; x++) // loop through all columns
{

// initialize column
wn1 = wn2 = wn3 = (double)(pget(x,y_start+2*rc,z)); 

// warmup forward pass
for (y=y_start+2*rc; y>y_start; y--)
{
wn = BB*(double)(pget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_start; y<y_end_new; y++) // real forward pass
{
wn = BB*(double)(pget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
tset(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_end_new-2; y>y_end_new-2-2*rc; y--) // chill down forward 
{
wn = BB*(double)(pget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
t2set(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

// warmup backward - no write
for (y=y_end_new-1-2*rc; y<y_end_new-1; y++) 
{
wn = BB*(double)(t2get(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_end_new-1; y>=y_start; y--) // backward pass
{
wn = BB*(double)(tget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
if (!doingProxy)
	t2set(x,y,z, (int)wn); // enlargement
else
	pset(x,y,z, (int)wn); // no enlargement

wn3=wn2; wn2=wn1; wn1=wn;
} // for y

if (updateProgress(2*(((z+1)*Y)+(z*X+x)), (X+Y)*zmax*7)) break;

} // for x

} // for z



////////// enlarge //////////

if (!doingProxy) {

  Xold = Xnew;
  Yold = Ynew;

  Xnew = (X*X)/(Xnew); // enlarge to original size
  Ynew = (Y*Y)/(Ynew);
  
  for (y=y_start; y<y_end; y++) { // rows

    //Calculate the coordinates in the old image for interpolation
    q1=(double)y*Y/Ynew;
    q=(int)(q1-0.5);

    //Calculate the fractional values
    fracy=q1-q;

  for (x=x_start; x<x_end; x++) { // columns

    //Calculate the coordinates in the old image for interpolation
    p1=(double)x*X/Xnew;
    p=(int)(p1-0.5);

    //Calculate the fractional values
    fracx=p1-p;

    pp = min(p+1,Xold-1); // prevent anti-aliasing at the right and bottom edges
    qp = min(q+1,Yold-1);

  for (z=0; z<zmax; z++) { // channels

    if (y<Ynew && x<Xnew) {

      if (fracx==0 && fracy==0) {
        //Nearest Neighbour
        Calc=t2get(p,q,z);

      } else if (fracx==0) {
        //Interpolation Formula
        CalcD= (1.0-fracy)*t2get(p,q,z) + fracy*t2get(p,qp,z) ;
        Calc= (int)CalcD;        

      } else if (fracy==0) {
        //Interpolation Formula
        CalcD= (1.0-fracx)*t2get(p,q,z)  + fracx*t2get(pp,q,z) ;
        Calc= (int)CalcD;
      } else {
        //Interpolation Formula
        CalcD=(1.0-fracx)*((1.0-fracy)*t2get(p,q,z) +
          fracy*t2get(p,qp,z)) + fracx*((1.0-fracy)*t2get(pp,q,z) + fracy*t2get(pp,qp,z));
        Calc= (int)CalcD;
      }

    }

    //else {
      //set the rest of the image to black
    //  Calc=0;
    //}

    pset(x, y, z, Calc);

  }
  }
  if (updateProgress(Y*zmax*3 + z*Y+y, Y*zmax*8)) break;
  }

}



////////// unsharp mask //////////

if (!ctl(3))
{

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{
for (z=0; z<zmax; z++)
{
  if (RGBOff)
    bg = (double)getArray(99,x,y,0);
  else
    bg = (double)src(x,y,z);
  fg = bg + alpha1*(bg - (double)pget(x,y,z));
  pset(x,y,z, (int)fg);
}
}
}

}

else if (RGBOff) {

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{
  pset(x,y,1, pget(x,y,0));
  pset(x,y,2, pget(x,y,0));
}
}

}

updateProgress(5, 8);



////////// darken/lighten mixer //////////

if ((ctl(4) != 0) && !ctl(3))
{

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{

if (isRGB)
  bg = (double)getArray(99,x,y,0);
else
  bg = (double)src(x,y,0);

if (RGBOn)
  fg = (double)pget(x,y,0)*0.299 + (double)pget(x,y,1)*0.587 + (double)pget(x,y,2)*0.114;
else
  fg = (double)pget(x,y,0);

for (z=0; z<zmax; z++)
{
  if (((ctl(4) < 0) && (bg < fg)) || ((ctl(4) > 0) && (bg > fg))) {
    if (RGBOff)
      bg2 = (double)getArray(99,x,y,0)/C;
    else
      bg2 = (double)src(x,y,z)/C;
    fg2 = (double)pget(x,y,z)/C;
    fg2 = alpha2*bg2 + (1-alpha2)*fg2;
    pset(x,y,z, (int)(fg2*C)); }
}

}
}

}

updateProgress(6, 8);



if ((ctl(37) || ctl(39)) && (ctl(31) + ctl(33) + ctl(34) + ctl(36) > 0) && !ctl(3)) { // begin mask



////////// shadow mask //////////

allocArray(0, X, Y, 1, 2); // the mask array

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{
  if (isRGB)
    fg = (double)pget(x,y,0)*0.299 + (double)pget(x,y,1)*0.587 + (double)pget(x,y,2)*0.114;
  else
    fg = (double)pget(x,y,0);

  if (fg < twS)
    putArray(0, x, y, 0, 0);
  else
    putArray(0, x, y, 0, C);
}
}



if (rcS > 0) { // can blur

////////// shadow blur //////////

qr = radiusS;
b0 = 1.57825 + 2.44413*qr + 1.4281*qr*qr + 0.422205*qr*qr*qr;
b1 = 2.44413*qr + 2.85619*qr*qr + 1.26661*qr*qr*qr;
b2 = -(1.4281*qr*qr + 1.26661*qr*qr*qr );
b3 = 0.422205*qr*qr*qr;

b1_b0 = b1/b0;
b2_b0 = b2/b0;
b3_b0 = b3/b0;
BB = 1.0 - (b1_b0+b2_b0+b3_b0);

z = 0;

for (y=y_start; y<y_end; y++) // loop through all rows
{

// initialize row
wn1 = wn2 = wn3 = (double)(getArray(0,x_start+2*rcS,y,z)); 

// warmup forward - no write
for (x=x_start+2*rcS; x>x_start; x--) 
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_start; x<x_end; x++) // real forward pass
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
tset(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_end-2; x>x_end-2-2*rcS; x--) // chill down forward
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
t2set(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

// warmup backward - no write
for (x=x_end-1-2*rcS; x<x_end-1; x++) 
{
wn = BB*(double)(t2get(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_end-1; x>=x_start; x--) // backward pass
{
wn = BB*(double)(tget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
putArray(0,x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

if (updateProgress(0,1)) break;

} // for y

// NOW DO COLUMNS

for (x=x_start; x<x_end; x++) // loop through all columns
{

// initialize column
wn1 = wn2 = wn3 = (double)(getArray(0,x,y_start+2*rcS,z));

// warmup forward pass
for (y=y_start+2*rcS; y>y_start; y--)
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_start; y<y_end; y++) // real forward pass
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
tset(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_end-2; y>y_end-2-2*rcS; y--) // chill down forward 
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
t2set(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

// warmup backward - no write
for (y=y_end-1-2*rcS; y<y_end-1; y++) 
{
wn = BB*(double)(t2get(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_end-1; y>=y_start; y--) // backward pass
{
wn = BB*(double)(tget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
putArray(0,x,y,z, (int)wn);

wn3=wn2; wn2=wn1; wn1=wn;
} // for y

if (updateProgress(0,1)) break;

} // for x

} // can blur



////////// apply shadow mask //////////

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{

  alpha4 = 1.0 - (1 - (double)getArray(0, x, y, 0)/C) * amtS;
  if (alpha4 < 0)
    alpha4 = 0.0;

if (!ctl(39)) {
for (z=0; z<zmax; z++)
{
  if (RGBOff)
    bg = (double)getArray(99,x,y,0)/C;
  else
    bg = (double)src(x,y,z)/C;
  fg = (double)pget(x,y,z)/C;
  fg = alpha4*fg + (1-alpha4)*bg;
  pset(x,y,z, (int)(fg*C));
}
}

else // to show the mask
  t3set(x,y,0, (int)(alpha4*C));

}
}

freeArray(0);



////////// highlight mask //////////

allocArray(0, X, Y, 1, 2); // the mask array

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{
  if (isRGB)
    fg = (double)pget(x,y,0)*0.299 + (double)pget(x,y,1)*0.587 + (double)pget(x,y,2)*0.114;
  else
    fg = (double)pget(x,y,0);

  if (fg > twH)
    putArray(0, x, y, 0, 0);
  else
    putArray(0, x, y, 0, C);
}
}



if (rcH > 0) { // can blur

////////// highlight blur //////////

qr = radiusH;
b0 = 1.57825 + 2.44413*qr + 1.4281*qr*qr + 0.422205*qr*qr*qr;
b1 = 2.44413*qr + 2.85619*qr*qr + 1.26661*qr*qr*qr;
b2 = -(1.4281*qr*qr + 1.26661*qr*qr*qr );
b3 = 0.422205*qr*qr*qr;

b1_b0 = b1/b0;
b2_b0 = b2/b0;
b3_b0 = b3/b0;
BB = 1.0 - (b1_b0+b2_b0+b3_b0);

z = 0;

for (y=y_start; y<y_end; y++) // loop through all rows
{

// initialize row
wn1 = wn2 = wn3 = (double)(getArray(0,x_start+2*rcH,y,z)); 

// warmup forward - no write
for (x=x_start+2*rcH; x>x_start; x--) 
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_start; x<x_end; x++) // real forward pass
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
tset(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_end-2; x>x_end-2-2*rcH; x--) // chill down forward
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
t2set(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

// warmup backward - no write
for (x=x_end-1-2*rcH; x<x_end-1; x++) 
{
wn = BB*(double)(t2get(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

for (x=x_end-1; x>=x_start; x--) // backward pass
{
wn = BB*(double)(tget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
putArray(0,x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for x

if (updateProgress(0,1)) break;

} // for y

// NOW DO COLUMNS

for (x=x_start; x<x_end; x++) // loop through all columns
{

// initialize column
wn1 = wn2 = wn3 = (double)(getArray(0,x,y_start+2*rcH,z));

// warmup forward pass
for (y=y_start+2*rcH; y>y_start; y--)
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_start; y<y_end; y++) // real forward pass
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
tset(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_end-2; y>y_end-2-2*rcH; y--) // chill down forward 
{
wn = BB*(double)(getArray(0,x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
t2set(x,y,z, (int)wn);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

// warmup backward - no write
for (y=y_end-1-2*rcH; y<y_end-1; y++) 
{
wn = BB*(double)(t2get(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
wn3=wn2; wn2=wn1; wn1=wn;
} // for y

for (y=y_end-1; y>=y_start; y--) // backward pass
{
wn = BB*(double)(tget(x,y,z)) + (b1_b0*wn1 + b2_b0*wn2 + b3_b0*wn3);
putArray(0,x,y,z, (int)wn);

wn3=wn2; wn2=wn1; wn1=wn;
} // for y

if (updateProgress(0,1)) break;

} // for x

} // can blur



////////// apply highlight mask //////////

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{

  alpha4 = 1.0 - (1 - (double)getArray(0, x, y, 0)/C) * amtH;
  if (alpha4 < 0)
    alpha4 = 0.0;

if (!ctl(39)) {
for (z=0; z<zmax; z++)
{
  if (RGBOff)
    bg = (double)getArray(99,x,y,0)/C;
  else
    bg = (double)src(x,y,z)/C;
  fg = (double)pget(x,y,z)/C;
  fg = alpha4*fg + (1-alpha4)*bg;
  pset(x,y,z, (int)(fg*C));
}
}

else { // to show the mask

  gg = (int)(alpha4*C);

  if (isRGB) {

    if (gg < t3get(x,y,0)) { // blend the highlight and shadow masks
      rr = C;
      bb = C - abs(t3get(x,y,0) - gg); }
    else if (t3get(x,y,0) < gg) {
      rr = C - abs(t3get(x,y,0) - gg);
      bb = C; }
    else {
      rr = C;
      bb = C; }

    gg = min(t3get(x,y,0), gg);

    pset(x,y,0,rr);
    pset(x,y,1,gg);
    pset(x,y,2,bb); }

  else
    pset(x,y,0, min(t3get(x,y,0), gg));

}

}
}

freeArray(0);



} // end mask



////////// color //////////

if (!ctl(3) && !ctl(39))
{

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{

clip = false;

// clipping, based on the fully processed image because of faint color tints

if (ctl(14) && doingProxy)
{

if (isRGB) // isRGB
{

if ( (RGBOff && (pget(x,y,0)<C1) ) || ( RGBOn && ((pget(x,y,0)<C1) && (pget(x,y,1)<C1) && (pget(x,y,2)<C1))) )
{
  pset(x,y,0, 0);
  pset(x,y,1, 0);
  pset(x,y,2, C);
  clip = true;
}

else if ( (RGBOff && (pget(x,y,0)>C254)) || (RGBOn && ((pget(x,y,0)>C254) && (pget(x,y,1)>C254) && (pget(x,y,2)>C254))) )
{
  pset(x,y,0, C);
  pset(x,y,1, 0);
  pset(x,y,2, 0);
  clip = true;
}

} // isRGB


else // grayscale
{

if (pget(x,y,0)<C1)
{
  pset(x,y,0, C);
}

else if (pget(x,y,0)>C254)
{
  pset(x,y,0, 0);
}

} // grayscale

}


if (isRGB && !clip)
{

  if (ctl(10)==1) // normal
    Yv = (double)pget(x,y,0)*0.299 + (double)pget(x,y,1)*0.587 + (double)pget(x,y,2)*0.114;
  else if (ctl(10)==2) // grayscale only
    Yv = (double)pget(x,y,0);
  else if (ctl(10)==0) // color only
    Yv = (double)getArray(99,x,y,0);

  if (alpha3==0.0) { // color components of the original image
    Cb = -(double)src(x,y,0)*0.168736 - (double)src(x,y,1)*0.331264 + (double)src(x,y,2)*0.500;
    Cr =  (double)src(x,y,0)*0.500 - (double)src(x,y,1)*0.418688 - (double)src(x,y,2)*0.081312; }
  else {
    Cb0 = -(double)src(x,y,0)*0.168736 - (double)src(x,y,1)*0.331264 + (double)src(x,y,2)*0.500;
    Cr0 =  (double)src(x,y,0)*0.500 - (double)src(x,y,1)*0.418688 - (double)src(x,y,2)*0.081312;
    Cb1 = -(double)pget(x,y,0)*0.168736 - (double)pget(x,y,1)*0.331264 + (double)pget(x,y,2)*0.500;
    Cr1 =  (double)pget(x,y,0)*0.500 - (double)pget(x,y,1)*0.418688 - (double)pget(x,y,2)*0.081312;
    Cb = alpha3*Cb1 + (1-alpha3)*Cb0;
    Cr = alpha3*Cr1 + (1-alpha3)*Cr0; }

  Rv = Yv + 1.402*Cr;
  Gv = Yv - 0.34414*Cb - 0.71414*Cr;
  Bv = Yv + 1.772*Cb;
  pset(x,y,0, (int)Rv);
  pset(x,y,1, (int)Gv);
  pset(x,y,2, (int)Bv);

}


}
}

}

updateProgress(8, 8);



////////// histograms //////////

allocArray(0, 2, 256, 0, 4);
fillArray (0, 0);

cWhite = RGB(255,255,255);
cLGray = RGB(200,200,200);
cGray = RGB(180,180,180);
cDGray = RGB(80,80,80);
cBlack = RGB(0,0,0);

for (x=x_start; x<x_end; x++)
{
for (y=y_start; y<y_end; y++)
{

  if (isRGB)
    Yv = (double)getArray(99,x,y,0);
  else
    Yv = (double)src(x,y,0);
  Yv *= to8bit; // convert to 8-bit
  xx = (int)Yv;
  n = getArray(0,0,xx,0) + 1;
  putArray(0,0,xx,0, n);
  if (n > nMax)
      nMax = n;

  if (!ctl(39)) { // show mask
  if (isRGB)
    Yv = (double)pget(x,y,0)*0.299 + (double)pget(x,y,1)*0.587 + (double)pget(x,y,2)*0.114;
  else
    Yv = (double)pget(x,y,0);
  Yv *= to8bit;
  xx = (int)Yv;
  n = getArray(0,1,xx,0) + 1;
  putArray(0,1,xx,0, n);
  if (n > nMax)
      nMax = n; }

}
}

if (!ctl(52)) {

if (nMax > (X*Y*4)/256)
  nMax = (X*Y*4)/256;

startSetPixel(38);

setRectFrame(0, 0, 258, 102, cGray);
setRectFill(1, 1, 257, 101, cLGray);

if (nMax > 0) {
for (i=0; i<256; i++) {

  xx0 = i+1;
  xx1 = i+1+1;
  hh0 = 100-min((int)((double)getArray(0,0,i,0)*(100.0/(double)nMax)),100)+1;
  hh1 = 100-min((int)((double)getArray(0,1,i,0)*(100.0/(double)nMax)),100)+1;
  hh2 = 100+1;

  if (hh0 < hh1) {
    setRectFill(xx0, hh0, xx1, hh1, cWhite);
    setRectFill(xx0, hh1, xx1, hh2, cDGray); }
  else {
    setRectFill(xx0, hh1, xx1, hh0, cBlack);
    setRectFill(xx0, hh0, xx1, hh2, cDGray); }
}
}

endSetPixel(38);

}


return true;
}


OnFilterEnd:
{
  updateProgress(0,1);
  return false;
}