12.6 External Storage or Memory-Local FFTs 535 na-mate[na]; fourew(file,&na,&nb,&nc,&nd); jk>>=1; if (jk 1)continue; mm=n; do if (nv ndim)jk=nn[++nv]; else free_vector(afc,1,KBF); free_vector(afb,1,KBF); http://www.nr. Permission is read able files free_vector(afa,1,KBF); return: 2 Jh11e(jk==1); .com or call (including this one) granted for 19881992 11-800-872 #include <stdio.h> #define SWAP(a,b)ftemp=(a);(a)=(b);(b)=ftemp n NUMERICAL RECIPES IN C: void fourew(FILE *file[5],int *na,int *nb,int *nc,int *nd) Utility used by fourfs.Rewinds and renumbers the four files. -7423(North America 电r:1s t users to make one paper by Cambridge University Press. THE int i; FILE *ftemp; 是 ART for (i=1;i<=4;i++)rewind(file[i]); Programs SWAP(file[2],file[4]); copy for thei strictly proh SwAP(f11e[1],f11e[3]); *na=3; *nb=4; Copyright (C) *nC=11 *nd=2; 1988-1992 by Numerical Recipes OF SCIENTIFIC COMPUTING(ISBN For one-dimensional data,Singleton's algorithm produces output in exactly the same order as a standard FFT(e.g.,four1).For multidimensional data,the output is the transpose of v@cam the conventional arrangement (e.g.,the output of fourn).This peculiarity,which is intrinsic to 10-:6211 the method,is generally only a minor inconvenience.For convolutions,one simply computes the component-by-component product of two transforms in their nonstandard arrangement. 43108 and then does an inverse transform on the result.Note that,if the lengths of the different dimensions are not all the same,then you must reverse the order of the values in nn [1.ndim] (thus giving the transpose dimensions)before performing the inverse transform.Note also (outside that,just like fourn,performing a transform and then an inverse results in multiplying the original data by the product of the lengths of all dimensions. North Software. We leave it as an exercise for the reader to figure out how to reorder fourfs's output into normal order,taking additional passes through the externally stored data.We doubt that Ame such reordering is ever really needed. You will likely want to modify fourfs to fit your particular application.For example, visit website machine as written,KBF =2 plays the dual role of being the size of the internal buffers,and the record size of the unformatted reads and writes.The latter role limits its size to that allowed by your machine's 1/O facility.It is a simple matter to perform multiple reads for a much larger KBF,thus reducing the number of passes by a few. Another modification of fourfs would be for the case where your virtual memory machine has sufficient address space,but not sufficient physical memory,to do an efficient FFT by the conventional algorithm (whose memory references are extremely nonlocal).In that case,you will need to replace the reads,writes,and rewinds by mappings of the arrays12.6 External Storage or Memory-Local FFTs 535 Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). } na=mate[na]; } fourew(file,&na,&nb,&nc,&nd); jk >>= 1; if (jk > 1) continue; mm=n; do { if (nv < ndim) jk=nn[++nv]; else { free_vector(afc,1,KBF); free_vector(afb,1,KBF); free_vector(afa,1,KBF); return; } } while (jk == 1); } } #include <stdio.h> #define SWAP(a,b) ftemp=(a);(a)=(b);(b)=ftemp void fourew(FILE *file[5], int *na, int *nb, int *nc, int *nd) Utility used by fourfs. Rewinds and renumbers the four files. { int i; FILE *ftemp; for (i=1;i<=4;i++) rewind(file[i]); SWAP(file[2],file[4]); SWAP(file[1],file[3]); *na=3; *nb=4; *nc=1; *nd=2; } For one-dimensional data, Singleton’s algorithm produces output in exactly the same order as a standard FFT (e.g., four1). For multidimensional data, the output is the transpose of the conventional arrangement (e.g., the output of fourn). This peculiarity, which is intrinsic to the method, is generally only a minor inconvenience. For convolutions, one simply computes the component-by-component product of two transforms in their nonstandard arrangement, and then does an inverse transform on the result. Note that, if the lengths of the different dimensions are not all the same, then you must reverse the order of the values in nn[1..ndim] (thus giving the transpose dimensions) before performing the inverse transform. Note also that, just like fourn, performing a transform and then an inverse results in multiplying the original data by the product of the lengths of all dimensions. We leave it as an exercise for the reader to figure out how to reorder fourfs’s output into normal order, taking additional passes through the externally stored data. We doubt that such reordering is ever really needed. You will likely want to modify fourfs to fit your particular application. For example, as written, KBF ≡ 2K plays the dual role of being the size of the internal buffers, and the record size of the unformatted reads and writes. The latter role limits its size to that allowed by your machine’s I/O facility. It is a simple matter to perform multiple reads for a much larger KBF, thus reducing the number of passes by a few. Another modification of fourfs would be for the case where your virtual memory machine has sufficient address space, but not sufficient physical memory, to do an efficient FFT by the conventional algorithm (whose memory references are extremely nonlocal). In that case, you will need to replace the reads, writes, and rewinds by mappings of the arrays