[R] Re: Rwave cgt plot time axis problem

zhu wang zhuw at mail.smu.edu
Fri Aug 15 16:48:41 CEST 2003

Thanks to Brandon Whitcher. I also would like to thank Andy Liaw for his
help through personal communication.

--
Zhu Wang

Statistical Science Department
Southern Methodist University
Phone: (214)768-2453
Fax: (214)768-4035
Email: zhuw at mail.smu.edu

On Fri, 2003-08-15 at 05:54, Brandon Whitcher wrote:
> On Mon, 4 Aug 2003, zhu wang wrote:
>
> > When I use the function cgt of library Rwave, the time axis of the plot
> > is always on the [0,1] scale regardless of the original time. This
> > happedn when I tried to reproduce some pictures of "Practical
> > Time-Frequency Analysis", e.g. Figure 3.5 and 3.6 using the codes
> > provided in the book. But the figures of the book display the real
> > sampling time using the same codes. It seems there is no argument for
> > the time scale for 'cgt'. Is this a bug or my misunderstanding?
>
> This is the result of how the function was originally written and how
> image() deals with its default x and y axis labels.  In R, they default to
> [0,1] whereas in S-plus (if I recall correclty) they default to 1:nrow(x)
> and 1:ncol(x) for a matrix x.
>
> I have been trying to produce text files for each chapter (using the
> examples from the book) with limited success.  That is, some of the
> example fail to work.  Luckily, for chapter 3 this doesn't appear to be
> the case.  The example is
>
> ## Example 3.4 Chirp
> ## Generate the chirp and compute the Gabor transform between frequencies
> ## 0 and 0.125 Hz:
> x <- 1:512
> chirp <- sin(2*pi*(x + 0.002*(x-256)*(x-256))/16)
> par(mfrow=c(3,1))
> plot(ts(chirp), xaxs="i", xlab="", ylab="")
> title("Chirp signal")
> cgtchirp <- cgt(chirp, 50, .005, 25)
> tmp <- cleanph(cgtchirp, .01)
> ## The result is displayed in Figure 3.5
>
> This works (at least on my SGI running R 1.6.2) but produces axis labels,
> for the two images, of [0,1]x[0,1].
>
> An easy work around, until I implement it myself in a new version of
> Rwave, is to rewrite the function (call it mycgt) using
>
> if(plot){
>   image(1:newsize, seq(0, nvoice*freqstep/2, length=nvoice),
>         Mod(output), xlab = "Time", ylab = "Frequency")
>   title("Gabor Transform Modulus")
> }
>
> ...for the plot statement near the bottom.  It appears to work.  A similar
> hack may be performed for cleanph().
>
> Brandon
>
> --------------------------------------------------------------------------
>  Senior Researcher in Imaging                             GlaxoSmithKline
>  Research Statistics Unit              New Frontiers Science Park (South)
>                                           Third Avenue, Harlow   CM19 5AW
>  phone:  +44 (0)127 963 1285                               United Kingdom
>  fax:  +44 (0)127 964 4004
> --------------------------------------------------------------------------