## powertf.R

powertf=function(T,coeff1, coeff0, coeff2){

## Input
## T = resolution in frequency (length of time series)
## coeff1 = vector of filter coefficients for PAST data
## coeff2 = FUTURE data 
## coeff0 == present
## Assume that length of coeff1 is the same as length of 
## coeff2

## Output 
## powtfout = spectrum of the ARMA process defined from (0, pi]

freq = 2*pi*seq(0, T-1)/T;
L = length(coeff1); 

powtf=c(1:T); 

for (k in 1:(T)){
theta.temp=coeff0;

for (m in 1:L)
{
theta.temp = theta.temp + coeff1[m]*complex(mod=1, arg=-1*m*freq[k]) + 
coeff2[m]*complex(mod=1, arg=1*m*freq[k])}

powtf[k] = (abs(theta.temp))^2
};

powtfout = powtf[2:(T/2+1)];
powtfout}