comp.lang.idl-pvwave archive
Messages from Usenet group comp.lang.idl-pvwave, compiled by Paulo Penteado

Home » Public Forums » archive » Re: mean() function
Show: Today's Messages :: Show Polls :: Message Navigator
E-mail to friend 
Switch to threaded view of this topic Create a new topic Submit Reply
Re: mean() function [message #46900] Thu, 12 January 2006 03:07 Go to next message
Maarten[1] is currently offline  Maarten[1]
Messages: 176
Registered: November 2005
Senior Member
And the reason you need that page, is in part because IDL uses the
moment routine described in Numerical Recipes (take total first, divide
later), instead of a proper running average, like the GNU scientific
library does.

However, since looping is slow in IDL, you don't want to implement that
in IDL, so the next best thing is to have that page.

Maarten
Re: mean() function [message #46924 is a reply to message #46900] Wed, 11 January 2006 06:19 Go to previous messageGo to next message
Paolo Grigis is currently offline  Paolo Grigis
Messages: 171
Registered: December 2003
Senior Member
>
> [...]
>
> I suppose math libraries or compilers changed between 6.1 and 6.2.
>
> Once the values you are adding differ by 6-7 orders of magnitude,
> precision is completely lost for single precision floats.
>
> IDL> print, total(replicate(1.0, 10^7))
> 1.00000e+07
> IDL> print, total(replicate(1.0, 10^8))
> 1.67772e+07
> IDL> print, total(replicate(1.0D0, 10^8))
> 1.0000000e+08

I have no problem with that, understanding the argument of finite
precision and all that, but there is one thing which has me puzzling
(all commands issued from the same linux machine)


IDL Version 5.4 (linux x86). (c) 2000, Research Systems, Inc.

IDL> help,total(replicate(1.,1d8))
<Expression> FLOAT = 1.67772e+07
IDL> exit


IDL Version 5.5a (linux x86). (c) 2001, Research Systems, Inc.

IDL> help,total(replicate(1.,1d8))
<Expression> FLOAT = 6.71089e+07
IDL> exit


IDL Version 5.6 (linux x86 m32). (c) 2002, Research Systems, Inc.

IDL> help,total(replicate(1.,1d8))
<Expression> FLOAT = 1.00000e+08
IDL> exit

I was thinking that the result of the operation would depend
on the hardware used, but I would have guessed that on the same
machine no difference would be seen between different versions
of IDL, since all of them should represent floats in the same way...

So what's happening here? Different compiler optimizations?

Ciao,
Paolo


>
> Cheers, Ken
Re: mean() function [message #46927 is a reply to message #46924] Wed, 11 January 2006 01:49 Go to previous messageGo to next message
Nigel Wade is currently offline  Nigel Wade
Messages: 286
Registered: March 1998
Senior Member
biocpu@yahoo.com wrote:

> { mipseb IRIX unix IRIX 6.0 Jun 27 2003 64 64}
>
>
> Thanks everybody!

IDL> print,!version
{ mipseb IRIX unix IRIX 6.0 Jun 27 2003 32 64}
IDL>
IDL> y = fltarr(1008879)+35
IDL> print, mean(y)
% Compiled module: MEAN.
% Compiled module: MOMENT.
35.0000

The only difference I see is that this binary is the n32 version whereas yours
looks like the 64bit version. Maybe that's the difference.

Do you have the n32 version you can test? I don't have the 64bit version
installed, and I can't locate it for download.

--
Nigel Wade, System Administrator, Space Plasma Physics Group,
University of Leicester, Leicester, LE1 7RH, UK
E-mail : nmw@ion.le.ac.uk
Phone : +44 (0)116 2523548, Fax : +44 (0)116 2523555
Re: mean() function [message #46929 is a reply to message #46927] Tue, 10 January 2006 19:22 Go to previous messageGo to next message
Kenneth P. Bowman is currently offline  Kenneth P. Bowman
Messages: 585
Registered: May 2000
Senior Member
In article <dq1gv6$v05$1@news.nems.noaa.gov>,
Paul Van Delst <Paul.vanDelst@noaa.gov> wrote:

> Kenneth Bowman wrote:
>> In article <1136932449.216202.42760@g44g2000cwa.googlegroups.com>,
>> biocpu@yahoo.com wrote:
>>
>>
>>> y = fltarr(1008879)+35
>>
>>
>> Looks like roundoff error to me
>>
>> IDL> y = fltarr(1008879)+35
>> IDL> print, mean(y)
>> 35.0497
>> IDL> print, mean(y, /double)
>> 35.000000
>> IDL> y = dblarr(1008879)+35
>> IDL> print, mean(y)
>> 35.000000
>
> Huh. I don't see this in single precision (see other post). What version of
> IDL did you use?
>
> paulv

The version I posted (quoted above) is

{ ppc darwin unix Mac OS X 6.2 Jun 20 2005 32 32}

If I run it on my PowerBook (now sadly obsolete ;-) ), which is running
6.1

{ ppc darwin unix Mac OS X 6.1 Jul 14 2004 32 32}

I get exactly what was in the original post (he was running 6.0 on IRIX)

IDL> y = fltarr(1008879) + 35
IDL> print, mean(y)
35.5249
IDL> print, mean(y, /double)
35.000000
IDL> print, version

I suppose math libraries or compilers changed between 6.1 and 6.2.

Once the values you are adding differ by 6-7 orders of magnitude,
precision is completely lost for single precision floats.

IDL> print, total(replicate(1.0, 10^7))
1.00000e+07
IDL> print, total(replicate(1.0, 10^8))
1.67772e+07
IDL> print, total(replicate(1.0D0, 10^8))
1.0000000e+08

Cheers, Ken
Re: mean() function [message #46930 is a reply to message #46929] Tue, 10 January 2006 16:18 Go to previous messageGo to next message
biocpu is currently offline  biocpu
Messages: 3
Registered: June 2004
Junior Member
{ mipseb IRIX unix IRIX 6.0 Jun 27 2003 64 64}


Thanks everybody!
Re: mean() function [message #46933 is a reply to message #46930] Tue, 10 January 2006 15:13 Go to previous messageGo to next message
Paul Van Delst[1] is currently offline  Paul Van Delst[1]
Messages: 1157
Registered: April 2002
Senior Member
Kenneth Bowman wrote:
> In article <1136932449.216202.42760@g44g2000cwa.googlegroups.com>,
> biocpu@yahoo.com wrote:
>
>
>> y = fltarr(1008879)+35
>
>
> Looks like roundoff error to me
>
> IDL> y = fltarr(1008879)+35
> IDL> print, mean(y)
> 35.0497
> IDL> print, mean(y, /double)
> 35.000000
> IDL> y = dblarr(1008879)+35
> IDL> print, mean(y)
> 35.000000

Huh. I don't see this in single precision (see other post). What version of IDL did you use?

paulv

--
Paul van Delst
CIMSS @ NOAA/NCEP/EMC
Re: mean() function [message #46934 is a reply to message #46933] Tue, 10 January 2006 15:09 Go to previous messageGo to next message
Paul Van Delst[1] is currently offline  Paul Van Delst[1]
Messages: 1157
Registered: April 2002
Senior Member
F�LDY Lajos wrote:
> Hi,
>
> have you tried mean(y,/double)?

Even regular old single worked for me:

IDL> y = fltarr(1008879)+35
IDL> id = where(y ne 35, cc)
IDL> print, cc
0
IDL> print, mean(y)
% Compiled module: MEAN.
% Compiled module: MOMENT.
35.0000
IDL> print, mean(y[0:400000])
35.0000

Maybe it's version related? I noticed that biocpu used () rather than [] for array
indexing. Maybe an earlier version of IDL had a MEAN() function that didn't use a
compensated summation algorithm? The OP mean from 0->400000 that worked suggests that's
not the case, but who knows? Anyway....

IDL> print, !version
{ x86 linux unix linux 6.0.3 Feb 26 2004 32 64}


paulv

>
> regards,
> lajos
>
>
> On Tue, 10 Jan 2006 biocpu@yahoo.com wrote:
>
>
>> The following looks very odd. Have any clues?
>>
>> IDL> y = fltarr(1008879)+35
>> IDL> id = where(y ne 35, cc)
>> IDL> print, cc
>> 0
>> ; so y is strictly 35.0 BUT
>> IDL> print, mean(y)
>> 35.5249
>>
>> IDL> print, mean(y(0:400000))
>> 35.0000
>>
>> Thanks,
>>
>>


--
Paul van Delst
CIMSS @ NOAA/NCEP/EMC
Re: mean() function [message #46935 is a reply to message #46934] Tue, 10 January 2006 15:07 Go to previous messageGo to next message
David Fanning is currently offline  David Fanning
Messages: 11724
Registered: August 2001
Senior Member
biocpu@yahoo.com writes:

> The following looks very odd. Have any clues?
>
> IDL> y = fltarr(1008879)+35
> IDL> id = where(y ne 35, cc)
> IDL> print, cc
> 0
> ; so y is strictly 35.0 BUT
> IDL> print, mean(y)
> 35.5249
>
> IDL> print, mean(y(0:400000))
> 35.0000

Goodness! It is the 10th already and this is the first
"The Sky is Falling" post of the year. I was beginning
to wonder what was going on! :-)

You might want to have a look at this:

http://www.dfanning.com/math_tips/sky_is_falling.html

Cheers,

David
--
David Fanning, Ph.D.
Fanning Software Consulting, Inc.
Coyote's Guide to IDL Programming: http://www.dfanning.com/
Re: mean() function [message #46936 is a reply to message #46935] Tue, 10 January 2006 15:00 Go to previous messageGo to next message
K. Bowman is currently offline  K. Bowman
Messages: 330
Registered: May 2000
Senior Member
In article <1136932449.216202.42760@g44g2000cwa.googlegroups.com>,
biocpu@yahoo.com wrote:

> y = fltarr(1008879)+35

Looks like roundoff error to me

IDL> y = fltarr(1008879)+35
IDL> print, mean(y)
35.0497
IDL> print, mean(y, /double)
35.000000
IDL> y = dblarr(1008879)+35
IDL> print, mean(y)
35.000000


Ken Bowman
Re: mean() function [message #46937 is a reply to message #46936] Tue, 10 January 2006 14:53 Go to previous messageGo to next message
Foldy Lajos is currently offline  Foldy Lajos
Messages: 268
Registered: October 2001
Senior Member
Hi,

have you tried mean(y,/double)?

regards,
lajos


On Tue, 10 Jan 2006 biocpu@yahoo.com wrote:

> The following looks very odd. Have any clues?
>
> IDL> y = fltarr(1008879)+35
> IDL> id = where(y ne 35, cc)
> IDL> print, cc
> 0
> ; so y is strictly 35.0 BUT
> IDL> print, mean(y)
> 35.5249
>
> IDL> print, mean(y(0:400000))
> 35.0000
>
> Thanks,
>
>
Re: mean() function [message #46957 is a reply to message #46935] Sun, 15 January 2006 00:58 Go to previous message
R.Bauer is currently offline  R.Bauer
Messages: 1424
Registered: November 1998
Senior Member
David Fanning wrote:

> biocpu@yahoo.com writes:
>
>> The following looks very odd. Have any clues?
>>
>> IDL> y = fltarr(1008879)+35
>> IDL> id = where(y ne 35, cc)
>> IDL> print, cc
>> 0
>> ; so y is strictly 35.0 BUT
>> IDL> print, mean(y)
>> 35.5249
>>
>> IDL> print, mean(y(0:400000))
>> 35.0000
>
> Goodness! It is the 10th already and this is the first
> "The Sky is Falling" post of the year. I was beginning
> to wonder what was going on! :-)
>
> You might want to have a look at this:
>
> http://www.dfanning.com/math_tips/sky_is_falling.html
>
> Cheers,
>
> David

mean is not useable if it results in this

IDL> print,mean( make_array(500000,val=35,/float) )
35.0413
IDL> print,mean( make_array(400000,val=35,/float) )
35.0000

I prefer a slower routine if this is right.

no one would accept 1.0 + 1.0 result = 1.5

The article is good to understand why mean should not be resolved this way!

Does one have implemented the gnu's version?

cheers
Reimar
Re: mean() function [message #46992 is a reply to message #46900] Thu, 12 January 2006 08:49 Go to previous message
savoie is currently offline  savoie
Messages: 68
Registered: September 1996
Member
"Maarten" <maarten.sneep@knmi.nl> writes:

> And the reason you need that page, is in part because IDL uses the
> moment routine described in Numerical Recipes (take total first, divide
> later), instead of a proper running average, like the GNU scientific
> library does.
>
> However, since looping is slow in IDL, you don't want to implement that
> in IDL, so the next best thing is to have that page.


Maarten,

Would you mind explaining this a bit for me? What's a proper running
average? And why is it better in general?

Thanks
Matt

--
Matthew Savoie - Scientific Programmer
National Snow and Ice Data Center
(303) 735-0785 http://nsidc.org
  Switch to threaded view of this topic Create a new topic Submit Reply
Previous Topic: Re: reading a ninary file
Next Topic: mean() function

-=] Back to Top [=-
[ Syndicate this forum (XML) ] [ RSS ] [ PDF ]

Current Time: Wed Oct 08 19:04:58 PDT 2025

Total time taken to generate the page: 0.00687 seconds