aovsub Module

Utility routines for Aov.
Note that only public routines are reliable and tested. Other ones are kept here for past/future development.


Uses


Contents


Interfaces

public interface sortx

  • private subroutine sortx_m_r4(dat, ina)

    (C) Alex Schwarzenberg-Czerny, 2005 alex@camk.edu.pl Stable index sorting by merge sort. Requires extra memory (inb)

    Arguments

    Type IntentOptional AttributesName
    real(kind=4), intent(in) :: dat(:)

    unsorted data

    integer, intent(out) :: ina(:)

    sort index

  • private subroutine sortx_m_r8(dat, ina)

    (C) Alex Schwarzenberg-Czerny, 2005 alex@camk.edu.pl Stable index sorting by merge sort. Requires extra memory (inb)

    Arguments

    Type IntentOptional AttributesName
    real(kind=8), intent(in) :: dat(:)

    unsorted data

    integer, intent(out) :: ina(:)

    sort index


Functions

public function dict_val(dict, isize, key)

Find key position in a keyword dictionary. Restriction: dictionary entries are of equal size, all kept in one string.

Arguments

Type IntentOptional AttributesName
character(len=*), intent(in) :: dict

Dictionary string

integer, intent(in) :: isize

size of keys

character(len=*), intent(in) :: key

key string

Return Value integer

Key position

public function tolower(str)

Change string case to lower

Arguments

Type IntentOptional AttributesName
character(len=*), intent(in) :: str

input string

Return Value character(len=Len(str))

result output string

public function toupper(str)

Change string case to upper

Arguments

Type IntentOptional AttributesName
character(len=*), intent(in) :: str

input string

Return Value character(len=Len(str))

result output string

public function givens(r, idf, detnth)

overdetermined linear system by accurate Givens rotations. Prior to GIVENS call m-times GIVENSA in order to build the triangle system R (C) Alex Schwarzenberg-Czerny,2000 alex@camk.edu.pl

Arguments

Type IntentOptional AttributesName
real(kind=SP), intent(inout) :: r(:,:)

r(n1,:) - unknowns (r(k,k),k=1,n) variances (r(:k,k),k=1,n) half of covariance matrix (r(k+1:,k),k=2,n) half of inverse triangle root Q r(n1,n1) - chi^2/idf (=RHS variance=std.dev^2) For idf<0: ((r(j,k),j=1,k),k=1,n) inverse triangle root Q

integer, intent(in) :: idf

idf - >0 solve normal equations for idf degrees of freedom <0 and r(n1,n1)=1: inverse triangle matrix R

real(kind=SP), intent(out) :: detnth

=Det(R)^(1/n) - n-th root of determinant

Return Value real(kind=SP)

result condition ratio (min eigenv/max eigenv), 0 signalls singular system


Subroutines

public subroutine givensacc(a, r)

overdetermined linear system by Givens rotations, a state-of-art algorithm. Routine takes current observation equation from a and accumulates it into r: (C) Alex Schwarzenberg-Czerny,2000 alex@camk.edu.pl

Arguments

Type IntentOptional AttributesName
real(kind=SP), intent(in) :: a(:)

a(n1) - an equation with its RHS as a(n1) element, multiplied by sqrt of weight. Note: sum of all weights must equall no.observations, lest GIVENS returns wrong covariance estimates.

real(kind=SP), intent(inout) :: r(:,:)