another fix for rectangular transforms

Author: MikaelSlevinsky

src/SphericalHarmonics/fastplan.jl

@@ -11,7 +11,7 @@ end
 function FastSphericalHarmonicPlan(A::Matrix{T}, L::Int; opts...) where T
     M, N = size(A)
     n = (N+1)÷2
-    RP = RotationPlan(T, n-1)
+    RP = RotationPlan(T, M-1)
     p1 = plan_normleg2cheb(A)
     p2 = plan_normleg12cheb2(A)
     p1inv = plan_cheb2normleg(A)

src/SphericalHarmonics/slowplan.jl

@@ -106,8 +106,8 @@ function Base.A_mul_B!(P::RotationPlan, A::AbstractMatrix)
     @stepthreads for m = M÷2:-1:2
         @inbounds for j = m:-2:2
             for l = N-j:-1:1
-                s = snm[l+(j-2)*(M+2-j)÷2]
-                c = cnm[l+(j-2)*(M+2-j)÷2]
+                s = snm[l+(j-2)*(2*N+3-j)÷2]
+                c = cnm[l+(j-2)*(2*N+3-j)÷2]
                 a1 = A[l+N*(2*m-1)]
                 a2 = A[l+2+N*(2*m-1)]
                 a3 = A[l+N*(2*m)]
@@ -129,8 +129,8 @@ function Base.At_mul_B!(P::RotationPlan, A::AbstractMatrix)
     @stepthreads for m = M÷2:-1:2
         @inbounds for j = reverse(m:-2:2)
             for l = 1:N-j
-                s = snm[l+(j-2)*(M+2-j)÷2]
-                c = cnm[l+(j-2)*(M+2-j)÷2]
+                s = snm[l+(j-2)*(2*N+3-j)÷2]
+                c = cnm[l+(j-2)*(2*N+3-j)÷2]
                 a1 = A[l+N*(2*m-1)]
                 a2 = A[l+2+N*(2*m-1)]
                 a3 = A[l+N*(2*m)]
@@ -199,8 +199,7 @@ end
 
 function SlowSphericalHarmonicPlan(A::Matrix{T}) where T
     M, N = size(A)
-    n = (N+1)÷2
-    RP = RotationPlan(T, n-1)
+    RP = RotationPlan(T, M-1)
     p1 = plan_normleg2cheb(A)
     p2 = plan_normleg12cheb2(A)
     p1inv = plan_cheb2normleg(A)

src/SphericalHarmonics/thinplan.jl

@@ -15,7 +15,7 @@ end
 function ThinSphericalHarmonicPlan(A::Matrix{T}, L::Int; opts...) where T
     M, N = size(A)
     n = (N+1)÷2
-    RP = RotationPlan(T, n-1)
+    RP = RotationPlan(T, M-1)
     p1 = plan_normleg2cheb(A)
     p2 = plan_normleg12cheb2(A)
     p1inv = plan_cheb2normleg(A)