subroutine letrfm1(n,p,s,a,wt,iflag)
*************************************
*** 1-D Legendre transform !! MUST call leinit.f to initialize !!
***Inputs:
***        p: input, array to be transformed: DESTROYED at output !!
***     a,wt: output from leinit.f
***  iflag=1: Spectral  to Physical
***  iflag.ne.1: Physical  to Spectral  
***
***        s: output, transformed array
***
* Rmk. This program is only correct if n is even !!
*      Can be easily modified to work for all n and m !
*
*     This program should replace leptos.f and lestop.f !
*************************************
      implicit double precision (a-h,o-z)
      dimension p(0:n), s(0:n), a(0:n/2,0:n), wt(0:n)
      data alpha/1.d0/,beta/0.d0/
      
      if (iflag.ne.1) then

      do j=0,n/2-1
         s(j)=(p(j)+p(n-j))*wt(j)
         s(n/2+j+1)=(p(j)-p(n-j))*wt(j)
      enddo
      s(n/2)=p(n/2)*wt(n/2)

      call  dgemv('n',n/2+1,n/2+1,alpha,a,n/2+1,s,1,beta,p,1)

      call  dgemv('n',n/2,n/2,alpha,a(0,n/2+1),n/2+1,s(n/2+1)
     1            ,1,beta,p(n/2+1),1)


      do i=0,n/2-1
         s(2*i)=p(i)*(2*i+.5d0)
         s(2*i+1)=p(n/2+i+1)*(2*i+1.5d0)
      enddo
      s(n)=p(n/2)*n*.5d0

      else

      do i=0,n/2-1
         s(i)=p(2*i)
         s(n/2+1+i)=p(2*i+1)
      enddo
      s(n/2)=p(n)

      call  dgemv('t',n/2+1,n/2+1,alpha,a,n/2+1,s,1,beta,p,1)

      call  dgemv('t',n/2,n/2,alpha,a(0,n/2+1),n/2+1,s(n/2+1)
     1            ,1,beta,p(n/2+1),1)

      do i=0,n/2-1
         s(i)=(p(i)+p(n/2+i+1))
         s(n-i)=(p(i)-p(n/2+i+1))
      enddo
      s(n/2)=p(n/2)


      endif

      return
      end