import Numeric
import pyclimate.LanczosFilter
import math

class NewLanczosFilter(pyclimate.LanczosFilter.LanczosFilter):
  """'LinearFilter' derived class for Lanczos filtering
  """
  def getcoefs(self):
    "Define Lanczos filter coefficients"
    n=self.length/2
    if self.filtertype=='lp':
        ocoefs = self.getcoefs_lp(self.fc2)
    elif self.filtertype=='hp':
        # Define this as 1 minus a low pass filter. The normalisation of the
        # low pass filter ensures coefficients sum to zero.
        ocoefs = -self.getcoefs_lp(self.fc1)
        ocoefs[self._place(0,n,1)] += 1.0
    elif self.filtertype=='bp':
        ocoefs = self.getcoefs_lp(self.fc2) - self.getcoefs_lp(self.fc1)
    return ocoefs
    
  def getcoefs_lp(self,freq):
    "Low pass filter coefficients with proper normalization"
    n=self.length/2
    thepi=math.acos(-1.)
    ocoefs=Numeric.zeros(self.length,Numeric.Float64)
    k=0
    ocoefs[self._place(k,n,1)]=2*freq
    for ik in xrange(n):
      k=ik+1
      sigma=math.sin(thepi*k/n)*n/thepi/k
      firstfactor=math.sin(2*thepi*freq*k)/thepi/k
      ocoefs[self._place(k,n,1)]=firstfactor*sigma
      ocoefs[self._place(k,n,-1)]=firstfactor*sigma
    # Normalise
    ocoefs = ocoefs/Numeric.sum(ocoefs)
    return ocoefs


# Instances of filters (as in pyclimatetest.py)
la = pyclimate.LanczosFilter
filters={}
filters["ll"]=la.LanczosFilter("lp",0.1,0.1,10)
filters["lh"]=la.LanczosFilter("hp",0.3,0.3,10)
filters["lb"]=la.LanczosFilter("bp",0.1,0.3,10)

newfilters={}
newfilters["ll"]=NewLanczosFilter("lp",0.1,0.1,10)
newfilters["lh"]=NewLanczosFilter("hp",0.3,0.3,10)
newfilters["lb"]=NewLanczosFilter("bp",0.1,0.3,10)

for f in ("ll", "lh", "lb"):
    c = filters[f].getcoefs()
    cx = newfilters[f].getcoefs()
    print "Coefficients for", f
    for k in range(len(c)):
        print "%2d  %8.4f  %8.4f" % (k, c[k], cx[k])
    print "-------------"
    print "Sum %8.4f  %8.4f" % (Numeric.sum(c), Numeric.sum(cx))