| filtering problem [message #22539] |
Thu, 16 November 2000 00:00  |
Dave Brennan
Messages: 28
Registered: December 1998
|
Junior Member |
|
|
Hi,
i don't know if anyone can help but it's worth a try!
I am trying to filter an array say (256x256) with a window of size 65x65
which scans across the array pixel by pixel. It should compare the
statistics of the area within the kernal with the global statistics of
the image to produce a correction image. (This is a particular type of
inhomogeneity correction)
In detail: 'the algorithm should correct the pixel value by a
multiplicative factor found by dividing the global mean by the window
mean'
A further problem is I want the ability to set a threshold where data
below the threshold are not included in the statistics and not corrected
by the algorithm.
At first I though I could just use convol to produce a correction map
but this does not allow me to set a threshold.
Does anyone have any ideas? It needs to be as fast as possible as it
will work on 128 images at a time.
Cheers
Dave Brennan
|
|
|
|
| Re: filtering problem [message #22603 is a reply to message #22539] |
Tue, 21 November 2000 00:00  |
John-David T. Smith
Messages: 384
Registered: January 2000
|
Senior Member |
|
|
Dave Brennan wrote:
>
> Once again thanks for the help.
>
> However, life is not quite as simple as I imagined. Previously I said :
>
>> In detail: 'the algorithm should correct the pixel value by a
>> multiplicative factor found by dividing the global mean by the window
>> mean'
>> A further problem is I want the ability to set a threshold where data
>> below the threshold are not included in the statistics and not corrected
>> by the algorithm.
>
> I have later found that this is not entirely correct.
>
> Although any value below the thershold should not be included in the
> statistics, pixels below the threshold should be corrected by the algorithm.
>
> Therefore is it possible to change the code:
>
> pro thresh, a, n, t
> m = a ge t
> wh = where(m,cnt)
> if cnt eq 0 then return
> a[wh] = a[wh] * mean(a[wh]) * (smooth(float(m),n,/EDGE) / $
> (smooth(a*m,n,/EDGE)+1.e-30))[wh]
> end
You can simply remove the [wh] from everywhere but the mean. By the
way, I rarely use constructs like adding 1.e-30 (I took that from the
previous poster's version to ensure comparable runtime penalties). If
you were ever using this on data arrays which were quite small, you'd be
in trouble. Dealing with overflows is why the !VALUES sysvar was
invented. And, on the other hand, if your sub-threshold pixels
dominated, such that it may be impossible to correct certain locations
(no window mean defined), you should explicitly test for this, and do
something sensible (either relax your threshold, warn the user, etc.)
JD
--
J.D. Smith | WORK: (607) 255-6263
Cornell Dept. of Astronomy | (607) 255-5842
304 Space Sciences Bldg. | FAX: (607) 255-5875
Ithaca, NY 14853 |
|
|
|
|
| Re: filtering problem [message #22619 is a reply to message #22539] |
Tue, 21 November 2000 00:00 |
Dave Brennan
Messages: 28
Registered: December 1998
|
Junior Member |
|
|
Thanks jaco. You just hit the nail on the head. I forgot to remove the [wh]
after the smooth operation (changed a[wh] to a and couldn't work out why it
wasn't working!!!), therefore I was getting array size errors. There were a
few further lines I had to add to do exactly what I wanted but that seems
to have solved my problem. Thanks folks. Now on to the next part of the
code.......
Cheers
Dave Brennan
>
> If I understand you correctly, you want to correct the below-threshold
> pixels by the factor calculated from the valid pixels surrounding it.
> Couln't you just leave out some of the specific subscripting [wh] in the
> final assignment?
> a = a * mean(a[wh]) * (smooth(float(m),n,/EDGE) / $
> (smooth(a*m,n,/EDGE)+1.e-30))
> Correct me if I'm wrong here...
>
> Jaco
>
> ----------------
> Jaco van Gorkom e-mail: gorkom@rijnh.nl
> FOM-Instituut voor Plasmafysica "Rijnhuizen", The Netherlands
|
|
|
|
| Re: filtering problem [message #22621 is a reply to message #22539] |
Tue, 21 November 2000 00:00 |
Jaco van Gorkom
Messages: 97
Registered: November 2000
|
Member |
|
|
"Dave Brennan" <9147261b@clinmed.gla.ac.uk> wrote in message
news:3A1A445B.14A912A3@clinmed.gla.ac.uk...
...
> Although any value below the thershold should not be included in the
> statistics, pixels below the threshold should be corrected by the
algorithm.
>
> Therefore is it possible to change the code:
>
> pro thresh, a, n, t
> m = a ge t
> wh = where(m,cnt)
> if cnt eq 0 then return
> a[wh] = a[wh] * mean(a[wh]) * (smooth(float(m),n,/EDGE) / $
> (smooth(a*m,n,/EDGE)+1.e-30))[wh]
> end
>
> to allow this.
If I understand you correctly, you want to correct the below-threshold
pixels by the factor calculated from the valid pixels surrounding it.
Couln't you just leave out some of the specific subscripting [wh] in the
final assignment?
a = a * mean(a[wh]) * (smooth(float(m),n,/EDGE) / $
(smooth(a*m,n,/EDGE)+1.e-30))
Correct me if I'm wrong here...
Jaco
----------------
Jaco van Gorkom e-mail: gorkom@rijnh.nl
FOM-Instituut voor Plasmafysica "Rijnhuizen", The Netherlands
|
|
|
|
| Re: filtering problem [message #22624 is a reply to message #22539] |
Tue, 21 November 2000 00:00 |
Dave Brennan
Messages: 28
Registered: December 1998
|
Junior Member |
|
|
Once again thanks for the help.
However, life is not quite as simple as I imagined. Previously I said :
> In detail: 'the algorithm should correct the pixel value by a
> multiplicative factor found by dividing the global mean by the window
> mean'
> A further problem is I want the ability to set a threshold where data
> below the threshold are not included in the statistics and not corrected
> by the algorithm.
I have later found that this is not entirely correct.
Although any value below the thershold should not be included in the
statistics, pixels below the threshold should be corrected by the algorithm.
Therefore is it possible to change the code:
pro thresh, a, n, t
m = a ge t
wh = where(m,cnt)
if cnt eq 0 then return
a[wh] = a[wh] * mean(a[wh]) * (smooth(float(m),n,/EDGE) / $
(smooth(a*m,n,/EDGE)+1.e-30))[wh]
end
to allow this.
After scratching my head for a while I can't think of a simple solution without
going back to my clunky method where I manually scan through the data and
threshold the window each time. There must be a better way!!
Any help would be greatly appreciated!!
Cheers
Dave Brennan
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
