Commit 7dfdde8e authored by thomas.forbriger's avatar thomas.forbriger

[FIX] (import_Seitosh): remove Fortran source code

The Fortran part of the vendor code is not used by DENISE. It can safely be
discarded. This way, no Fortran compiler is required when compiling and
installing the code in subdirectory contrib, except for code in aff/tests,
which however is not part of the install target.
parent 3ae84559
c this is <fcommand.f>
c ----------------------------------------------------------------------------
c ($Id$)
c
c Copyright (c) 2002 by Thomas Forbriger (IMG Frankfurt)
c
c evaluate filters specified in text form
c
c ----
c libfourier is free software; you can redistribute it and/or modify
c it under the terms of the GNU General Public License as published by
c the Free Software Foundation; either version 2 of the License, or
c (at your option) any later version.
c
c This program is distributed in the hope that it will be useful,
c but WITHOUT ANY WARRANTY; without even the implied warranty of
c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
c GNU General Public License for more details.
c
c You should have received a copy of the GNU General Public License
c along with this program; if not, write to the Free Software
c Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
c ----
c
c
c REVISIONS and CHANGES
c 11/11/2002 V1.0 Thomas Forbriger
c
cS
c ============================================================================
c
subroutine foucmd_revision
c
c print library code revision
c
cE
print *,'$Id$'
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foucmd_help(verbose)
c
logical verbose
c
c print command definitions
c
c verbose: tell library code revision
c
cE
c
print *,'filter library commands:'
print 50,'rem',' ', 'any comment (ignored)'
print 50,'# ',' ', 'any comment (ignored)'
print 50,'dif',' ', 'differentiate'
print 50,'int',' ', 'integrate'
print 50,'lp2','in,h',
& 'second order low-pass (see foufil_lp2)'
print 50,'hp2','in,h',
& 'second order high-pass (see foufil_hp2)'
print 50,'lpb','in,ord',
& 'Butterworth low-pass of order ord (see foufil_lpb)'
print 50,'hpb','in,ord',
& 'Butterworth high-pass of order ord (see foufil_hpb)'
print 50,'fac','factor', 'mutiply by factor'
print 50,'mod','normal',
& 'switch back to normal filters (default)'
print 50,'mod','inverse', 'switch to inverse filters'
print 50,'mod','period',
& 'switch to specification by period (default)'
print 50,'mod','frequency','switch to specification by frequency '
print 50,'end',' ', 'terminate file reading'
c
if (verbose) then
print *,' '
call fou_revision
call foufil_revision
call foucmd_revision
endif
c
return
50 format(1x,a3,1x,a,t18,1x,a)
end
c
cS
c----------------------------------------------------------------------
c
subroutine foucmd_command(line)
c
c evaluate filter specification given as a line of text
c
c rem any comment (ignored)
c # any comment (ignored)
c dif differentiate
c int integrate
c lp2 in,h second order low-pass (see foufil_lp2)
c hp2 in,h second order high-pass (see foufil_hp2)
c lpb in,ord Butterworth low-pass of order ord (see foufil_lpb)
c hpb in,ord Butterworth high-pass of order ord (see foufil_hpb)
c fac factor mutiply by factor
c mod normal switch back to normal filters (default)
c mod inverse switch to inverse filters
c mod period switch to specification by period (default)
c mod frequency switch to specification by frequency
c
character*(*) line
c
cE
double precision v1,v2
integer i1
c
if ((line(1:3).ne.'rem').and.
& (line(1:1).ne.'#')) then
if (line(1:4).eq.'mod ') then
if (line(5:8).eq.'norm') then
call fou_normal
elseif (line(5:8).eq.'inve') then
call fou_inverse
elseif (line(5:8).eq.'freq') then
call foufil_frequency
elseif (line(5:8).eq.'peri') then
call foufil_period
else
print *,'WARNING (foucmd_command): unknown mode'
print *,line
endif
elseif (line(1:3).eq.'dif') then
call foufil_dif
elseif (line(1:3).eq.'int') then
call foufil_int
elseif (line(1:4).eq.'fac ') then
read(line(5:), *, err=99,end=98) v1
call fou_numer(v1)
elseif (line(1:4).eq.'lp2 ') then
read(line(5:), *, err=99,end=98) v1,v2
call foufil_lp2(v1,v2)
elseif (line(1:4).eq.'hp2 ') then
read(line(5:), *, err=99,end=98) v1,v2
call foufil_hp2(v1,v2)
elseif (line(1:4).eq.'lpb ') then
read(line(5:), *, err=99,end=98) v1,i1
call foufil_lpb(v1,i1)
elseif (line(1:4).eq.'hpb ') then
read(line(5:), *, err=99,end=98) v1,i1
call foufil_hpb(v1,i1)
else
print *,'WARNING (foucmd_command): unknown command'
print *,line
endif
endif
c
return
99 stop 'ERROR (foucmd_command): reading arguments'
98 stop 'ERROR (foucmd_command): too few arguments'
end
c
cS
c----------------------------------------------------------------------
c
subroutine foucmd_unit(lu, verbose)
c
c read filter specifications from logical file unit lu
c
c end terminates reading
c
c NOTICE: this subroutine does NOT call any of the initialization routines
c
integer lu
logical verbose
c
cE
character*79 line
line=' '
do while (line(1:3).ne.'end')
read(lu, '(a79)', err=99, end=98) line
if (verbose) print *,line
if (line(1:3).ne.'end') call foucmd_command(line)
enddo
c
return
99 stop 'ERROR (foucmd_unit): reading line'
98 stop 'ERROR (foucmd_unit): unexpected end'
end
c
cS
c----------------------------------------------------------------------
c
subroutine foucmd_file(filename, verbose)
c
c read filter specifications from file
c
c NOTICE: this subroutine does NOT call any of the initialization routines
c
character*(*) filename
logical verbose
c
cE
integer lu
parameter(lu=21)
if (verbose) print *,'read filter commands from ',
& filename(1:index(filename,' ')-1)
open(lu, file=filename, err=99)
call foucmd_unit(lu, verbose)
close(lu, err=98)
c
return
99 stop 'ERROR (foucmd_file): opening file'
98 stop 'ERROR (foucmd_file): closing file'
end
c
c ----- END OF fcommand.f -----
c this is <filters.f>
c ----------------------------------------------------------------------------
c ($Id$)
c
c Copyright (c) 2002 by Thomas Forbriger (IMG Frankfurt)
c
c ready made filters
c
c ----
c libfourier is free software; you can redistribute it and/or modify
c it under the terms of the GNU General Public License as published by
c the Free Software Foundation; either version 2 of the License, or
c (at your option) any later version.
c
c This program is distributed in the hope that it will be useful,
c but WITHOUT ANY WARRANTY; without even the implied warranty of
c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
c GNU General Public License for more details.
c
c You should have received a copy of the GNU General Public License
c along with this program; if not, write to the Free Software
c Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
c ----
c
c
c REVISIONS and CHANGES
c 11/11/2002 V1.0 Thomas Forbriger
c
cS
c ============================================================================
c
subroutine foufil_revision
c
c print library code revision
c
cE
print *,'$Id$'
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_clear
c
c initialize fourier filters
c set to period mode to be compatible with seife
c
include 'filters.inc'
c
cE
c
call fou_clear
call foufil_period
c
return
end
c
cS
c----------------------------------------------------------------------
c
double precision function foufil_omega(in)
c
c calculate angular frequency
c in: will be treated as a frequency given in Hz, if fourier_frequency is true
c otherwise it will be treated as a period given in sec
c
include 'filters.inc'
c
double precision in
c
cE
c
double precision pi2, value
parameter(pi2=2.d0*3.14159265358979d0)
c
if (fourier_frequency) then
value=pi2*in
else
value=pi2/in
endif
c
foufil_omega=value
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_frequency
c
c set frequency mode
c
include 'filters.inc'
c
cE
fourier_frequency=.true.
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_period
c
c set period mode
c
include 'filters.inc'
c
cE
fourier_frequency=.false.
c
return
end
c
c======================================================================
c define filters
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_int
c
c integrate once
c
cE
double complex value
value=(0.d0,0.d0)
call fou_pole(value)
value=(0.d0,1.d0)
call fou_cdenom(value)
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_dif
c
c differentiate once
c
cE
double complex value
value=(0.d0,0.d0)
call fou_zero(value)
value=(0.d0,1.d0)
call fou_cnumer(value)
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_hpb(in, ord)
c
c Butterworth high pass filter
c in: frequency or period (cf. foufil_omega)
c ord: order of filter
c
double precision in
integer ord
c
cE
double precision foufil_omega, omega, pi
integer i
double complex ime, pole, zero
parameter (ime=(0.d0,1.d0), pi=3.1415926535897931159d0)
c
if (ord.lt.1) stop 'ERROR (foufil_hpb): illegal order'
omega=foufil_omega(in)
c
zero=(0.d0,0.d0)
do i=1,ord
call fou_zero(zero)
pole=omega*exp(ime*pi*(2.d0*i-1.d0)/(2.d0*ord))
call fou_pole(pole)
enddo
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_lpb(in, ord)
c
c Butterworth low pass filter
c in: frequency or period (cf. foufil_omega)
c ord: order of filter
c
double precision in
integer ord
c
cE
double precision foufil_omega, omega, pi
integer i
double complex ime, pole, cfac
parameter (ime=(0.d0,1.d0), pi=3.1415926535897931159d0)
c
if (ord.lt.1) stop 'ERROR (foufil_hpb): illegal order'
omega=foufil_omega(in)
c
cfac=-ime*omega
do i=1,ord
call fou_cnumer(cfac)
pole=omega*exp(ime*pi*(2.d0*i-1.d0)/(2.d0*ord))
call fou_pole(pole)
enddo
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_lp2(in, h)
c
c second order low pass filter
c in: frequency or period (cf. foufil_omega)
c h: damping as a fraction of critical
c
double precision in, h
c
cE
double precision foufil_omega, omega
double complex ime, pole, cfac
parameter (ime=(0.d0,1.d0))
c
omega=foufil_omega(in)
c
cfac=-ime*omega
call fou_cnumer(cfac)
call fou_cnumer(cfac)
pole=omega*(ime*h+sqrt(1.d0-h*h))
call fou_pole(pole)
pole=omega*(ime*h-sqrt(1.d0-h*h))
call fou_pole(pole)
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine foufil_hp2(in, h)
c
c second order high pass filter
c in: frequency or period (cf. foufil_omega)
c h: damping as a fraction of critical
c
double precision in, h
c
cE
double precision foufil_omega, omega
double complex ime, pole, zero
parameter (ime=(0.d0,1.d0))
c
omega=foufil_omega(in)
c
zero=(0.d0,0.d0)
call fou_zero(zero)
call fou_zero(zero)
pole=omega*(ime*h+sqrt(1.d0-h*h))
call fou_pole(pole)
pole=omega*(ime*h-sqrt(1.d0-h*h))
call fou_pole(pole)
c
return
end
c
c
c ----- END OF filters.f -----
c this is <polesnzeroes.f>
c ----------------------------------------------------------------------------
c ($Id$)
c
c Copyright (c) 2002 by Thomas Forbriger (IMG Frankfurt)
c
c pole and zero access subroutines
c
c ----
c libfourier is free software; you can redistribute it and/or modify
c it under the terms of the GNU General Public License as published by
c the Free Software Foundation; either version 2 of the License, or
c (at your option) any later version.
c
c This program is distributed in the hope that it will be useful,
c but WITHOUT ANY WARRANTY; without even the implied warranty of
c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
c GNU General Public License for more details.
c
c You should have received a copy of the GNU General Public License
c along with this program; if not, write to the Free Software
c Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
c ----
c
c
c REVISIONS and CHANGES
c 11/11/2002 V1.0 Thomas Forbriger
c 13/11/2002 V1.1 introduce waterlevel to omega
c
cS
c ============================================================================
c
subroutine fou_revision
c
c print library code revision
c
cE
print *,
& '$Id$'
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine fou_clear
c
c initialize poles and zeros database
c
include 'polesnzeros.inc'
c
cE
call fou_normal
fourier_numer=(1.d0,0.d0)
fourier_denom=(1.d0,0.d0)
fourier_npoles=0
fourier_nzeros=0
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine fou_inverse
c
c set inverse flag
c all following coefficients will be treated as inverse filters
c
include 'polesnzeros.inc'
c
cE
c
fourier_normal=.false.
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine fou_normal
c
c reset inverse flag
c all following coefficients will be treated as normal filters
c
include 'polesnzeros.inc'
c
cE
c
fourier_normal=.true.
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine fou_pole(pole)
c
c set one pole
c
double complex pole
c
include 'polesnzeros.inc'
c
cE
c
if (fourier_normal) then
if (fourier_npoles.lt.fourier_nmax) then
fourier_npoles=fourier_npoles+1
fourier_poles(fourier_npoles)=pole
else
stop 'ERROR (fourier_pole): too many poles'
endif
else
if (fourier_nzeros.lt.fourier_nmax) then
fourier_nzeros=fourier_nzeros+1
fourier_zeros(fourier_nzeros)=pole
else
stop 'ERROR (fourier_pole): too many zeros'
endif
endif
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine fou_zero(zero)
c
c set one zero
c
double complex zero
c
include 'polesnzeros.inc'
c
cE
c
if (fourier_normal) then
if (fourier_nzeros.lt.fourier_nmax) then
fourier_nzeros=fourier_nzeros+1
fourier_zeros(fourier_nzeros)=zero
else
stop 'ERROR (fourier_zero): too many zeros'
endif
else
if (fourier_npoles.lt.fourier_nmax) then
fourier_npoles=fourier_npoles+1
fourier_poles(fourier_npoles)=zero
else
stop 'ERROR (fourier_zero): too many poles'
endif
endif
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine fou_numer(factor)
c
c set one numerator factor
c
double precision factor
c
include 'polesnzeros.inc'
c
cE
c
if (fourier_normal) then
fourier_numer=factor*fourier_numer
else
fourier_denom=factor*fourier_denom
endif
c
return
end
c
cS
c----------------------------------------------------------------------
c
subroutine fou_denom(factor)
c
c set one denominator factor
c
double precision factor
c
include 'polesnzeros.inc'
c
cE
c
if (fourier_normal) then
fourier_denom=factor*fourier_denom
else
fourier_numer=factor*fourier_numer
endif
c
return
end
c
cS
c----------------------------------------------------------------------