Cumulative standard deviation: Difference between revisions
Content added Content deleted
(→Old style, four ways: Agh! More mistypes! And cut&pasted too.) |
|||
Line 1,118: | Line 1,118: | ||
=={{header|Fortran}}== |
=={{header|Fortran}}== |
||
{{works with|Fortran|2003 and later}} |
|||
{{incorrect|Fortran|It does not return the ''running'' standard deviation of the series.}} |
|||
<lang fortran> |
|||
{{trans|C}} |
|||
program standard_deviation |
|||
implicit none |
|||
integer(kind=4), parameter :: dp = kind(0.0d0) |
|||
real(kind=dp), dimension(:), allocatable :: vals |
|||
{{works with|Fortran|95 and later}} |
|||
integer(kind=4) :: i |
|||
real(kind=dp), dimension(8) :: sample_data = (/ 2, 4, 4, 4, 5, 5, 7, 9 /) |
|||
This one imitates C and suffers the same problems: the function is not thread-safe and must be used to compute the stddev for one set per time. |
|||
do i = lbound(sample_data, 1), ubound(sample_data, 1) |
|||
<lang fortran>program Test_Stddev |
|||
call sample_add(vals, sample_data(i)) |
|||
implicit none |
|||
write(*, fmt='(''#'',I1,1X,''value = '',F3.1,1X,''stddev ='',1X,F10.8)') & |
|||
i, sample_data(i), stddev(vals) |
|||
real, dimension(8) :: v = (/ 2,4,4,4,5,5,7,9 /) |
|||
integer :: i |
|||
real :: sd |
|||
do i = 1, size(v) |
|||
sd = stat_object(v(i)) |
|||
end do |
end do |
||
if (allocated(vals)) deallocate(vals) |
|||
print *, "std dev = ", sd |
|||
contains |
contains |
||
! Adds value :val: to array :population: dynamically resizing array |
|||
subroutine sample_add(population, val) |
|||
real(kind=dp), dimension(:), allocatable, intent (inout) :: population |
|||
real(kind=dp), intent (in) :: val |
|||
real(kind=dp), dimension(:), allocatable :: tmp |
|||
recursive function stat_object(a, cmd) result(stddev) |
|||
integer(kind=4) :: n |
|||
real, intent(in) :: a |
|||
character(len=*), intent(in), optional :: cmd |
|||
if (.not. allocated(population)) then |
|||
real, save :: summa = 0.0, summa2 = 0.0 |
|||
allocate(population(1)) |
|||
integer, save :: num = 0 |
|||
population(1) = val |
|||
real :: m |
|||
if ( .not. present(cmd) ) then |
|||
num = num + 1 |
|||
summa = summa + a |
|||
summa2 = summa2 + a*a |
|||
stddev = stat_object(0.0, "stddev") |
|||
else |
else |
||
n = size(population) |
|||
call move_alloc(population, tmp) |
|||
stddev = sqrt(stat_object(0.0, "variance")) |
|||
case("variance") |
|||
m = stat_object(0.0, "mean") |
|||
if ( num > 0 ) then |
|||
stddev = summa2/real(num) - m*m |
|||
else |
|||
stddev = 0.0 |
|||
end if |
|||
case("count") |
|||
stddev = real(num) |
|||
case("mean") |
|||
if ( num > 0 ) then |
|||
stddev = summa/real(num) |
|||
else |
|||
stddev = 0.0 |
|||
end if |
|||
case("reset") |
|||
summa = 0.0 |
|||
summa2 = 0.0 |
|||
num = 0 |
|||
case default |
|||
stddev = 0.0 |
|||
end select |
|||
end if |
|||
end function stat_object |
|||
end program Test_Stddev</lang> |
|||
===Using built-in array awareness=== |
|||
{{incorrect|Fortran|Function takes array of data at once instead of single values at a time.}} |
|||
This uses Fortran's built-in array features (which aren't available in C) |
|||
{{works with|Fortran|95 and later}} |
|||
<lang fortran> |
|||
program stats |
|||
implicit none |
|||
integer, parameter :: N = 8 |
|||
integer :: data(N) |
|||
real(8) :: mean |
|||
real(8) :: std_dev1, std_dev2 |
|||
! Set the data |
|||
data = [2,4,4,4,5,5,7,9] ! Strictly this is a Fortran 2003 construct |
|||
! Use intrinsic function 'sum' to calculate the mean |
|||
mean = sum(data)/N |
|||
! Method1: |
|||
! Calculate the standard deviation directly from the definition |
|||
std_dev1 = sqrt(sum((data - mean)**2)/N) |
|||
allocate(population(n + 1)) |
|||
! Method 2: |
|||
population(1:n) = tmp |
|||
! Use the alternative version that is less susceptible to rounding error |
|||
population(n + 1) = val |
|||
std_dev2 = sqrt(sum(data**2)/N - mean**2) |
|||
endif |
|||
end subroutine sample_add |
|||
! Calculates standard deviation for given set of values |
|||
write(*,'(a,8i2)') 'Data = ',data |
|||
real(kind=dp) function stddev(vals) |
|||
write(*,'(a,f3.1)') 'Mean = ',mean |
|||
real(kind=dp), dimension(:), intent(in) :: vals |
|||
write(*,'(a,f3.1)') 'Standard deviation (method 1) = ',std_dev1 |
|||
real(kind=dp) :: mean |
|||
write(*,'(a,f3.1)') 'Standard deviation (method 2) = ',std_dev2 |
|||
integer(kind=4) :: n |
|||
n = size(vals) |
|||
end program stats |
|||
mean = sum(vals)/n |
|||
stddev = sqrt(sum((vals - mean)**2)/n) |
|||
end function stddev |
|||
end program standard_deviation |
|||
</lang> |
</lang> |
||