1 c===========================================================================
3 c This file is part of TISEAN
5 c Copyright (c) 1998-2007 Rainer Hegger, Holger Kantz, Thomas Schreiber
7 c TISEAN is free software; you can redistribute it and/or modify
8 c it under the terms of the GNU General Public License as published by
9 c the Free Software Foundation; either version 2 of the License, or
10 c (at your option) any later version.
12 c TISEAN is distributed in the hope that it will be useful,
13 c but WITHOUT ANY WARRANTY; without even the implied warranty of
14 c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 c GNU General Public License for more details.
17 c You should have received a copy of the GNU General Public License
18 c along with TISEAN; if not, write to the Free Software
19 c Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
21 c===========================================================================
22 c space time separation plot
23 c see H. Kantz, T. Schreiber, Nonlinear Time Series Analysis, Cambridge
24 c University Press (1997,2004)
25 c author T. Schreiber (1998) based on earlier version
26 c==========================================================================
27 parameter(nx=1000000,mdt=500,mfrac=100)
28 dimension x(nx), stp(mfrac,mdt)
29 character*72 file, fout
30 data idt/1/, perc/0.05/, ndt/100/
33 call whatido("space-time separation plot",iverb)
37 ndt=min(ican("t",ndt),mdt)
39 nfrac=min(mfrac,int(1/perc))
44 isout=igetout(fout,iverb)
45 if(iv_io(iverb).eq.1) write(istderr(),*) "computing ", nfrac,
46 . " levels at fractions ", perc, 2*perc, "..."
48 call nthstring(1,file)
49 call readfile(nmax,x,nexcl,jcol,file,iverb)
50 call minmax(nmax,x,xmin,xmax)
51 call stplot(nmax,x,id,m,xmax-xmin,stp,nfrac,ndt,idt)
52 if(isout.eq.1) call addsuff(fout,file,"_stp")
53 call outfile(fout,iunit,iverb)
56 20 write(iunit,*) it*idt, stp(iper,it)
64 . " -d# -m# [-## -t# -%# -o outfile -l# -x# -c# -V# -h] file")
65 call popt("d","delay")
66 call popt("m","embedding dimension")
67 call popt("#","time resolution (1)")
68 call popt("t","time steps (100, <500)")
69 call popt("%","fraction at wich to create levels (0.05, >0.01)")
70 call popt("l","number of values to be read (all)")
71 call popt("x","number of values to be skipped (0)")
72 call popt("c","column to be read (1 or file,#)")
78 subroutine stplot(nmax,y,id,m,epsmax,stp,nfrac,mdt,idt)
79 parameter(meps=1000,mfrac=100)
80 dimension y(nmax),stp(mfrac,mdt),ihist(meps)
85 do 30 n=it*idt+(m-1)*id+1,nmax
86 dis=0 ! compute distance in m dimensions
88 40 dis=max(dis,abs(y(n-me*id)-y(n-me*id-it*idt)))
89 ih=min(int(meps*dis/epsmax)+1,meps)
90 30 ihist(ih)=ihist(ih)+1
92 need=(nmax-it*idt-(m-1)*id)*ifrac/real(nfrac)
96 50 if(is.ge.need) goto 1
97 1 stp(ifrac,it)=ieps*epsmax/meps