R devel - Apr 1999 - mvfft (PR#176)

This message is in MIME format
--_=XFMail.1.3.p0.Linux:990422185944:4782=_
Content-Type: text/plain; charset=us-ascii

R-0.64.0 on RedHat Linux 5.1

Some problems with multivariate fast fourier transform. I
have attached the dump of a 30 x 3 matrix that seems to
reliably reproduce these problems

1) mvfft doesn't like vectors. It complains about not having
   enough memory:

R> mvfft(xx[,1])
Error: heap memory (2048 Kb) exhausted [needed 1024 Kb more]
       See "help(Memory)" on how to increase the heap size.

What in fact is happening here is that do_mvfft tries to
get the dimensions of the argument (fourier.c: 125), but 
gets junk instead, since the argument is not an array.

2) mvfft causes a segmentation fault or gets stuck in a loop

R> mvfft(xx)

Program received signal SIGSEGV, Segmentation fault.
0x80d36cb in fftmx (a=0x40617eb0, b=0x40617eb8, ntot=30, n=30, nspan=30, 
    isn=-2, m=3, kt=0, at=0x407c8f58, ck=0x407c8f80, bt=0x407c8fa8, 
    sk=0x407c8fd0, np=0x407c8ee4, nfac=0x8115a78) at fft.c:558
558         bt[j] = b[k1] - b[k2];
(gdb) print bt[j]
Cannot access memory at address 0x407cb000.
(gdb) print kspan
$4 = 0

We're in a loop here where kspan is the increment for k1 and decrement
for k2. The loop breaks when k1 < k2, but somehow kspan is zero.  How it
got there, I can't say, as we're in goto hell at this point.

Martyn

--_=XFMail.1.3.p0.Linux:990422185944:4782=_
Content-Disposition: attachment; filename="dumpdata"
Content-Transfer-Encoding: 7bit
Content-Description: dumpdata
Content-Type: text/plain; charset=us-ascii; name=dumpdata; SizeOnDisk=946

"xx" <-
structure(c(2.80791, 2.92504, 3.02911, 2.85911, 3.28071, 3.43015, 
2.3898, 1.67055, 3.81164, 2.43832, 2.92065, 1.37547, 3.23811, 
2.96203, 2.83998, 2.36596, 3.51832, 2.59301, 3.08633, 2.87675, 
3.16977, 2.81822, 2.48206, 3.28322, 2.51132, 3.01153, 3.25603, 
2.98599, 2.93238, 2.12555, 1.15452, 0.625684, 0.708469, 0.984351, 
0.576846, 0.228625, 1.40307, 1.21012, 0.143615, 1.4046, 1.32799, 
0.0254974, 1.17005, 1.28567, 1.51562, 0.842245, 0.852613, 0.780749, 
0.967637, 0.810536, 0.858167, 0.960625, 0.613509, 1.38389, 1.22199, 
0.93934, 0.479022, 0.684411, 0.907639, 0.694256, 0.637226, 0.698318, 
0.404093, 0.979977, 1.42097, 0.761441, 1.0802, 2.09035, 2.85473, 
2.44305, 1.81988, 1.78082, 1.69385, 0.643862, 1.5311, 0.960032, 
0.69531, 0.821253, 0.49187, 0.394229, 0.806075, 0.557793, 0.694918, 
1.46079, 1.56847, 0.823876, 0.73726, 0.698766, 0.800476, 1.45114
), .Dim = c(30, 3), .Dimnames = list(NULL, c("alpha",
"beta",
"sigma")))

--_=XFMail.1.3.p0.Linux:990422185944:4782=_--
End of MIME message

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To:
r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._

This message is in MIME format
--_=XFMail.1.3.p0.Linux:990427135205:436=_
Content-Type: text/plain; charset=us-ascii

This should fix it.
1) Check that the dimensions are not null.
2) Call fft_factor before each call to fft_work.  The array nfac
   is set by fft_factor and overwritten by fftmx (in src/appl/fft.c,
   called by fft_work). Since nfac is not an argument of fft_work,
   nor of fft_factor it is not easy to spot the bug!

--_=XFMail.1.3.p0.Linux:990427135205:436=_
Content-Disposition: attachment; filename="fourier.diff"
Content-Transfer-Encoding: 7bit
Content-Description: fourier.diff
Content-Type: text/plain; charset=us-ascii; name=fourier.diff; SizeOnDisk=833

diff -uNr R-0.64.0/src/main/fourier.c R-fixed/src/main/fourier.c
--- R-0.64.0/src/main/fourier.c	Fri Feb  5 07:44:44 1999
+++ R-fixed/src/main/fourier.c	Tue Apr 27 12:44:48 1999
@@ -123,7 +123,7 @@
     z = CAR(args);
 
     d = getAttrib(z, R_DimSymbol);
-    if (length(d) > 2)
+    if (d == R_NilValue || length(d) > 2)
 	errorcall(call, "vector-valued series required\n");
     n = INTEGER(d)[0];
     p = INTEGER(d)[1];
@@ -156,9 +156,11 @@
 	    errorcall(call, "fft factorization error\n");
 	work = (double*)R_alloc(4 * maxf, sizeof(double));
 	iwork = (int*)R_alloc(maxp, sizeof(int));
-	for (i = 0; i < p; i++)
+	for (i = 0; i < p; i++) {
+	    fft_factor(n, &maxf, &maxp);
 	    fft_work(&(COMPLEX(z)[i*n].r), &(COMPLEX(z)[i*n].i),
 		     1, n, 1, inv, work, iwork);
+        }
 	vmaxset(vmax);
     }
     UNPROTECT(1);

--_=XFMail.1.3.p0.Linux:990427135205:436=_--
End of MIME message
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To:
r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._