feat: implemented encoding up to 7th order

Alexander Nguyen

18 hours ago d6bb375789821f4e4e175465f50e77e8a6b1aaf7

 sources/classes/sox-hoaencode.lisp

@@ -1,11 +1,13 @@
;Authors: A. Nguyen, 2025.
;Authors: 
; A. Nguyen, 2025
; M. Schumacher, 2025

; Design limitations:
; 1) Sounds must have identical sample rate; otherwise, SoX fails silently, literally.

; Questions:
; 1) Is there a built-in (flatten lst)-function?
; 2) Implement up to which order? <=> CLI issues to be expected (e.g., maximum command character length)?
; 1) Is there a built-in (flatten lst)-function? Yes, om:flat
; 2) Implement up to which order? <=> CLI issues to be expected (e.g., maximum command character length)? max character length for CLI should be changeable 

(in-package :om)

@@ -34,29 +36,74 @@
    )
)

; Ambisonics
;  ### Ambisonics ###

(defun sox-hoaencode-sn3d-factor (order degree)
"This is a placeholder for documentation"
    (ecase order
        (0 1)
        (1 1)
        (2 
            (ecase (abs degree)
                (2 (/ (sqrt 3) 6))
                (1 (/ (sqrt 3) 3))
                (0 1)
                (1 (/ (sqrt 3) 3))
                (2 (/ (sqrt 3) 6))
            ))
        (3 
            (ecase (abs degree)
                (3 (/ (sqrt 10) 60))
                (2 (/ (sqrt 15) 30))
                (1 (/ (sqrt 6) 6))
                (0 1)
                (1 (/ (sqrt 6) 6))
                (2 (/ (sqrt 15) 30))
                (3 (/ (sqrt 10) 60))
            )
        )
        (4 
            (ecase (abs degree)
                (0 1)
                (1 (/ (sqrt 10) 10))
                (2 (/ (sqrt 5) 30))
                (3 (/ (sqrt 70) 420))
                (4 (/ (sqrt 35) 840))
            )
        )
        (5
            (ecase (abs degree)
                (0 1)
                (1 (/ (sqrt 15) 15))
                (2 (/ (sqrt 105) 210))
                (3 (/ (sqrt 70) 840))
                (4 (/ (sqrt 35) 2520))
                (5 (/ (sqrt 14) 5040))
            )
        )
        (6
            (ecase (abs degree)
                (0 1)
                (1 (/ (sqrt 21) 21))
                (2 (/ (sqrt 210) 420))
                (3 (/ (sqrt 210) 2520))
                (4 (/ (sqrt 7) 2520))
                (5 (/ (sqrt 154) 55440))
                (6 (/ (sqrt 462) 332640))
            )
        )
        (7
            (ecase (abs degree)
                (0 1)
                (1 (/ (sqrt 7) 14))
                (2 (/ (sqrt 42) 252))
                (3 (/ (sqrt 21) 1260))
                (4 (/ (sqrt 231) 27720))
                (5 (/ (sqrt 231) 166320))
                (6 (/ (sqrt 6006) 4324320))
                (7 (/ (sqrt 429) 4324320))
            )
        )
    )
)

(defun sox-hoaencode-azimuth-factor (degree theta_deg)
"This is a placeholder for documentation"
    (let* 
        (
            (theta (sox-hoaencode-deg-to-rad theta_deg))
@@ -69,6 +116,7 @@
)

(defun sox-hoaencode-elevation-factor (order degree phi_deg)
"This is a placeholder for documentation"
    (let* 
        (
            (phi (sox-hoaencode-deg-to-rad phi_deg))
@@ -83,7 +131,7 @@
            (2 
                (ecase (abs degree)
                    (0 (- (/ (* 3 (expt (sin phi) 2)) 2) (/ 1 2)))
                    (1 (/ (* 3 (expt (sin phi) 2)) 2))
                    (1 (/ (* 3 (* (sin phi) 2)) 2))
                    (2 (* 3 (expt (cos phi) 2)))
                ))
            (3 
@@ -93,45 +141,54 @@
                    (2 (* -15 (sin phi) (- (expt (sin phi) 2) 1)))
                    (3 (* 15 (expt (cos phi) 3)))
                ))
            (4 
                (ecase (abs degree)
                    (0 (+ (/ (* 35 (expt (cos phi) 4)) 8) (- (* 5 (expt (cos phi) 2))) 1))
                    (1 (- (/ (* (sin (* 2 phi)) (- (/ (* 35 (expt (cos phi) 2)) 2) 10)) 2)))
                    (2 (- (/ (* (* 15 (expt (cos phi) 2)) (- (* 7 (expt (cos phi) 2)) 6)) 2)))
                    (3 (* 105 (expt (cos phi) 3) (sin phi) ))
                    (4 (* 105 (expt (cos phi) 4)))
                ))
            (5 
                (ecase (abs degree)
                    (0    (/ (*     (sin phi) (+ (- (* 63 (expt (sin phi) 4)) (* 70 (expt (sin phi) 2))) 15)) 8))
                    (1    (/ (*  15 (sin phi) (+ (- (* 21 (expt (cos phi) 4)) (* 28 (expt (cos phi) 2)))  8)) 8))
                    (2 (- (/ (* 105 (sin phi) (+ (- (*  3 (expt (sin phi) 4)) (*  4 (expt (sin phi) 2)))  1)) 2)))
                    (3 (- (/ (* 105 (expt (cos phi) 3) (- (* 9 (expt (cos phi) 2)) 8))                        2)))
                    (4       (* 945 (sin phi) (expt (- (expt (sin phi) 2) 1) 2)))
                    (5       (* 945 (expt (cos phi) 5)))
                ))
            (6 
                (ecase (abs degree)
                    (0 (+ (- (/ (* 231 (expt (cos phi) 6)) 16)) (/ (* 189 (expt (cos phi) 4)) 8) (- (/ (* 21 (expt (cos phi) 2)) 2)) 1))
                    (1 (/ (* 21 (sin (* 2 phi)) (+ (* 33 (expt (sin phi) 4)) (- (* 30 (expt (sin phi) 2))) 5)) 16))
                    (2 (/ (* 105 (expt (cos phi) 2) (+ (* 33 (expt (cos phi) 4)) (- (* 48 (expt (cos phi) 2))) 16)) 8))
                    (3 (- (/ (* 315 (expt (cos phi) 3) (sin phi)) (- (* 11 (expt (cos phi) 2)) 8)) 2))
                    (4 (- (/ (* 945 (expt (cos phi) 4) (- (* 11 (expt (cos phi) 2)) 10)) 2)))
                    (5 (* 10395 (expt (cos phi) 5) (sin phi)))
                    (6 (* 10395 (expt (cos phi) 6)))
                ))
            (7
                (ecase (abs degree)
                    (0    (/ (* (sin phi) (+ (* 429 (expt (sin phi) 6)) (- (* 693 (expt (sin phi) 4))) (* 315 (expt (sin phi) 2)) (- 35))) 16))
                    (1 (- (/ (* (* 7 (cos phi)) (+ (* 429 (expt (cos phi) 6)) (- (* 792 (expt (cos phi) 4))) (* 432 (expt (cos phi) 2)) (- 64))) 16)))
                    (2 (- (/ (* 63 (sin phi) (- (expt (sin phi) 2) 1) (+ (* 143 (expt (sin phi) 4)) (- (* 110 (expt (sin phi) 2))) 15)) 8)))
                    (3    (/ (* (* 315 (expt (cos phi) 3)) (+ (* 143 (expt (cos phi) 4)) (- (* 220 (expt (cos phi) 2))) 80)) 8))
                    (4    (/ (* 3465 (sin phi) (- (* 13 (expt (sin phi) 2)) 3) (expt (- (expt (sin phi) 2) 1) 2)) 2))
                    (5 (- (/ (* 10395 (expt (cos phi) 5) (- (* 13 (expt (cos phi) 2)) 12)) 2)))
                    (6 (* -135135 (sin phi) (expt (- (expt (sin phi) 2) 1) 3)))
                    (7 (* 135135 (expt (cos phi) 7)))
                ))
        )
    ))

; Convenience functions

; Ensures that positions are given as '((azimuth elevation) ...) (in degrees), i.e.
; 1) if list of list of three values, input is assumed to be '(x y z) coordinates and will be transformed to '(azimuth elevation) coordinates (in the navigational spherical coordinate system).
; 2) if list of list of two values, input is assumed to be '(azimuth elevation) coordinates and won't be transformed any further.
; 3) if 3dc, input is transformed to '(azimuth elevation) coordinates.
(defun sox-hoaencode-auto-convert-positions (positions)
    (cond 
        ((subtypep (type-of positions) '3dc)  
            (progn ; convert xyz->aed, keep azimuth and elevation only.
                (mat-trans 
                    (butlast 
                        (mat-trans 
                            (mapcar 
                                (lambda (xyz) (multiple-value-list (om:xyz->aed (first xyz) (second xyz) (third xyz)))) 
                                (mat-trans (list (x-points positions) (y-points positions) (z-points positions)))
                            ))))))
        ((subtypep (type-of positions) 'list) 
            (cond
                ((= 2 (length (first positions)))
                    positions)
                ((= 3 (length (first positions))) 
                    (mat-trans 
                        (butlast 
                            (mat-trans 
                                (mapcar 
                                    (lambda (xyz) (multiple-value-list (om:xyz->aed (first xyz) (second xyz) (third xyz)))) 
                                    positions
                                )))))
                (t (error "Positions must consist of lists of length 2 (ae) or 3 (xyz)."))
            ))
        (t (error "Positions must be of type 3dc or list."))))

; High-level API
; Returns the gain value (linear, -1 to 1) for a single ACN-channel
; ####### High-level API ########## 


(defun sox-hoaencode-gain-single-component (order degree azimuth_deg elevation_deg)
  "Returns the gain value (linear, -1 to 1) for a single ACN-channel"
    (let
        (
            ; It is assumed that azimuth_deg follows the implementation details of SpatDIF, 
@@ -146,13 +203,18 @@
            (sox-hoaencode-elevation-factor order degree elevation_deg)))
    )

; Returns the gain values for all components at a specific order
; (sox-hoaencode-gain-single-component 1 1 45 0)
; @AN: what is the difference between "order" and "degree"?


(defun sox-hoaencode-gains-by-order (order azimuth_deg elevation_deg)
  "Returns the gain values for all components at a specific order"
    (loop for degree from (* -1 order) to order collect 
        (sox-hoaencode-gain-single-component order degree azimuth_deg elevation_deg)))

; Returns the gain values for all components up to a specific order

(defun sox-hoaencode-gains-up-to-order (order azimuth_deg elevation_deg)
  "Returns the gain values for all components up to a specific order"
    (sox-hoaencode-double-to-float 
        (flatten
            (loop for ord from 0 to order collect 
@@ -161,11 +223,13 @@

(defclass! sox-hoaencode (sox-input)
    (
        (positions :accessor positions :initarg :positions :initform nil :documentation *sox-hoaencode-positions-doc*)
        (order :accessor order :initarg :order :initform nil :documentation *sox-hoaencode-order-doc*)
     (gains :accessor gains :initarg :gains :initform nil :documentation *sox-gain-doc*)
     (azimuth :accessor azimuth :initarg :azimuth :initform 0 :documentation "azimuth angle in degrees")
     (elevation :accessor elevation :initarg :elevation :initform 0 :documentation "elevation angle in degrees")
     (order :accessor order :initarg :order :initform 3 :documentation *sox-hoaencode-order-doc*)
    )
    (:icon 100)
    (:documentation "Sox-hoaencode transforms <sound> into a <order>-th ambisonic (HOA) signal at <positions>.
    (:documentation "Sox-hoaencode encodes <sound> into a <order>-th ambisonic (HOA) signal at <positions>.

    The signal follows the ambiX convention, i.e. it uses SN3D normalization and ACN channel ordering. The resulting file has (<order>+1)^2 channels (order 0: 1 channel, order 1: 4 channels, order 2: 9 channels, order 3: 16 channels, ...).")
)
@@ -173,6 +237,8 @@
(defmethod initialize-instance :after ((self sox-hoaencode) &rest l)
    (declare (ignore l))
    (when (sound self)
        (sox-init-sound self 'list)
        (sox-init-sound self 'atom)
    )
)