-- aₙ, bₙ, cₙ, dₙ in Kuhl+Giardina
λxs.λys.λN.
  { sum ← [(+)/x]
  ; tieSelf ← [({.x)⊳x]; Δ ← ((-)\~) ∴ tieSelf
  ; dxs ⟜ Δ xs; dys ⟜ Δ ys
  ; dts ⟜ [√(x^2+y^2)]`dxs dys
  ; dxss ⟜ ((%)`dxs dts); dyss ⟜ ((%)`dys dts)
  ; pts ⟜ (+)Λₒ 0 dts; T ⟜}. pts
  ; coeffs ← λn.
    { n ⟜ ℝn; k ⟜ 2*n*𝜋%T
    ; cosDiffs ⟜ (-)\~((cos. ∴ (k*))'pts)
    ; sinDiffs ⟜ (-)\~((sin. ∴ (k*))'pts)
    ; c ⟜ T%(2*n^2*𝜋^2)
    ; aₙ ← c*sum ((*)`dxss cosDiffs)
    ; bₙ ← c*sum ((*)`dxss sinDiffs)
    ; cₙ ← c*sum ((*)`dyss cosDiffs)
    ; dₙ ← c*sum ((*)`dyss sinDiffs)
    ; (aₙ,bₙ,cₙ,dₙ)
    }
  ; coeffs'(irange 0 N 1)
  }