A double precision version PDA_DRFFTB of the routine has been added.
******************************************************************
subroutine pda_rfftb(n,r,wsave)
******************************************************************
subroutine pda_rfftb computes the real perodic sequence from its
fourier coefficients (fourier synthesis). the transform is defined
below at output parameter r.
input parameters
n the length of the array r to be transformed. the method
is most efficient when n is a product of small primes.
n may change so long as different work arrays are provided
r a real array of length n which contains the sequence
to be transformed
wsave a work array which must be dimensioned at least 2*n+15.
in the program that calls pda_rfftb. the wsave array must be
initialized by calling subroutine pda_rffti(n,wsave) and a
different wsave array must be used for each different
value of n. this initialization does not have to be
repeated so long as n remains unchanged thus subsequent
transforms can be obtained faster than the first.
the same wsave array can be used by pda_rfftf and pda_rfftb.
output parameters
r for n even and for i = 1,...,n
r(i) = r(1)+(-1)**(i-1)*r(n)
plus the sum from k=2 to k=n/2 of
2.*r(2*k-2)*cos((k-1)*(i-1)*2*pi/n)
-2.*r(2*k-1)*sin((k-1)*(i-1)*2*pi/n)
for n odd and for i = 1,...,n
r(i) = r(1) plus the sum from k=2 to k=(n+1)/2 of
2.*r(2*k-2)*cos((k-1)*(i-1)*2*pi/n)
-2.*r(2*k-1)*sin((k-1)*(i-1)*2*pi/n)
***** note
this transform is unnormalized since a call of pda_rfftf
followed by a call of pda_rfftb will multiply the input
sequence by n.
wsave contains results which must not be destroyed between
calls of pda_rfftb or pda_rfftf.