digitalFAQ.com Forums [Archives]

digitalFAQ.com Forums [Archives] (http://www.digitalfaq.com/archives/)
-   Avisynth Scripting (http://www.digitalfaq.com/archives/avisynth/)
-   -   Avisynth: Wavelet denoising? (http://www.digitalfaq.com/archives/avisynth/1519-avisynth-wavelet-denoising.html)

GFR 11-08-2002 01:20 PM

Avisynth: Wavelet denoising?
 
After someone (blackprince?) posted a link about it I studied some references about wavelet noise reduction and it seems really promising.

The idea is the following (very simplified):

At each level of a Wavelet transform, you split the image in two components: a "blured" component (low-pass) and a "detail" component (high-pass).

If you look at the "detail" component, most of it should be zero, except where you have, well, details: edges, textures. When you have noise, instead of mostly zero, you have random small amplitude points (the noise) with the higher amplitude edges and textures.

If you "gate" this "detail" component so that these small points due to noise are zeroed and the higher amplitude edges and textures are kept untouched, and then combine it with the "blured" component, the noise is killed but the edges and textures are as sharp as before!

This is unlike traditional smoothing where all high frequencies are attenuated, and by the time noise is killed the edges and textures are also smoothed...

Besides removing normal "white" noise, it is also impressive at removing block noise if your source is a low bitrate AVI or MPEG file :)

One nice feature is that you can try to estimate an "optimal" threshold automatically from the variance of the "detail" component.

I have programmed an implementation of it in Delphi (that works with still images), I'll try to translate it to C and compile as a AVISynth filter. (the original filter documentation is in Japanese and it won't work with my old Pentium II :(, it needs SSE).

kwag 11-08-2002 04:14 PM

Hi GFR,

Is there a reference model or a definition in C, C++ ?
I'm interested in looking at that 8)

-kwag

GFR 11-11-2002 06:45 AM

You can find some C/C++ or Matlab libraries on the web (try it in Google).

While there are many different wavelets you can try (each one with its pros and cons) the wavelet I'm trying to code is a very simple algorythm called "lifting".

It's like this:

1) Split the image in two sub-images, one with the even colums, the other with the odd columns (the first column is column 0).

2) Make the odd image = odd - even. The odd image is now a high pass in the horizontal direction (note: it can be negative).

3) Make the even image = even + odd div 2 This makes the even image a low pass in the horizontal direction.

4) Split each of the two subimages into even and odd lines, and now calculate odd-odd as a high pass in horizontal and vertical directions ("corners" detail), odd-even is high pass in H and low pass in V (vertical edges), even-odd is low pass H high pass V (horizontal edges) and the even-even is a lower resolution version of the original image.

5) If you want you can run the even-even image through the algorythm again for another level of the transform and so on.

You can estimate the variance of the noise using the highest resolution odd-odd image, because it is problably mostly noise. The median of the aboslute deviation to the median is said to be a good estimative.

Now "gate" all the detail images with
x=0 if abs(x)<thresh
x=x-sign(x)*thresh if abs(x)>=thresh

Or any other gating function you want (this one avoids creation of false oscillations).

Use the thresh=k*(variance of the noise), k from 2 to 3 should be good enough (the optimal k is subject of many researches). You can also try a different thresh for each subimage.

The inverse transform is straightforward, just make even=even-odd div 2, odd = odd+even and interleave the even with odd.

Good things:
1) very fast, few operations, can be all integer.
2) no need for extra memory, the calculations can be in-place (you don't actually split the images just keep track of the pointers)

Bad things:
it's not the smoothest wavelet around and it is not the one that best decorrelates the low freq coeffs from the high freq ones.

I hope the Delphi code is working today, I can post it to you if you want it.

kwag 11-11-2002 07:41 AM

Hi GFR,

Sure!, post it. I would be nice if it can be ported to C. My Pascal (Dephi) programming is very (VERY) rusty, but I should be able to read the code and do an un-optimized working port to C, and then worry about optimizations. It would be very interesting if this would complement SansGrip "Blockbuster" filter :wink:

-kwag

GFR 11-11-2002 10:42 AM

Just the forward transform - I'm debugging the backward transform.
It seems long but it's repetitive and can problably be optimized.
Nothing special, only you have to be careful so that you know which index is pointing to what!

Code:

procedure TForm1.Button1Click(Sender: TObject);
var
  BitMap1,BitMap2:TBitMap;
  i,j: Integer;
  tempint,tempint2:Integer;
  Line1,Line2:PByteArray; // Delphi pointer to an array of Bytes
  sigma,lambda,t: real;
begin
  BitMap1:=Image1.Picture.Bitmap;
  BitMap2:=Image4.Picture.Bitmap;

  // I copy Image1 to (invisible) Image4
  // so I have the original Image to compare to
  // You can do it all in place
 
  For i:=0 to BitMap1.Height-1 do
    begin
      Line1:=BitMap1.ScanLine[i]; // Line1 points to line i of Image1
      Line2:=BitMap2.ScanLine[i]; // Line2 points to line i of Image4
      For j:=0 to BitMap1.Width-1 do
        Line2[j]:=Line1[j]; // Get each pixel i,j
    end;

  {Take even,odd COLUMNS
  odds:=err (odd-even)
  even:=mean (even+1/2*err)
  }
  for j:=0 to BitMap2.Height-1 do // Scan the j lines
    begin
      Line2:=BitMap2.ScanLine[j];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          tempint2:=Max(Min(Line2[2*i+1]-Line2[2*i],127),-127);
          // column 2*i+1 are the odd columns
          // I clip the values to the range -127..127
          tempint:=Line2[2*i] + (tempint2 div 2);
          // column 2*i are even columns
          // I used tempint and tempint2 because the Bitmaps only handle unsigned
          // bytes in the range 0..255
          {high-pass - odds}
          Line2[2*i+1]:=128+tempint2;
          // Add a DC level 128 so that it's in the appropriate range
          // And you can store it in the bitmap and preview it
          {passa-baixas - pares}
          Line2[2*i]:=tempint;
        end;
    end;

  {draws the bitmap}
  // I use Image2 to draw because Image4 is "scrambled"
  // This draws a low pass, squashed picture on the left and the high pass on the right
  {        L | H        }
  Bitmap1:=Image2.Picture.Bitmap;
  for j:=0 to BitMap2.Height-1 do
    begin
      Line2:=Bitmap2.ScanLine[j];
      Line1:=Bitmap1.ScanLine[j];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i]:=Line2[2*i];
          Line1[i+(Image1.Width div 2)]:=Line2[2*i+1];
        end;
    end;

  Image2.Invalidate; // forces the Image to redraw


  {Take even,odd LINES
  odds:=err (odd-even)
  even:=mean (even+1/2*err)
  }
  for j:=0 to (BitMap2.Height div 2)-1 do // scan the lines
    begin
      Line1:=BitMap2.ScanLine[2*j]; // even lines
      Line2:=BitMap2.ScanLine[2*j+1]; // odd lines
      for i:=0 to BitMap2.Width-1 do
        begin
          tempint2:=Max(Min(Line2[i]-Line1[i],127),-127);
          tempint:=Line1[i] + (tempint2 div 2);
          {passa-altas - impares}
          Line2[i]:=128+tempint2;
          {passa-baixas - pares}
          Line1[i]:=tempint;
        end;
    end;

  {draws}
  {don't care about sigma it's an estimate of the variance but it's not correct}
  {the picture is arranged like this
 
    LL | HL
    ---+---
    LH | HH
   
  }
  {You don't need to draw this but it is very useful for debugging
  and for learning about the images and wavelets}
  sigma:=0;
  Bitmap1:=Image2.Picture.Bitmap;
  for j:=0 to (BitMap2.Height div 2)-1 do
    begin
      Line2:=Bitmap2.ScanLine[2*j];
      Line1:=Bitmap1.ScanLine[j];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i]:=Line2[2*i];
          Line1[i+(Image1.Width div 2)]:=Line2[2*i+1];
        end;
      Line2:=Bitmap2.ScanLine[2*j+1];
      Line1:=Bitmap1.ScanLine[j+(BitMap2.Height div 2)];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i]:=Line2[2*i];
          Line1[i+(Image1.Width div 2)]:=Line2[2*i+1];
          sigma:=sigma+{sqr(Line2[2*i+1]-128)}abs(Line2[2*i+1]-128); {MAD?}
        end;
    end;

  sigma:=(sigma/((BitMap2.Width/2)*(BitMap2.Height/2)))/0.6745{sqrt(sigma/((BitMap2.Width/2)*(BitMap2.Height/2)))};
  lambda:=sqrt(2*Log10((BitMap2.Width)*(BitMap2.Height)));
  t:=sigma*lambda;
  // write to the screen to debug the variance estimate
  Label2.Caption:=FloatToStr(sigma);
  Edit1.Text:=FloatToStr(lambda);
  Label5.Caption:=FloatToStr(t);

  Image2.Invalidate;

  {The nest level or SCALE}

  {Takes COLUMNS *4 ("evens"), *4+2 ("odds")
  odds:=err(odds-evens)
  evens:=mean (ebens+1/2*err)
  }
  for j:=0 to (BitMap2.Height div 2)-1 do // scan 1 column out of 4
    begin
      Line2:=BitMap2.ScanLine[2*j];
      for i:=0 to (BitMap2.Width div 4)-1 do
        begin
          tempint2:=Max(Min(Line2[4*i+2]-Line2[4*i],127),-127);
          tempint:=Line2[4*i] + (tempint2 div 2);
          {passa-altas - impares}
          Line2[4*i+2]:=128+tempint2;
          {passa-baixas - pares}
          Line2[4*i]:=tempint;
        end;
    end;

  {draws}
  {    LLL | LLH |  HL     
          |    |
    -------+-----+-----------
                |
          LH    |  HH
                |                }               
  Bitmap1:=Image2.Picture.Bitmap;
  for j:=0 to (BitMap2.Height div 2)-1 do
    begin
      Line2:=Bitmap2.ScanLine[2*j];
      Line1:=Bitmap1.ScanLine[j];
      for i:=0 to (BitMap2.Width div 4)-1 do
        begin
          Line1[i]:=Line2[4*i];
          Line1[i+(Image1.Width div 4)]:=Line2[4*i+2];
        end;
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i+(Image1.Width div 2)]:=Line2[2*i+1];
        end;
      Line2:=Bitmap2.ScanLine[2*j+1];
      Line1:=Bitmap1.ScanLine[j+(BitMap2.Height div 2)];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i]:=Line2[2*i];
          Line1[i+(Image1.Width div 2)]:=Line2[2*i+1];
        end;
    end;

  Image2.Invalidate;

  {Takes LINES *4 ("evens"), *4+2 ("odds")
  odds:=err(odds-evens)
  evens:=mean (ebens+1/2*err)
  }
  for j:=0 to (BitMap2.Height div 4)-1 do // each line out of 4
    begin
      Line1:=BitMap2.ScanLine[4*j];
      Line2:=BitMap2.ScanLine[4*j+2];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          tempint2:=Max(Min(Line2[2*i]-Line1[2*i],127),-127);
          tempint:=Line1[2*i] + (tempint2 div 2);
          {passa-altas - impares}
          Line2[2*i]:=128+tempint2;
          {passa-baixas - pares}
          Line1[2*i]:=tempint;
        end;
    end;

  {draws}
  {    LLLL | LLHL |  HL     
    --------+------+       
      LLLH | LLHH |
    --------+------+-----------
                  |
          LH      |  HH
                  |                }               
  Bitmap1:=Image2.Picture.Bitmap;
  for j:=0 to (BitMap2.Height div 4)-1 do
    begin
      Line2:=Bitmap2.ScanLine[4*j];
      Line1:=Bitmap1.ScanLine[j];
      for i:=0 to (BitMap2.Width div 4)-1 do
        begin
          Line1[i]:=Line2[4*i];
          Line1[i+(Image1.Width div 4)]:=Line2[4*i+2];
        end;
      Line2:=Bitmap2.ScanLine[4*j+2];
      Line1:=Bitmap1.ScanLine[j+(BitMap2.Height div 4)];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i]:=Line2[4*i];
          Line1[i+(Image1.Width div 4)]:=Line2[4*i+2];
        end;
    end;
  for j:=0 to (BitMap2.Height div 2)-1 do
    begin
      Line2:=Bitmap2.ScanLine[2*j];
      Line1:=Bitmap1.ScanLine[j];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i+(Image1.Width div 2)]:=Line2[2*i+1];
        end;
      Line2:=Bitmap2.ScanLine[2*j+1];
      Line1:=Bitmap1.ScanLine[j+(BitMap2.Height div 2)];
      for i:=0 to (BitMap2.Width div 2)-1 do
        begin
          Line1[i]:=Line2[2*i];
          Line1[i+(Image1.Width div 2)]:=Line2[2*i+1];
        end;
    end;

  Image2.Invalidate;

  // enables the threshold setting
  TrackBar1.Enabled:=true;
  // enable the gating and inverse transform buttons
  Button2.Enabled:=true;
  Button4.Enabled:=false;
end;


GFR 11-12-2002 09:44 AM

Here are some tests I ran:

http://www.geocities.com/gfr.geo/wavelets.html

GFR 11-12-2002 01:19 PM

Hi kwag

Below:

http://www.cs.kuleuven.ac.be/˜wavelets/

you can find a C++ wavelet library and some nice articles.

I'm afraid you need some background in signal/image processing to understand these.

kwag 11-12-2002 01:49 PM

Quote:

Originally Posted by GFR
Hi kwag

Below:

http://www.cs.kuleuven.ac.be/˜wavelets/

you can find a C++ wavelet library and some nice articles.

I'm afraid you need some background in signal/image processing to understand these.

Hi GFR,

You're right :) I can read and write code, but I've never gone into signal processing ( DSP's, etc. ) algorithm's that deep :roll:
To many things I want to do, not enough time :D
The most I've done in digital coding is generating a Golay (23,12) encoder in software, using a 6522 PIA (Peripheral Interface Adapter ), for test encoding on paging equipment ( circa 1985 8O ) . Some info related to what I had to deal with here: http://www.math.uic.edu/~fields/Deco...roduction.html
And that was ~15 years ago 8O . After that, my programming has been in the communications field ( more of I/O control, some embedded stuff ,etc ) but not on video. So I have a lot to learn, before I can start applying my programming concepts in this field ;)

-kwag

SansGrip 11-13-2002 02:12 PM

Re: wavelet denoising
 
Quote:

Originally Posted by GFR
After someone (blackprince?) posted a link about it I studied some references about wavelet noise reduction and it seems really promising.

This does sound very interesting, and I'm intrigued about the possibility of including it (or a variation) in NoMoSmooth. Could you post the references?

black prince 11-13-2002 04:41 PM

Hi SansGrip,

SansGrip wrote:
Quote:

SansGrip Posted: Wed Nov 13, 2002 3:12 pm Post subject: Re: wavelet denoising

--------------------------------------------------------------------------------

GFR wrote:
After someone (blackprince?) posted a link about it I studied some references about wavelet noise reduction and it seems really promising.


This does sound very interesting, and I'm intrigued about the possibility of including it (or a variation) in NoMoSmooth. Could you post the references?
Follow this link for the refrences to wavelet denoising filter. I hope this
helps:


http://www.kvcd.net/forum/viewtopic.php?t=1440

-black prince

GFR 11-14-2002 05:19 AM

References:

http://www.cs.kuleuven.ac.be/˜wavelets/

Some publications, a C++ package.

You can find the articles below with google:

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 9, NO. 9, SEPTEMBER 2000
Adaptive Wavelet Thresholding for Image Denoising
and Compression
S. Grace Chang, Student Member, IEEE, Bin Yu, Senior Member, IEEE, and Martin Vetterli, Fellow, IEEE

On Denoising and Best Signal Representation
Hamid Krim, Senior Member, IEEE, Dewey Tucker, St´ephane Mallat, Member, IEEE, and David Donoho

DE-NOISING BY
SOFT-THRESHOLDING
David L. Donoho
Department of Statistics
Stanford University

Multimedia Applications
of the Wavelet Transform
(PHD thesis)
Claudia Kerstin Schremmer

Maarten Jansen
pieflab - matlab package, some articles

GFR 11-14-2002 06:53 AM

Copy this and paste in the browser's address bar (without http:// )

www.geocities.com/gfr.geo/smooth-lift.gif

This looks horrible :( but it is very clear to illustrate the concept.

I'm using the integer lifting wavelet (lazy wavelet) , then CLIPPING the high level details instead of killing the low-level details. Since this particular wavelet produces blocks when you eliminate details, it's very easy to see which portions of the the picture are "smoothed" (in fact "blocked") and which are kept intact.

With a better choice of wavelet this can be a good selective edge blur filter.

SansGrip 11-15-2002 11:43 AM

Experimental Avisynth wavelet filter
 
I threw together a quick filter using the algorithm posted above. I hope I got it right (it seems to reconstruct the image fine, but I'd be grateful if someone could check the source against the description of the algorithm). Currently it only does one pass, and only operates on luma.

You can specify a threshold for each of the detail components. See the readme for usage.

Doesn't seem very effective to me, with low thresholds not really touching the noise and higher ones making the image blocky.

Is this because it only does one pass? How many passes would be good to be effective against noise? Or is it just not a very good type of wavelet? If so, what would be more suitable for denoising?

Incidentally, what would happen if one were to apply this to even/odd frames instead of pixels?

kwag 11-15-2002 01:28 PM

Re: Experimental Avisynth wavelet filter
 
Quote:

Originally Posted by SansGrip
I threw together a quick filter using the algorithm posted above. I hope I got it right (it seems to reconstruct the image fine, but I'd be grateful if someone could check the source against the description of the algorithm). Currently it only does one pass, and only operates on luma.

You can specify a threshold for each of the detail components. See the readme for usage.

Doesn't seem very effective to me, with low thresholds not really touching the noise and higher ones making the image blocky.

Is this because it only does one pass? How many passes would be good to be effective against noise? Or is it just not a very good type of wavelet? If so, what would be more suitable for denoising?

Incidentally, what would happen if one were to apply this to even/odd frames instead of pixels?

Hi SansGrip,

I tested the wavelet filter. But as you said, the lower the values, the lower the filtering. But I still can see the artifacts. If I increase the values, they just get blured, but they are still present. I believe because it is applied on the input stage of the encoder, it won't be as usefull as if it was applied on a output stage of an encoder. Which in this case ( TMPEG, CCE, etc. ) , it's impossible to do. If the encoders had a "Pre filter" hook stage and a "Post filter" hook stage to plug filters, then I believe it would work on the output stage, just as it does on the still images samples presented earlier on this thread. Just my thought :roll:

-kwag

SansGrip 11-15-2002 04:28 PM

Re: Experimental Avisynth wavelet filter
 
@kwag, GFR, all:

Quote:

Originally Posted by kwag
I tested the wavelet filter. But as you said, the lower the values, the lower the filtering. But I still can see the artifacts. If I increase the values, they just get blured, but they are still present.

Edit: The most likely reason that the output from this filter is so poor is because of the "lifting" approach, which basically involves halving the resolution and then reconstructing the entire image from that. "Proper" wavelets should perform much better, but of course are slower and less easy to code.

Theoretically speaking wavelets should be a very good technique for denoising. There are several key features of noise that I've identified:

1) Low amplitude
2) High frequency
3) Random

The key features of the details we want to keep are:

1) Low amplitude
2) Mostly low frequency
3) Non-random

Because dumb spatial/temporal softeners usually only pay attention to the first feature of noise (the characteristic low amplitude), they usually annihilate details too (because they're also low amplitude).

In order for a filter to remove noise and nothing else, it needs to answer the following questions:

Low amplitude? No -> Skip
Low frequency? Yes -> Skip
Random? No -> Skip
Replace with (possibly weighted) average

Only if all tests pass should the pixel be considered influenced by noise.

Step 1 is easy, because low amplitude just means "how much does the actual value vary from the average?" If a lot, it's high amplitude. If very little, it's low.

Step 2 is tricky in that it involves some kind of transformation of the signal. Generally to perform frequency analysis on a signal one must transform it into the frequency domain, i.e. using Fourier techniques. These are a pain to code and slow to run.

The nice thing about wavelets is that they let you (with varying effectiveness based on the particular wavelet used) isolate the high-frequency and low-frequency parts of the signal without jumping through too many hoops. They're both (usually) easier to code and faster than Fourier methods.

Step 3 is either very easy or very hard, I've not really thought about it enough yet ;). Randomness is the defining factor of noise, though, and it's becoming my opinion that any truly non-destructive denoiser must measure the (temporal) randomness of what it thinks is noise. Simply being within a variance threshold is not sufficient.

Another complicating factor is motion. To properly account for motion it's necessary to use motion compensation (aka motion estimation), which is a very tricky thing to code. Motion adaption -- as used in NoMoSmooth -- is a step in the right direction, though it might be rendered redundant if the three steps outlined above are effective.

Anyway, that's a summary of my thoughts over the last few days trying to improve NoMoSmooth. Differentiating between noise and detail is very hard, but wavelets might open an avenue to solving the puzzle.

I'd be very interested to hear about various wavelets and their strengths/weaknesses, preferably without too many mathematical details, at least initially ;).

SansGrip 11-15-2002 04:37 PM

Re: wavelet denoising
 
Quote:

Originally Posted by GFR
If you "gate" this "detail" component so that these small points due to noise are zeroed and the higher amplitude edges and textures are kept untouched, and then combine it with the "blured" component, the noise is killed but the edges and textures are as sharp as before!

Note that the gating must be done only on low-amplitude high-frequency parts of signal. If we could find a wavelet that produces only this component we'd have a much easier time denoising, since we wouldn't have to worry about those low-amplitude low-frequency details that we want to keep.

snowmoon 11-15-2002 10:53 PM

Couldn't a FFT filter be used to get rid of those?

SansGrip 11-16-2002 11:35 AM

Quote:

Originally Posted by snowmoon
Couldn't a FFT filter be used to get rid of those?

FFT is a very effective way of targeting specific frequencies, but it's also a little tricky to code and slower than wavelets (since you have to transform from the spatial domain to the frequency domain, do your processing, then transform back again). That said I'm looking into FFT as a possibility.

kwag 11-16-2002 12:14 PM

Hi SansGrip,

This is not related to wavelets, but it might be worth a look.
I found this in the file readpic.c, after studying some sources originally by the MPEG Group, and modified by John Schlichther as used in his program "AVI2MPG1".

Here's the code:

/*
Settings for the softfilter: Maybe these should be adjustable from the
commandline, but Tolearance 10 and Filtersize 6 are good defaults.

- Tolerance: An absolute value (on a scale from 0 to 255), which
determines, if the colour diffence should be exposed to filtering or
not. This value needs to be bigger than the image noise peak to peak
value. All small stuctures in the image, that have a contrast of less
that this value are filtered down and might get lost.

- Filtersize: The array used for every pixel for filtering. The larger
this value is, the stronger the filter is, but also the longer will
filtering take. doubling this value will take four time longer for
filtering.
*/

const int Tolerance = 10;
const int Filtersize = 6;

#define MAX(a, b) ((a)>(b)?(a):(b))
#define MIN(a, b) ((a)<(b)?(a):(b))

/*
Softfilter is an advanced blur filter, that will not blindly soften everything
in the image, but will look at all pixels around the pixel to be modified, and
will only use these pixels, that have almost the same value.
By this simple rule, it will not soften down sharp edges, but only adjacent
pixels of almost the same color. The effect of this is dramatic. It can
eliminate almost all noise in large single colored areas of the image, without
significantly decreasing sharpness in other areas, where there is high
contrast.
The code is not really highly optimised, but it produces good results. If you
need speed badly, you can either not use it, or use some inline assembler.
I assume, you should be able to get quite a speedup from using assembler, but
I guess there is not more than 20% speedup possible by using C techniques.
Thomas Hieber (thieber@gmx.net)
*/

void softfilter(unsigned char* dest, unsigned char* src, int width, int height)
{
int refval, aktval, upperval, lowerval, numvalues, sum, rowindex;
int x, y, fx, fy, fy1, fy2, fx1, fx2;

for(y = 0; y < height; y++)
{
for(x = 0; x < width; x++)
{
refval = src[x+y*width];
upperval = refval + Tolerance;
lowerval = refval - Tolerance;

numvalues = 1;
sum = refval;

fy1 = MAX(y-Filtersize, 0);
fy2 = MIN(y+Filtersize+1, height);

for (fy = fy1; fy<fy2; fy++)
{
rowindex = fy*width;
fx1 = MAX(x-Filtersize, 0) + rowindex;
fx2 = MIN(x+Filtersize+1, width) + rowindex;

for (fx = fx1; fx<fx2; fx++)
{
aktval = src[fx];
if ((lowerval-aktval)*(upperval-aktval)<0)
{
numvalues ++;
sum += aktval;
}
} //for fx
} //for fy

dest[x+y*width] = sum/numvalues;
}
}
}


Regards,
-kwag

SansGrip 11-16-2002 01:48 PM

Quote:

Originally Posted by kwag
This is not related to wavelets, but it might be worth a look.

NoMoSmooth already does this (as do TemporalSoften and SpatialSoften). This technique is especially useful in temporal smoothing because it completely eliminates ghosting.

Unfortunately it's also significantly slower than just doing a straight average because you have to check each pixel within the radius to see if it's inside the threshold.

GFR 11-18-2002 07:35 AM

Hi

The "lazy" wavelet from the code I've posted is not good for this application - when you remove the details, the reconstruction gets blocky.

Sansgrip, you definetely need another stage to get usable results, the highest frequency bands are of too high frequency, you should have very little energy in these bands (with a good wavelet) so that they make very little difference unless you are too agressive to them :)

Some papers I've read state that you get a significant improvement in the smoothness of the reconstructed image if you use redundant wavelets, that is, if each band is not critically sampled and you use several "shifts" for each band. This of course involves LOTS more computation and memory.

I think we could try an intermediary thing like this:

1) smooth the image so you get rid of the high frequencies (but dont subsample the blured image, keep it at full resolution). You could try any smoothing filter you want, a moving average filter, splines, etc.

You could use separable filters (filter horizontal, then vertical) or non separable filters (that work with a "radius" around a given pixel).

2) Generate a detail image that is the difference of the original image and the blured image. This image is not subsampled too.

3) You can go on and smooth the smoothed image even more, with a more drastic version of your smoothing filter, and generate another difference image.

4) Apply the thresholding to the difference images (low level gate for noise removing, clipping for selective edge smoothing)

5) and then reconstruct the image adding the most smoothed image to its difference, then to the previous difference, etc.

You can easily change the smoothing filter without changing the rest of the program, and try different filters until you're satisfied with the smoothness. I think you could even implement the full algorythm in a AVISynth script :)

Even if maybe it's not exactly a Wavelet (depends on the filter you use) it will be a subband decomposition with perfect reconstruction, and the fact that we don't subsample the subbands may help the visual quality of the denoised image.

I think that eliminating the upsampling + filtering stage is the key to get rid of most of the spurious oscillations on the denoised image.

SansGrip 11-18-2002 10:09 AM

Quote:

Originally Posted by GFR
I think we could try an intermediary thing like this

This process you described (smoothing, threshold the difference, recompose) is basically the same as:

1) Find average of pixel values in n-by-n block
2) If centre pixel is out by more than threshold, average it with neighours

which is pretty much what all smoothers do already, and has the drawback that fine details are smoothed away along with noise. I don't think smoothing the smoothed (or even the difference) image again would help, because noise is generally so low-amplitude that it would be taken away in the first pass.

I'm intruiged by (non-lazy) wavelets, though, and would love to find one specifically targeted at denoising. Targeted at video denoising would be even better, but I can dream ;).

So if you hear of any wavelet algorithms (even very computationally expensive ones) you think would be particularly suitable, it would be interesting to code one of them up with your mathematical help to see how effective it really is.

GFR 11-18-2002 04:27 PM

1) The problem: we have a noisy image and we want to recover the best approximation of the clean image we can, from this noisy image.

2) The clean image (the signal) is not available but must be guessed from the noisy image. We can make some assumptions about the clean image.

The first is that each pixel in the image is highly correlated to its neighbours. This means that the information is spread all over the image.

The second assumption we can make is that we can apply a "transform" to the image, that is, represent the image in an alternate way such that most of the information is "compacted" in a few coefficients, and the other, many, coefficients are not so important. One of such transforms is the DCT, another is the Wavelet Transform (WT). If you look at the DCT of an image or of a block of an image you'll see that most of the energy is compacted near the DC coefficient and it rapidly vanishes as you move to the highest frequencies. In a WT most of the energy is at the lowest resolution subband, and each detail band has less energy as you move to the highest frequency subbands.

This is because most images have most of their energy concentrated in the lowest frequencies, and have very little high frequency information. This is, of course, a generalization, but this is what lossy image compressing schemes like JPEG or MPEG rely on.

3) A practical denoising technique must also make some assumptions on the noise.

In the so called white noise one pixel of noise is not correlated to other pixels of noise; This means it doesn't have a DC level (zero mean) and that if you take a DCT or WT there won't be any concentration of significant coefficients - the white noise has equal energy at all frequencies. This doesn't hold if the noise is "coloured".

White noise is also not related to the clean signal in any way. Some kinds of noise, like JPEG artifacts, are highly correlated to the clean signal.

The noise is considered additive, this meaning that the noisy image is a "sum" of the clean image with some amount of pure noise.

To be able to make a nice approximation of the clean image using the noisy image, the amount of noise should not be "too much". It's not only much harder to guess what's signal in the middle of an ocean of noise than to remove a little noise in an almost good image, but the more noise you have the worse will be your results. As the level of noise increases, the more subtle information is masked by the noise and you can only extract the more strong features of the clean signal.

4) The Wavelet denoising.

Since the noise is supposed to be additive, and the WT is linear, the WT of the noisy image is equal to the WT of the clean signal added to the WT of the pure noise.

The WT of the clean image has most energy compacted at the low frequency subbands. The WT of the noise has its energy spread evenly all over the subbands, and its energy level is hopefully very much smaller than the energy level of the lowest frequencies of the clean signal and hopefully smaller than the most important signal features at the highest frequency bands of the clean signal.

Since the energy level of the noise is hopefully very much smaller than the energy level of the lowest frequencies of the clean signal, you can leave the low frequencies alone as the noise is masked by the signal and wont bother you :)

Traditional denoising kills the noise by attenuating the high frequencies (with some smoothing filter), but this also attenuates the high frequencies of the signal by the same amount.

The WT thresholding technique uses a "kill or keep" strategy. If the high frequency coefficient amplitude is below a certain level, it is noise or some detail you can't see anyway because of the noise; so you kill it. If it is above a certain threshold it problably ain't noise, but it is some detail you can see even with the noise; so you keep it. The result is that the noise is "killed" but the details that stand above the noise level are kept intact. The higher the noise level, the more details you have to sacrifice, just like in traditional denoising using "smoothing", but the difference is that even if you're sacrificing lots of details, the most proeminent details are still kept intact.

================================================== ====================================

GFR 11-19-2002 06:39 AM

Hi,

take a look at

http://www.geocities.com/gfr.geo/teste2.html

Looks a lot better than the lazy Wavelet!.


here's a nice tutorial (but you need some math):

http://cas.ensmp.fr/~chaplais/Waveto...tation_US.html

This one is easier to understand, very good:

http://engineering.rowan.edu/~polikar/homepage.html

In this thesis you've got a chapter about still images and another chapter on video encoding:

http://www.informatik.uni-mannheim.d...mmer2002b.html
Multimedia Applications of the Wavelet Transform
Claudia Kerstin Schremmer

There are many papers on still images denoising you can find with google.

This book:

Signal and Image processing with neural networks (A C++ sourcebook)
Timothy Masters
John Wiley & Son Inc

Has the most intuitive wavelet explanation I've read, and also has got C code for image manipulation...

Some interesting stuff:

http://www.informatik.uni-mannheim.d...mmer2001d.html

http://www.informatik.uni-mannheim.d...mmer2001g.html

kwag 11-19-2002 08:27 AM

Hi GFR,

Now that's a lot of good information :D
Some of it, I can understand, and some I need an ice bag on my head after reading. 8O
Excelent references.
This one is really nice and graphical: http://engineering.rowan.edu/%7epoli...Ttutorial.html

Thanks,
-kwag

GFR 11-19-2002 11:53 AM

More to torture your neurons :)

I've found a C source for Wavelets & images.

Note: this Mallat guy is THE guru for wavelets!

Quote:

From: Stephane Mallat
Date: Tue, 22 Dec 92 12:56:39 -0500
Subject: Wavelet Softwares Availables

/************************************************** **************/
/* Wave1 Software */
/* Authors: Emmanuel Bacry, Wen-Liang Hwang, */
/* Stephane Mallat and Sifen Zhong */
/* */
/* (c) Copyright 1992 New York University */
/* All rights reserved */
/************************************************** **************/


/************************************************** **************/
/* Wave2 Software */
/* Authors: Wen-Liang Hwang, Stephane Mallat, Sifen Zhong */
/* */
/* (c) Copyright 1992 New York University */
/* All rights reserved */
/************************************************** **************/

Wave1 and Wave2 are two softwares written in C
to process signals with the wavelet transform.
These software run under Unix and X11 windows.
They provide shell environements that include
image displays, plotting and postscript facilities,
as well as general signal processing subroutines.
They can be particularly useful for researchers or engineers who
want to make experiments with the wavelet tranform.
The wavelet transform algorithms implemented in these softwares are
mainly based on two published papers:

[1] "Characterization of signals from multiscale edges" by Stephane Mallat
and Sifen Zhong,
IEEE Transactions on Pattern Analysis and Machine Intelligence,
Vol. 14, No. 7, p. 710-732, July 1992,

[2] "Singularity detection and processing with wavelets", by Stephane Mallat
and Wen Liang Hwang, IEEE Transactions on Information Theory,
Vol. 38, No. 2, p. 617-643, March 1992.

The Wave1 software processes one-dimensional signals whereas Wave2
processes images. Dyadic wavelet tranform decompositions and reconstruction
are implemented. Signal singularities and edges are detected from the
local maxima of the wavelet transform. We implemented the
subroutines that reconstruct signals from the local maxima of their wavelet
transform, with the algorithm described in [1]. One-dimensional signals
as well as images can thus be processed through the local
maxima of their wavelet transform (multiscale edges). Many subroutines
to process these local maxima are available. For images, applications to noise
removal and optical flow detection are implemented. An on-line
documentation is available.

A nonexclusive license to use this software and its documentation
for scholarly and research purposes only is hereby granted by
New York University. There is no fee to use this software for these
purposes. In using the software and documentation for scholarly and
research purposes you may copy or modify them but (a) you must include
the NYU copyright notice and the original authors' names on each copy
or modification you make of the software and of the supporting documentation;
and (b) you may not use, rewrite, or adapt the software or the documentation
(or any part of either) as the basis for any commercial software or hardware
product, without in each instance obtaining the prior written consent of
NYU or an appropriate license from NYU. You may not use this software
and its documentation (or any part of either) for any other purpose
without obtaining an appropriate license from NYU. To obtain any
license or consent from NYU, please contact: Patrick Franc, Office of
Industrial Liaison, New York University, 251 Mercer Street, New York,
NY 10012.

The source codes of Wave1 and Wave2 are available via anonymous ftp on the
machine cs.nyu.edu (IP number 128.122.140.24). After login with login name
"anonymous", do "cd pub/wave". The wave1 and wave2 softwares are
respectively in the compressed tarfiles "wave1.tar.Z" and "wave2.tar.Z".
These are unix-compressed tarfiles, so don't forget to put ftp under "binary"
mode before getting them. You must use the "uncompress" and "tar" unix
commands to restore the software.
Instructions for compilation are in the README files.

GFR 11-21-2002 11:05 AM

This is great!

http://www.informatik.uni-mannheim.d...d/animationen/

It has java applets for Wavelets, DCT, etc. very instructive!

GFR 01-08-2003 12:26 PM

http://perso.wanadoo.fr/polyvalens/clemens/clemens.html

The above site has an intuitive wavelet tutorial (downloadable in pdf) that covers noise reduction among other things, and it has pascal and c codes.

It also has a link to a page with a collection!!! of transforms :)

Gaudi 01-08-2003 11:16 PM

Don´t know much about wavelets ( :oops: in fact know nothing), but used Matlab and know it has Wavelets and even an Image Processing Toolbox.

Matlab Homepage

These guy at Matlab are very good, so perhaps going through the code of the toolbox functions will be a good start. As far as I know, the coda can be compiled to C++.


Hope it helps a little.


Gaudi


All times are GMT -5. The time now is 07:28 PM  —  vBulletin © Jelsoft Enterprises Ltd

Site design, images and content © 2002-2024 The Digital FAQ, www.digitalFAQ.com
Forum Software by vBulletin · Copyright © 2024 Jelsoft Enterprises Ltd.