#### RemiAudio
#### Copyright (C) 2022-2024 Remilia Scarlet <remilia@posteo.jp>
####   Based on Zita-Rev1
####   Copyright (C) 2003-2017 Fons Adriaensen <fons@linuxaudio.org>
####
#### This program is free software; you can redistribute it and/or modify it
#### under the terms of the GNU General Public License as published by the Free
#### Software Foundation; either version 3 of the License, or (at your option)
#### any later version.
####
#### This program is distributed in the hope that it will be useful, but WITHOUT
#### ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
#### FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License for
#### more details.
####
#### You should have received a copy of the GNU General Public License along
#### with this program.  If not, see <http://www.gnu.org/licenses/>.
require "./reverb"

module RemiAudio::DSP
  class ZitaReverb < Reverb
    # A parametric equalizer specifically for `ZitaReverb` and
    # `ZitaReverbSimple`.
    private class ParametricEQ
      enum Mode
        Bypass = 1
        Static
        Smooth
        MaxCh
      end

      @touch0 : Int16 = 0
      @touch1 : Int16 = 0
      @state  : Mode = Mode::Bypass
      @bypass : Bool = false
      @sampleRate : Float64 = 0.0
      @invSampleRate : Float64 = 0.0
      getter g0 : Float64 = 1.0
      getter g1 : Float64 = 1.0
      getter f0 : Float64 = 1.0e3
      getter f1 : Float64 = 1.0e3
      @c1 : Float64 = 0.0
      @dc1 : Float64 = 0.0
      @c2 : Float64 = 0.0
      @dc2 : Float64 = 0.0
      @gg : Float64 = 0.0
      @dgg : Float64 = 0.0
      @z1 : Array(Float64) = Array.new(PARA_EQ_MAX_CH, 0.0f64)
      @z2 : Array(Float64) = Array.new(PARA_EQ_MAX_CH, 0.0f64)

      def initialize
        self.sampleRate = 44100.0
      end

      @[AlwaysInline]
      def sampleRate=(@sampleRate : Float64) : Nil
        @invSampleRate = 1.0 / @sampleRate
        reset
      end

      @[AlwaysInline]
      def setParam(freq : Float64, gain : Float64) : Nil
        @f0 = freq
        @g0 = 10.0f64 ** (0.05 * gain)
        @touch0 += 1
      end

      @[AlwaysInline]
      def reset : Nil
        @z1.fill(0.0f64)
        @z2.fill(0.0f64)
      end

      def prepare(numFrames : Int) : Nil
        update : Bool = false
        g : Float64 = 0.0
        f : Float64 = 0.0

        unless @touch0 == @touch1
          g = @g0
          f = @f0

          unless g == @g1
            update = true
            if g > (@g1 + @g1)
              @g1 *= 2.0
            elsif @g1 > (g + g)
              @g1 *= 0.5
            else
              @g1 = g
            end
          end

          unless f == @f1
            update = true
            if f > (@f1 + @f1)
              @f1 *= 2.0
            elsif @f1 > (f + f)
              @f1 *= 0.5
            else
              @f1 = f
            end
          end

          if update
            if @state.bypass? && @g1 == 1.0
              calculate(0, @g1, @f1)
            else
              @state = Mode::Smooth
              calculate(numFrames, @g1, @f1)
            end
          else
            @touch1 = @touch0
            if (@g1 - 1).abs < 0.001
              @state = Mode::Bypass
              reset
            else
              @state = Mode::Static
            end
          end
        end
      end

      @[AlwaysInline]
      def calculate(numFrames : Int, gain : Float64, freq : Float64) : Nil
        freq = freq * Math::PI * @invSampleRate
        b : Float64 = (2 * freq) / Math.sqrt(gain)
        gg : Float64 = 0.5 * (gain - 1.0)
        c1 : Float64 = -RemiMath.fastCos(2 * freq)
        c2 : Float64 = (1 - b) / (1 + b)

        if numFrames > 0
          @dc1 = (c1 - @c1) / numFrames + 1.0e-30f64
          @dc2 = (c2 - @c2) / numFrames + 1.0e-30f64
          @dgg = (gg - @gg) / numFrames + 1.0e-30f64
        else
          @c1 = c1
          @c2 = c2
          @gg = gg
        end
      end

      private macro defineProcessFuncs(numType, fromDataFn, toDataFn)
        {% if (toDataFn && !fromDataFn) || (!toDataFn && fromDataFn) %}
          {% raise "fromDataFn and toDataFn must both be specified, or neither" %}
        {% end %}

        private def processSmooth(dataLeft : Array({{numType}})|Slice({{numType}}),
                                  dataRight : Array({{numType}})|Slice({{numType}})) : Nil
          c1 : Float64 = 0.0
          c2 : Float64 = 0.0
          gg : Float64 = 0.0
          z1 : Float64 = 0.0
          z2 : Float64 = 0.0
          x  : Float64 = 0.0
          y  : Float64 = 0.0

          # Left...
          z1 = @z1.unsafe_fetch(0)
          z2 = @z2.unsafe_fetch(0)
          c1 = @c1
          c2 = @c2
          gg = @gg

          dataLeft.size.times do |j|
            c1 += @dc1
            c2 += @dc2
            gg += @dgg

            {% if fromDataFn %}
              x = dataLeft.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = dataLeft.unsafe_fetch(j)
            {% end %}

            y = x - (c2 * z2)

            {% if toDataFn %}
              dataLeft.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              dataLeft.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= c1 * z1
            z2 = z1 + c1 * y
            z1 = y + 1.0e-20f64
          end
          @z1.unsafe_put(0, z1)
          @z2.unsafe_put(0, z2)

          # Right...
          z1 = @z1.unsafe_fetch(1)
          z2 = @z2.unsafe_fetch(1)
          c1 = @c1
          c2 = @c2
          gg = @gg
          dataRight.size.times do |j|
            c1 += @dc1
            c2 += @dc2
            gg += @dgg

            {% if fromDataFn %}
              x = dataRight.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = dataRight.unsafe_fetch(j)
            {% end %}

            x = dataRight[j]
            y = x - c2 * z2

            {% if toDataFn %}
              dataRight.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              dataRight.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= c1 * z1
            z2 = z1 + c1 * y
            z1 = y + 1.0e-20f64
          end
          @z1.unsafe_put(1, z1)
          @z2.unsafe_put(1, z2)

          # Finish up
          @c1 = c1
          @c2 = c2
          @gg = gg
        end

        private def processNonSmooth(dataLeft : Array({{numType}})|Slice({{numType}}),
                                     dataRight : Array({{numType}})|Slice({{numType}})) : Nil
          c1 : Float64 = 0.0
          c2 : Float64 = 0.0
          gg : Float64 = 0.0
          z1 : Float64 = 0.0
          z2 : Float64 = 0.0
          x  : Float64 = 0.0
          y  : Float64 = 0.0

          # Left...
          z1 = @z1.unsafe_fetch(0)
          z2 = @z2.unsafe_fetch(0)

          dataLeft.size.times do |j|
            {% if fromDataFn %}
              x = dataLeft.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = dataLeft.unsafe_fetch(j)
            {% end %}

            y = x - (c2 * z2)

            {% if toDataFn %}
              dataLeft.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              dataLeft.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= c1 * z1
            z2 = z1 + (c1 * y)
            z1 = y + 1.0e-20f64
          end
          @z1.unsafe_put(0, z1)
          @z2.unsafe_put(0, z2)

          # Right...
          z1 = @z1.unsafe_fetch(1)
          z2 = @z2.unsafe_fetch(1)
          dataRight.size.times do |j|
            {% if fromDataFn %}
              x = dataRight.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = dataRight.unsafe_fetch(j)
            {% end %}

            y = x - (c2 * z2)
            {% if toDataFn %}
              dataRight.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              dataRight.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= (c1 * z1)
            z2 = z1 + (c1 * y)
            z1 = y + 1.0e-20f64
          end
          @z1.unsafe_put(1, z1)
          @z2.unsafe_put(1, z2)
        end

        # Interleaved version
        private def processSmooth(data : Array({{numType}})|Slice({{numType}})) : Nil
          c1 : Float64 = 0.0
          c2 : Float64 = 0.0
          gg : Float64 = 0.0
          z1 : Float64 = 0.0
          z2 : Float64 = 0.0
          x  : Float64 = 0.0
          y  : Float64 = 0.0
          j = 0

          # Left...
          z1 = @z1.unsafe_fetch(0)
          z2 = @z2.unsafe_fetch(0)
          c1 = @c1
          c2 = @c2
          gg = @gg
          while j < data.size
            c1 += @dc1
            c2 += @dc2
            gg += @dgg

            {% if fromDataFn %}
              x = data.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = data.unsafe_fetch(j)
            {% end %}

            y = x - (c2 * z2)

            {% if toDataFn %}
              data.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              data.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= c1 * z1
            z2 = z1 + c1 * y
            z1 = y + 1.0e-20f64
            j += 2
          end
          @z1.unsafe_put(0, z1)
          @z2.unsafe_put(0, z2)

          # Right...
          z1 = @z1.unsafe_fetch(1)
          z2 = @z2.unsafe_fetch(1)
          c1 = @c1
          c2 = @c2
          gg = @gg
          j = 1
          while j < data.size
            c1 += @dc1
            c2 += @dc2
            gg += @dgg

            {% if fromDataFn %}
              x = data.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = data.unsafe_fetch(j)
            {% end %}

            y = x - c2 * z2

            {% if toDataFn %}
              data.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              data.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= c1 * z1
            z2 = z1 + c1 * y
            z1 = y + 1.0e-20f64
            j += 2
          end
          @z1.unsafe_put(1, z1)
          @z2.unsafe_put(1, z2)

          # Finish up
          @c1 = c1
          @c2 = c2
          @gg = gg
        end

        # Interleaved version
        private def processNonSmooth(data : Array({{numType}})|Slice({{numType}})) : Nil
          c1 : Float64 = 0.0
          c2 : Float64 = 0.0
          gg : Float64 = 0.0
          z1 : Float64 = 0.0
          z2 : Float64 = 0.0
          x  : Float64 = 0.0
          y  : Float64 = 0.0
          j = 0

          # Left...
          z1 = @z1.unsafe_fetch(0)
          z2 = @z2.unsafe_fetch(0)
          while j < data.size
            {% if fromDataFn %}
              x = data.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = data.unsafe_fetch(j)
            {% end %}

            y = x - (c2 * z2)

            {% if toDataFn %}
              data.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              data.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= c1 * z1
            z2 = z1 + (c1 * y)
            z1 = y + 1.0e-20f64
            j += 2
          end
          @z1.unsafe_put(0, z1)
          @z2.unsafe_put(0, z2)

          # Right...
          z1 = @z1.unsafe_fetch(1)
          z2 = @z2.unsafe_fetch(1)
          j = 1
          while j < data.size
            {% if fromDataFn %}
              x = data.unsafe_fetch(j).{{fromDataFn.id}}
            {% else %}
              x = data.unsafe_fetch(j)
            {% end %}

            y = x - (c2 * z2)

            {% if toDataFn %}
              data.unsafe_put(j, (x - gg * (z2 + c2 * y - x)).{{toDataFn.id}})
            {% else %}
              data.unsafe_put(j, x - gg * (z2 + c2 * y - x))
            {% end %}

            y -= (c1 * z1)
            z2 = z1 + (c1 * y)
            z1 = y + 1.0e-20f64
            j += 2
          end
          @z1.unsafe_put(1, z1)
          @z2.unsafe_put(1, z2)
        end

      end

      defineProcessFuncs(Float64, nil, nil)
      defineProcessFuncs(Float32, :to_f64!, :to_f32!)

      {% begin %}
        {% for typ in [:Float64, :Float32] %}
          @[AlwaysInline]
          def process(dataLeft : Array({{typ.id}})|Slice({{typ.id}}),
                      dataRight : Array({{typ.id}})|Slice({{typ.id}})) : Nil
            if @state.smooth?
              processSmooth(dataLeft, dataRight)
            else
              processNonSmooth(dataLeft, dataRight)
            end
          end

          @[AlwaysInline]
          def process(data : Array({{typ.id}})|Slice({{typ.id}})) : Nil
            if @state.smooth?
              processSmooth(data)
            else
              processNonSmooth(data)
            end
          end

          @[AlwaysInline]
          def process(dataLeft : {{typ.id}}, dataRight : {{typ.id}}) : Tuple({{typ.id}}, {{typ.id}})
            left = Slice({{typ.id}}).new(1, dataLeft)
            right = Slice({{typ.id}}).new(1, dataRight)

            if @state.smooth?
              processSmooth(dataLeft, left)
            else
              processNonSmooth(dataLeft, right)
            end

            {left.unsafe_fetch(0), right.unsafe_fetch(0)}
          end
        {% end %}
      {% end %}
    end
  end
end
