Commit 5f771d7e authored by thomas.forbriger's avatar thomas.forbriger Committed by thomas.forbriger
Browse files

not much better

This is a legacy commit from before 2015-03-01.
It may be incomplete as well as inconsistent.
See COPYING.legacy and README.history for details.


SVN Path:     http://gpitrsvn.gpi.uni-karlsruhe.de/repos/TFSoftware/trunk
SVN Revision: 513
SVN UUID:     67feda4a-a26e-11df-9d6e-31afc202ad0c
parent 5204261a
c this is <siggen.f>
c------------------------------------------------------------------------------
c ($Id: siggen.f,v 1.1 2001-10-24 16:03:24 forbrig Exp $)
c ($Id: siggen.f,v 1.2 2001-10-24 21:32:49 forbrig Exp $)
c
c 24/10/2001 by Thomas Forbriger (IMGF Frankfurt)
c
......@@ -16,16 +16,17 @@ c
character*(*) version
parameter(version='SIGGEN V1.0 SIGnal GENerator')
character*(*) SIGGEN_CVS_ID
parameter(SIGGEN_CVS_ID='$Id: siggen.f,v 1.1 2001-10-24 16:03:24 forbrig Exp $')
parameter(SIGGEN_CVS_ID='$Id: siggen.f,v 1.2 2001-10-24 21:32:49 forbrig Exp $')
c
c parameters
integer nsig,ncyc
real f,t,a,d,ta,te,td,tm,f1,f2
double precision f,t,a,d,ta,te,td,tm,f1,f2,expo
logical overwrite
character*80 filename
c internal parameters
integer nsamples,i
double precision ti,freq,tend,b,fh1,fh2
double precision myexpo
real pi
parameter(pi= 3.1415926535897)
integer lu
......@@ -37,7 +38,7 @@ c data space
c commandline
integer maxopt, lastarg, iargc
character*80 argument
parameter(maxopt=14)
parameter(maxopt=15)
character*3 optid(maxopt)
character*40 optarg(maxopt)
logical optset(maxopt), opthasarg(maxopt)
......@@ -45,9 +46,10 @@ c debugging
logical debug, verbose
c here are the keys to our commandline options
data optid/2h-D, 2h-v, 2h-o, 2h-f, 2h-T, 2h-a, 2h-d,
& 3h-Ta,3h-Te,2h-n,3h-Td,3h-Tm,3h-f1,3h-f2/
data opthasarg/3*.FALSE.,11*.TRUE./
data optarg/3*1h-,3h20.,3*2h1.,2h0.,3h20.,1h5,5h1.e20,3h30.,2*3h20./
& 3h-Ta,3h-Te,2h-n,3h-Td,3h-Tm,3h-f1,3h-f2,2h-e/
data opthasarg/3*.FALSE.,12*.TRUE./
data optarg/3*1h-,3h20.,3*2h1.,2h0.,3h20.,1h5,5h1.e20,3h30.,
& 2*3h20.,1h1/
c
c------------------------------------------------------------------------------
c basic information
......@@ -94,6 +96,8 @@ c
print *,' (default: ',optarg(12)(1:3),')'
print *,'-n n set cycle parameter n to ''n'' cycles'
print *,' (default: ',optarg(10)(1:3),')'
print *,'-e e set exponent parameter e to ''e'' '
print *,' (default: ',optarg(15)(1:3),')'
print *,'-f1 f1 set frequency parameter f1 to ''f1''Hz'
print *,' (default: ',optarg(13)(1:3),')'
print *,'-f2 f2 set frequency parameter f2 to ''f2''Hz'
......@@ -166,6 +170,7 @@ c
read(optarg(12), *) tm
read(optarg(13), *) f1
read(optarg(14), *) f2
read(optarg(15), *) expo
c
c------------------------------------------------------------------------------
......@@ -232,34 +237,35 @@ c
print 50,'onset of wavelet','Ta',ta,'ms'
print 50,'damping time constant','Td',td,'ms'
print 51,'number of cycles','n',ncyc,' '
print 50,'frequency variation exponent','e',expo,' '
print 50,'initial frequency','f1',f1,'Hz'
print 50,'final frequency','f2',f2,'Hz'
endif
nsamples=int(1.e3*t/d)
c b=((2*f2-f1)**2-f1**2)/(4*ncyc)
if (f1.gt.f2) then
b=-b
endif
tend=4.*float(ncyc)*(f2-f1)/((2.*f2-f1)**2-f1**2)
b=(f2-f1)/tend
tend=ta+1.e3*tend
c tend=1.e3*(1.e-3*ta-f1/(2.*b)+sqrt(f1**2+4.*b*float(ncyc))/(2.*b))
tend=float(ncyc*(expo+1))/(f2+float(expo)*f1)
b=((f2-f1)/float(expo+1))*myexpo(tend,-expo)
tend=ta+1.e3*tend
if (verbose) then
print *,'derived parameters:'
print 52,'number of samples',nsamples,' '
print 53,'end of signal',tend,'ms'
print 53,'frequency modulation slope',b,'Hz/s'
print 53,'initial frequency',f1,'Hz/s'
print 53,'final frequency',f1+b*1.e-3*(tend-ta),'Hz/s'
print 53,'initial frequency',f1,'Hz'
print 53,'final frequency',
& f1+(expo+1)*b*myexpo(dble(1.e-3*(tend-ta)),expo),'Hz'
print 53,'total cycles',
& (f1+b*myexpo(dble(1.e-3*(tend-ta)),expo))*1.e-3*(tend-ta),' '
endif
do i=1,nsamples
ti=d*1.e-3*(i-1)
ti=d*1.e-3*float(i-1)
if (ti.lt.(1.e-3*ta)) then
data(i)=0.
elseif(ti.le.(1.e-3*tend)) then
freq=f1+b*ti
data(i)=a*(sin(2*pi*freq*(ti-1.e-3*Ta))*
data(i)=a*(sin(2*pi*(f1+b*myexpo(dble(ti-1.e-3*Ta),expo))*(ti-1.e-3*Ta))*
& exp(-1.e3*(ti-1.e-3*Ta)/Td))
data(i)=f1+(expo+1.)*b*myexpo(dble(ti-1.e-3*Ta),expo)
else
data(i)=0.
endif
......@@ -281,5 +287,19 @@ c tend=1.e3*(1.e-3*ta-f1/(2.*b)+sqrt(f1**2+4.*b*float(ncyc))/(2.*b))
52 format(3x,a30,1x,3x,1x,i10,a4)
53 format(3x,a30,1x,3x,1x,f10.3,a4)
end
c
c----------------------------------------------------------------------
double precision function myexpo(x,y)
double precision x,y,result
if (x.lt.0.d0) stop 'ERROR: argument of myexpo'
if (x.lt.1.d-2) then
result=0.d0
else
result=exp(y*log(x))
endif
myexpo=result
return
end
c
c ----- END OF siggen.f -----
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment