[AT3P] Add GHA processing in to the bitstream encoding.

This comit bring some glue to extract residual.
Also we need additional one frame delay after pqf.
Temporary use ffmpeg code for gha synth.

9 files changed, 553 insertions, 27 deletions

diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index eba4332..47c8cc9 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -93,6 +93,7 @@ set(SOURCE_ATRACDENC_IMPL
     atrac/atrac3.cpp
     atrac/atrac3_bitstream.cpp
     atrac/atrac3plus_pqf/atrac3plus_pqf.c
+    atrac/at3p/ff/atrac3plusdsp.c
     atrac/at3p/at3p.cpp
     atrac/at3p/at3p_bitstream.cpp
     atrac/at3p/at3p_gha.cpp
diff --git a/src/atrac/at3p/at3p.cpp b/src/atrac/at3p/at3p.cpp
index 514fa96..278d362 100644
--- a/src/atrac/at3p/at3p.cpp
+++ b/src/atrac/at3p/at3p.cpp
@@ -39,7 +39,9 @@ public:
         : BitStream(out, 2048)
         , ChannelCtx(channels)
         , GhaProcessor(MakeGhaProcessor0(channels == 2))
-    {}
+    {
+        delay.NumToneBands = 0;
+    }
 
     TPCMEngine::EProcessResult EncodeFrame(const float* data, int channels);
 private:
@@ -59,6 +61,7 @@ private:
         float* CurBuf = nullptr;
         float Buf1[TAt3PEnc::NumSamples];
         float Buf2[TAt3PEnc::NumSamples];
+        float PrevBuf[TAt3PEnc::NumSamples];
         TAt3pMDCT::THistBuf MdctBuf;
         std::vector<float> Specs;
     };
@@ -68,6 +71,7 @@ private:
     TAt3PBitStream BitStream;
     vector<TChannelCtx> ChannelCtx;
     std::unique_ptr<IGhaProcessor> GhaProcessor;
+    TAt3PGhaData delay;
 };
 
 TPCMEngine::EProcessResult TAt3PEnc::TImpl::
@@ -95,21 +99,30 @@ EncodeFrame(const float* data, int channels)
 
     assert(needMore == 0);
 
+    float* b1Prev = ChannelCtx[0].PrevBuf;
     const float* b1Cur = ChannelCtx[0].CurBuf;
     const float* b1Next = ChannelCtx[0].NextBuf;
+    float* b2Prev = (channels == 2) ? ChannelCtx[1].PrevBuf : nullptr;
     const float* b2Cur = (channels == 2) ? ChannelCtx[1].CurBuf : nullptr;
     const float* b2Next = (channels == 2) ? ChannelCtx[1].NextBuf : nullptr;
 
-    auto tonalBlock = GhaProcessor->DoAnalize({b1Cur, b1Next}, {b2Cur, b2Next});
+
+    const TAt3PGhaData* p = nullptr;
+    if (delay.NumToneBands) {
+        p = &delay;
+    }
+
+    const TAt3PGhaData* tonalBlock = GhaProcessor->DoAnalize({b1Cur, b1Next}, {b2Cur, b2Next}, b1Prev, b2Prev);
 
     std::vector<std::vector<TScaledBlock>> scaledBlocks;
     for (int ch = 0; ch < channels; ch++) {
+        float* x = (ch == 0) ? b1Prev : b2Prev;
         auto& c = ChannelCtx[ch];
         TAt3pMDCT::TPcmBandsData p;
         float tmp[2048];
         //TODO: scale window
         for (size_t i = 0; i < 2048; i++) {
-            tmp[i] = c.CurBuf[i] / 32768.0;
+            tmp[i] = x[i] / 32768.0;
         }
         for (size_t b = 0; b < 16; b++) {
             p[b] = tmp + b * 128;
@@ -120,11 +133,17 @@ EncodeFrame(const float* data, int channels)
         scaledBlocks.push_back(block);
     }
 
-    BitStream.WriteFrame(channels, tonalBlock, scaledBlocks);
+    BitStream.WriteFrame(channels, p, scaledBlocks);
 
     for (int ch = 0; ch < channels; ch++) {
+        memcpy(ChannelCtx[ch].PrevBuf, ChannelCtx[ch].CurBuf, sizeof(float) * TAt3PEnc::NumSamples);
         std::swap(ChannelCtx[ch].NextBuf, ChannelCtx[ch].CurBuf);
     }
+    if (tonalBlock) {
+        delay = *tonalBlock;
+    } else {
+        delay.NumToneBands = 0;
+    }
 
     return TPCMEngine::EProcessResult::PROCESSED;
 }
diff --git a/src/atrac/at3p/at3p_bitstream.cpp b/src/atrac/at3p/at3p_bitstream.cpp
index 444ace2..f00fe31 100644
--- a/src/atrac/at3p/at3p_bitstream.cpp
+++ b/src/atrac/at3p/at3p_bitstream.cpp
@@ -537,7 +537,7 @@ static void WriteTonalBlock(NBitStream::TBitStream& bs, int channels, const TAt3
     }
 }
 
-void TAt3PBitStream::WriteFrame(int channels, const TAt3PGhaData* /*tonalBlock*/, const std::vector<std::vector<TScaledBlock>>& scaledBlocks)
+void TAt3PBitStream::WriteFrame(int channels, const TAt3PGhaData* tonalBlock, const std::vector<std::vector<TScaledBlock>>& scaledBlocks)
 {
     NBitStream::TBitStream bitStream;
     // First bit must be zero
@@ -564,7 +564,13 @@ void TAt3PBitStream::WriteFrame(int channels, const TAt3PGhaData* /*tonalBlock*/
         bitStream.Write(0, 1); //gain comp
     }
 
-    bitStream.Write(0, 1);
+    if (tonalBlock && tonalBlock->NumToneBands) {
+        bitStream.Write(1, 1);
+        WriteTonalBlock(bitStream, channels, tonalBlock);
+    } else {
+        bitStream.Write(0, 1);
+    }
+
     bitStream.Write(0, 1); // no noise info
     // Terminator
     bitStream.Write(3, 2);
diff --git a/src/atrac/at3p/at3p_bitstream_impl.h b/src/atrac/at3p/at3p_bitstream_impl.h
index 4b4ea6b..52a86a5 100644
--- a/src/atrac/at3p/at3p_bitstream_impl.h
+++ b/src/atrac/at3p/at3p_bitstream_impl.h
@@ -48,7 +48,6 @@ struct TSpecFrame {
             : ScaledBlocks(scaledBlock)
         {}
         const std::vector<TScaledBlock>& ScaledBlocks;
-        int Mant[2048];
     };
 
     std::vector<TChannel> Chs;
diff --git a/src/atrac/at3p/at3p_gha.cpp b/src/atrac/at3p/at3p_gha.cpp
index 96b3eba..2f20940 100644
--- a/src/atrac/at3p/at3p_gha.cpp
+++ b/src/atrac/at3p/at3p_gha.cpp
@@ -17,6 +17,7 @@
  */
 
 #include "at3p_gha.h"
+#include "ff/atrac3plus.h"
 
 #include <assert.h>
 #include <atrac/atrac_psy_common.h>
@@ -115,6 +116,8 @@ public:
         gha_set_max_magnitude(LibGhaCtx, 32768);
 
         if (!StaticInited) {
+            ff_atrac3p_init_dsp_static();
+
             FillSubbandAth(&SubbandAth[0]);
 
             AmpSfTab = CreateAmpSfTab();
@@ -124,6 +127,16 @@ public:
 
             StaticInited = true;
         }
+
+        for (size_t ch = 0; ch < 2; ch++) {
+           ChUnit.channels[ch].tones_info      = &ChUnit.channels[ch].tones_info_hist[0][0];
+           ChUnit.channels[ch].tones_info_prev = &ChUnit.channels[ch].tones_info_hist[1][0];
+        }
+
+        ChUnit.waves_info      = &ChUnit.wave_synth_hist[0];
+        ChUnit.waves_info_prev = &ChUnit.wave_synth_hist[1];
+        ChUnit.waves_info->tones_present = false;
+        ChUnit.waves_info_prev->tones_present = false;
     }
 
     ~TGhaProcessor()
@@ -131,8 +144,10 @@ public:
         gha_free_ctx(LibGhaCtx);
     }
 
-    const TAt3PGhaData* DoAnalize(TBufPtr b1, TBufPtr b2) override;
+    const TAt3PGhaData* DoAnalize(TBufPtr b1, TBufPtr b2, float *w1, float *w2) override;
+
 private:
+    void ApplyFilter(const TAt3PGhaData*, float *b1, float *b2);
     static void FillSubbandAth(float* out);
     static TAmpSfTab CreateAmpSfTab();
     static void CheckResuidalAndApply(float* resuidal, size_t size, void* self) noexcept;
@@ -151,12 +166,15 @@ private:
     gha_ctx_t LibGhaCtx;
     TAt3PGhaData ResultBuf;
     TAt3PGhaData ResultBufHistory;
+
     const bool Stereo;
 
     static float SubbandAth[SUBBANDS];
     static float SineTab[2048];
     static bool StaticInited;
     static TAmpSfTab AmpSfTab;
+
+    Atrac3pChanUnitCtx ChUnit;
 };
 
 bool TGhaProcessor::StaticInited = false;
@@ -292,7 +310,115 @@ void TGhaProcessor::CheckResuidalAndApply(float* resuidal, size_t size, void* d)
     memcpy(b, resuidal, sizeof(float) * SAMPLES_PER_SUBBAND);
 }
 
-const TAt3PGhaData* TGhaProcessor::DoAnalize(TBufPtr b1, TBufPtr b2)
+void TGhaProcessor::ApplyFilter(const TAt3PGhaData* d, float* b1, float* b2)
+{
+    for (size_t ch_num = 0; ch_num < 2; ch_num++)
+        memset(ChUnit.channels[ch_num].tones_info, 0,
+           sizeof(*ChUnit.channels[ch_num].tones_info) * ATRAC3P_SUBBANDS);
+
+    if (d) {
+        memset(ChUnit.waves_info->waves, 0, sizeof(ChUnit.waves_info->waves));
+        ChUnit.waves_info->num_tone_bands = d->NumToneBands;
+        ChUnit.waves_info->tones_present = true;
+        ChUnit.waves_info->amplitude_mode = 1;
+        ChUnit.waves_info->tones_index = 0;
+    } else {
+        ChUnit.waves_info->tones_present = false;
+    }
+
+    if (d) {
+        for (size_t ch = 0; ch <= (size_t)Stereo; ch++) {
+            for (int i = 0; i < ChUnit.waves_info->num_tone_bands; i++) {
+                if (ch && d->ToneSharing[i]) {
+                    continue;
+                }
+
+                ChUnit.channels[ch].tones_info[i].num_wavs = d->GetNumWaves(ch, i);
+
+                const auto envelope = d->GetEnvelope(ch, i);
+                if (envelope.first != TAt3PGhaData::EMPTY_POINT) {
+                    // start point present
+                    ChUnit.channels[ch].tones_info[i].pend_env.has_start_point = true;
+                    ChUnit.channels[ch].tones_info[i].pend_env.start_pos = envelope.first;
+                } else {
+                    ChUnit.channels[ch].tones_info[i].pend_env.has_start_point = false;
+                    ChUnit.channels[ch].tones_info[i].pend_env.start_pos = -1;
+                }
+
+                if (envelope.second != TAt3PGhaData::EMPTY_POINT) {
+                    // stop point present
+                    ChUnit.channels[ch].tones_info[i].pend_env.has_stop_point  = true;
+                    ChUnit.channels[ch].tones_info[i].pend_env.stop_pos = envelope.second;
+                } else {
+                    ChUnit.channels[ch].tones_info[i].pend_env.has_stop_point  = false;
+                    ChUnit.channels[ch].tones_info[i].pend_env.stop_pos = 32;
+                }
+            }
+
+            for (int sb = 0; sb < ChUnit.waves_info->num_tone_bands; sb++) {
+                if (d->GetNumWaves(ch, sb)) {
+                    if (ChUnit.waves_info->tones_index + ChUnit.channels[ch].tones_info[sb].num_wavs > 48) {
+                        std::cerr << "too many tones: " << ChUnit.waves_info->tones_index + ChUnit.channels[ch].tones_info[sb].num_wavs << std::endl;
+                        abort();
+                    }
+                    ChUnit.channels[ch].tones_info[sb].start_index           = ChUnit.waves_info->tones_index;
+                    ChUnit.waves_info->tones_index += ChUnit.channels[ch].tones_info[sb].num_wavs;
+                }
+            }
+
+            Atrac3pWaveParam *iwav;
+            for (int sb = 0; sb < ChUnit.waves_info->num_tone_bands; sb++) {
+                if (d->GetNumWaves(ch, sb)) {
+                    iwav = &ChUnit.waves_info->waves[ChUnit.channels[ch].tones_info[sb].start_index];
+                    auto w = d->GetWaves(ch, sb);
+                    ChUnit.channels[ch].tones_info[sb].num_wavs = w.second;
+                    for (size_t j = 0; j < w.second; j++) {
+                        iwav[j].freq_index = w.first[j].FreqIndex;
+                        iwav[j].amp_index = w.first[j].AmpIndex;
+                        iwav[j].amp_sf = w.first[j].AmpSf;
+                        iwav[j].phase_index = w.first[j].PhaseIndex;
+                    }
+                }
+            }
+        }
+
+        if (Stereo) {
+            for (int i = 0; i < ChUnit.waves_info->num_tone_bands; i++) {
+                if (d->ToneSharing[i]) {
+                    ChUnit.channels[1].tones_info[i] = ChUnit.channels[0].tones_info[i];
+                }
+
+                if (d->SecondIsLeader) {
+                    std::swap(ChUnit.channels[0].tones_info[i],
+                        ChUnit.channels[1].tones_info[i]);
+                }
+            }
+        }
+    }
+
+    for (size_t ch = 0; ch <= (size_t)Stereo; ch++) {
+        float* x = (ch == 0) ? b1 : b2;
+        if (ChUnit.waves_info->tones_present ||
+            ChUnit.waves_info_prev->tones_present) {
+            for (size_t sb = 0; sb < SUBBANDS; sb++) {
+                if (ChUnit.channels[ch].tones_info[sb].num_wavs ||
+                    ChUnit.channels[ch].tones_info_prev[sb].num_wavs) {
+
+                    ff_atrac3p_generate_tones(&ChUnit, ch, sb,
+                        (x + sb * 128));
+                }
+            }
+        }
+    }
+
+    for (size_t ch = 0; ch <= (size_t)Stereo; ch++) {
+        std::swap(ChUnit.channels[ch].tones_info, ChUnit.channels[ch].tones_info_prev);
+    }
+
+    std::swap(ChUnit.waves_info, ChUnit.waves_info_prev);
+}
+
+const TAt3PGhaData* TGhaProcessor::DoAnalize(TBufPtr b1, TBufPtr b2, float* w1, float* w2)
 {
     vector<TChannelData> data((size_t)Stereo + 1);
     bool progress[2] = {false};
@@ -322,13 +448,17 @@ const TAt3PGhaData* TGhaProcessor::DoAnalize(TBufPtr b1, TBufPtr b2)
     } while ((progress[0] || progress[1]) && totalTones < 48);
 
     if (totalTones == 0) {
+        ApplyFilter(nullptr, w1, w2);
         return nullptr;
     }
 
     FillResultBuf(data);
+
     ResultBufHistory = ResultBuf;
 
-    return &ResultBuf;
+    ApplyFilter(&ResultBuf, w1, w2);
+
+    return  &ResultBuf;
 }
 
 bool TGhaProcessor::CheckNextFrame(const float* nextSrc, const vector<gha_info>& ghaInfos) const
diff --git a/src/atrac/at3p/at3p_gha.h b/src/atrac/at3p/at3p_gha.h
index 45f3cbc..ff708f4 100644
--- a/src/atrac/at3p/at3p_gha.h
+++ b/src/atrac/at3p/at3p_gha.h
@@ -68,7 +68,7 @@ class IGhaProcessor {
 public:
     using TBufPtr = std::array<const float*, 2>;
     virtual ~IGhaProcessor() {}
-    virtual const TAt3PGhaData* DoAnalize(TBufPtr b1, TBufPtr b2) = 0;
+    virtual const TAt3PGhaData* DoAnalize(TBufPtr b1, TBufPtr b2, float* w1, float* w2) = 0;
 };
 
 std::unique_ptr<IGhaProcessor> MakeGhaProcessor0(bool stereo);
diff --git a/src/atrac/at3p/at3p_gha_ut.cpp b/src/atrac/at3p/at3p_gha_ut.cpp
index fcf5095..769c1a0 100644
--- a/src/atrac/at3p/at3p_gha_ut.cpp
+++ b/src/atrac/at3p/at3p_gha_ut.cpp
@@ -29,6 +29,12 @@ static void __attribute__ ((noinline)) Gen(const TTestParam& p, vector<float>& o
     }
 }
 
+static TAt3PGhaData DoAnalize(IGhaProcessor* p, IGhaProcessor::TBufPtr b1, IGhaProcessor::TBufPtr b2) {
+    float w1[2048] = {0};
+    float w2[2048] = {0};
+    return *p->DoAnalize(b1, b2, w1, w2);
+}
+
 static const TAt3PGhaData __attribute__ ((noinline)) GenAndRunGha(vector<TTestParam> p1, vector<TTestParam> p2)
 {
     vector<float> buf1(2048 * 2);
@@ -51,7 +57,7 @@ static const TAt3PGhaData __attribute__ ((noinline)) GenAndRunGha(vector<TTestPa
     const float* b1 = buf1.data();
     const float* b2 = buf2.empty() ? nullptr : buf2.data();
 
-    return *(processor->DoAnalize({b1, b1 + 2048}, {b2, b2 + 2048}));
+    return DoAnalize(processor.get(), {b1, b1 + 2048}, {b2, b2 + 2048});
 }
 
 static class TDumper {
@@ -470,7 +476,7 @@ TEST(AT3PGHA, 100hz__two_frames_mono) {
     auto processor = MakeGhaProcessor0(false);
 
     {
-        const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+        const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
         EXPECT_EQ(res.NumToneBands, 1);
         EXPECT_EQ(res.Waves[0].WaveParams.size(), 1);
@@ -486,7 +492,7 @@ TEST(AT3PGHA, 100hz__two_frames_mono) {
 
     {
         memset(&buf[0], 0, sizeof(float) * 128);
-        const auto& res = *processor->DoAnalize({&buf[2048], &buf[0]}, {nullptr, nullptr});
+        const auto res = DoAnalize(processor.get(), {&buf[2048], &buf[0]}, {nullptr, nullptr});
 
         EXPECT_EQ(res.NumToneBands, 1);
         EXPECT_EQ(res.Waves[0].WaveParams.size(), 1);
@@ -515,7 +521,7 @@ TEST(AT3PGHA, 100hz_than_500hz_than_100hz__3_frames_mono) {
     auto processor = MakeGhaProcessor0(false);
 
     {
-        const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+        const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
         EXPECT_EQ(res.NumToneBands, 1);
         EXPECT_EQ(res.Waves[0].WaveParams.size(), 1);
@@ -532,7 +538,7 @@ TEST(AT3PGHA, 100hz_than_500hz_than_100hz__3_frames_mono) {
     {
         memset(&buf[0], 0, sizeof(float) * 128);
         Gen({100.0f, 0, 32768, 0, 128}, buf);
-        const auto& res = *processor->DoAnalize({&buf[2048], &buf[0]}, {nullptr, nullptr});
+        const auto res = DoAnalize(processor.get(), {&buf[2048], &buf[0]}, {nullptr, nullptr});
 
         EXPECT_EQ(res.NumToneBands, 1);
         EXPECT_EQ(res.Waves[0].WaveParams.size(), 1);
@@ -547,7 +553,7 @@ TEST(AT3PGHA, 100hz_than_500hz_than_100hz__3_frames_mono) {
     }
     {
         memset(&buf[2048], 0, sizeof(float) * 128);
-        const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+        const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
         EXPECT_EQ(res.NumToneBands, 1);
         EXPECT_EQ(res.Waves[0].WaveParams.size(), 1);
@@ -564,7 +570,6 @@ TEST(AT3PGHA, 100hz_than_500hz_than_100hz__3_frames_mono) {
     Dumper.Dump(resBuf.data(), 1, 3);
 }
 
-
 TEST(AT3PGHA, 100hz__phase_two_frames_mono) {
     vector<float> buf(2048 * 2);
 
@@ -574,7 +579,7 @@ TEST(AT3PGHA, 100hz__phase_two_frames_mono) {
     memset(&buf[128], 0, sizeof(float) * 128);
 
     auto processor = MakeGhaProcessor0(false);
-    const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+    const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
     EXPECT_EQ(res.NumToneBands, 1);
     EXPECT_EQ(res.Waves[0].WaveParams.size(), 1);
@@ -586,7 +591,6 @@ TEST(AT3PGHA, 100hz__phase_two_frames_mono) {
     EXPECT_EQ(res.GetWaves(0, 0).first->PhaseIndex, 4);
 }
 
-
 TEST(AT3PGHA, 689hz0625__two_frames_mono) {
     vector<float> buf(2048 * 2);
 
@@ -596,7 +600,7 @@ TEST(AT3PGHA, 689hz0625__two_frames_mono) {
     memset(&buf[128], 0, sizeof(float) * 128);
 
     auto processor = MakeGhaProcessor0(false);
-    const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+    const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
     EXPECT_EQ(res.NumToneBands, 1);
     EXPECT_EQ(res.Waves[0].WaveParams.size(), 1);
@@ -617,7 +621,7 @@ TEST(AT3PGHA, 689hz0625_1000hz__two_frames_mono) {
     memset(&buf[128], 0, sizeof(float) * 128);
 
     auto processor = MakeGhaProcessor0(false);
-    const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+    const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
     EXPECT_EQ(res.NumToneBands, 1);
     EXPECT_EQ(res.Waves[0].WaveParams.size(), 2);
@@ -640,7 +644,7 @@ TEST(AT3PGHA, 500hz_1000hz__two_frames_mono) {
     memset(&buf[128], 0, sizeof(float) * 128);
 
     auto processor = MakeGhaProcessor0(false);
-    const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+    const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
     EXPECT_EQ(res.NumToneBands, 1);
     EXPECT_EQ(res.Waves[0].WaveParams.size(), 2);
@@ -663,7 +667,7 @@ TEST(AT3PGHA, 500hz_1000hz__phase_two_frames_mono) {
     memset(&buf[128], 0, sizeof(float) * 128);
 
     auto processor = MakeGhaProcessor0(false);
-    const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+    const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
     EXPECT_EQ(res.NumToneBands, 1);
     EXPECT_EQ(res.Waves[0].WaveParams.size(), 2);
@@ -677,7 +681,6 @@ TEST(AT3PGHA, 500hz_1000hz__phase_two_frames_mono) {
     Dumper.Dump(&res, 1, 1);
 }
 
-
 TEST(AT3PGHA, 250hz_500hz_1000hz__two_frames_mono) {
     vector<float> buf(2048 * 2);
 
@@ -689,7 +692,7 @@ TEST(AT3PGHA, 250hz_500hz_1000hz__two_frames_mono) {
     memset(&buf[128], 0, sizeof(float) * 128);
 
     auto processor = MakeGhaProcessor0(false);
-    const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+    const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
     EXPECT_EQ(res.NumToneBands, 1);
     EXPECT_EQ(res.Waves[0].WaveParams.size(), 3);
@@ -716,7 +719,7 @@ TEST(AT3PGHA, 250hz_500hz_1000hz_1200hz__two_frames_mono) {
     memset(&buf[128], 0, sizeof(float) * 128);
 
     auto processor = MakeGhaProcessor0(false);
-    const auto& res = *processor->DoAnalize({&buf[0], &buf[2048]}, {nullptr, nullptr});
+    const auto res = DoAnalize(processor.get(), {&buf[0], &buf[2048]}, {nullptr, nullptr});
 
     EXPECT_EQ(res.NumToneBands, 1);
     EXPECT_EQ(res.Waves[0].WaveParams.size(), 4);
diff --git a/src/atrac/at3p/ff/atrac3plus.h b/src/atrac/at3p/ff/atrac3plus.h
new file mode 100644
index 0000000..3aeb110
--- /dev/null
+++ b/src/atrac/at3p/ff/atrac3plus.h
@@ -0,0 +1,165 @@
+/*
+ * ATRAC3+ compatible decoder
+ *
+ * Copyright (c) 2010-2013 Maxim Poliakovski
+ *
+ * This file is part of FFmpeg.
+ *
+ * FFmpeg is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public
+ * License as published by the Free Software Foundation; either
+ * version 2.1 of the License, or (at your option) any later version.
+ *
+ * FFmpeg 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
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with FFmpeg; if not, write to the Free Software
+ * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
+ */
+
+/**
+ * @file
+ * Global structures, constants and data for ATRAC3+ decoder.
+ */
+
+#ifndef AVCODEC_ATRAC3PLUS_H
+#define AVCODEC_ATRAC3PLUS_H
+
+#include <stdint.h>
+
+/** Global unit sizes */
+#define ATRAC3P_SUBBANDS        16  ///< number of PQF subbands
+#define ATRAC3P_SUBBAND_SAMPLES 128 ///< number of samples per subband
+#define ATRAC3P_FRAME_SAMPLES   (ATRAC3P_SUBBAND_SAMPLES * ATRAC3P_SUBBANDS)
+
+#define ATRAC3P_PQF_FIR_LEN     12  ///< length of the prototype FIR of the PQF
+
+/** Global constants */
+#define ATRAC3P_POWER_COMP_OFF  15  ///< disable power compensation
+
+#ifdef __cplusplus
+extern "C"
+{
+#endif /* __cplusplus */
+
+/** ATRAC3+ channel unit types */
+enum Atrac3pChannelUnitTypes {
+    CH_UNIT_MONO       = 0, ///< unit containing one coded channel
+    CH_UNIT_STEREO     = 1, ///< unit containing two jointly-coded channels
+    CH_UNIT_EXTENSION  = 2, ///< unit containing extension information
+    CH_UNIT_TERMINATOR = 3  ///< unit sequence terminator
+};
+
+/** Amplitude envelope of a group of sine waves */
+typedef struct Atrac3pWaveEnvelope {
+    int has_start_point;    ///< indicates start point within the GHA window
+    int has_stop_point;     ///< indicates stop point within the GHA window
+    int start_pos;          ///< start position expressed in n*4 samples
+    int stop_pos;           ///< stop  position expressed in n*4 samples
+} Atrac3pWaveEnvelope;
+
+/** Parameters of a group of sine waves */
+typedef struct Atrac3pWavesData {
+    Atrac3pWaveEnvelope pend_env;   ///< pending envelope from the previous frame
+    Atrac3pWaveEnvelope curr_env;   ///< group envelope from the current frame
+    int num_wavs;           ///< number of sine waves in the group
+    int start_index;        ///< start index into global tones table for that subband
+} Atrac3pWavesData;
+
+/** Parameters of a single sine wave */
+typedef struct Atrac3pWaveParam {
+    int   freq_index;   ///< wave frequency index
+    int   amp_sf;       ///< quantized amplitude scale factor
+    int   amp_index;    ///< quantized amplitude index
+    int   phase_index;  ///< quantized phase index
+} Atrac3pWaveParam;
+
+/** Sound channel parameters */
+typedef struct Atrac3pChanParams {
+    //int ch_num;
+    //int num_coded_vals;         ///< number of transmitted quant unit values
+    //int fill_mode;
+    //int split_point;
+    //int table_type;             ///< table type: 0 - tone?, 1- noise?
+    //int qu_wordlen[32];         ///< array of word lengths for each quant unit
+    //int qu_sf_idx[32];          ///< array of scale factor indexes for each quant unit
+    //int qu_tab_idx[32];         ///< array of code table indexes for each quant unit
+    //int16_t spectrum[2048];     ///< decoded IMDCT spectrum
+    //uint8_t power_levs[5];      ///< power compensation levels
+
+    /* imdct window shape history (2 frames) for overlapping. */
+    //uint8_t wnd_shape_hist[2][ATRAC3P_SUBBANDS];    ///< IMDCT window shape, 0=sine/1=steep
+    //uint8_t *wnd_shape;         ///< IMDCT window shape for current frame
+    //uint8_t *wnd_shape_prev;    ///< IMDCT window shape for previous frame
+
+    //int num_gain_subbands;      ///< number of subbands with gain control data
+
+    /* tones data history (2 frames) for overlapping. */
+    Atrac3pWavesData tones_info_hist[2][ATRAC3P_SUBBANDS];
+    Atrac3pWavesData *tones_info;
+    Atrac3pWavesData *tones_info_prev;
+} Atrac3pChanParams;
+
+/* Per-unit sine wave parameters */
+typedef struct Atrac3pWaveSynthParams {
+    int tones_present;                      ///< 1 - tones info present
+    int amplitude_mode;                     ///< 1 - low range, 0 - high range
+    int num_tone_bands;                     ///< number of PQF bands with tones
+    //uint8_t tone_sharing[ATRAC3P_SUBBANDS]; ///< 1 - subband-wise tone sharing flags
+    //uint8_t tone_master[ATRAC3P_SUBBANDS];  ///< 1 - subband-wise tone channel swapping
+    uint8_t invert_phase[ATRAC3P_SUBBANDS]; ///< 1 - subband-wise phase inversion
+    int tones_index;                        ///< total sum of tones in this unit
+    Atrac3pWaveParam waves[48];
+} Atrac3pWaveSynthParams;
+
+/** Channel unit parameters */
+typedef struct Atrac3pChanUnitCtx {
+    /* channel unit variables */
+    //int unit_type;                          ///< unit type (mono/stereo)
+    //int num_quant_units;
+    //int num_subbands;
+    //int used_quant_units;                   ///< number of quant units with coded spectrum
+    //int num_coded_subbands;                 ///< number of subbands with coded spectrum
+    //int mute_flag;                          ///< mute flag
+    //int use_full_table;                     ///< 1 - full table list, 0 - restricted one
+    //int noise_present;                      ///< 1 - global noise info present
+    //int noise_level_index;                  ///< global noise level index
+    //int noise_table_index;                  ///< global noise RNG table index
+    //uint8_t swap_channels[ATRAC3P_SUBBANDS];    ///< 1 - perform subband-wise channel swapping
+    //uint8_t negate_coeffs[ATRAC3P_SUBBANDS];    ///< 1 - subband-wise IMDCT coefficients negation
+    Atrac3pChanParams channels[2];
+
+    /* Variables related to GHA tones */
+    Atrac3pWaveSynthParams wave_synth_hist[2];     ///< waves synth history for two frames
+    Atrac3pWaveSynthParams *waves_info;
+    Atrac3pWaveSynthParams *waves_info_prev;
+
+    //Atrac3pIPQFChannelCtx ipqf_ctx[2];
+    //DECLARE_ALIGNED(32, float, prev_buf)[2][ATRAC3P_FRAME_SAMPLES]; ///< overlapping buffer
+} Atrac3pChanUnitCtx;
+
+/**
+ * Initialize sine waves synthesizer and ff_sine_* tables.
+ */
+void ff_atrac3p_init_dsp_static(void);
+
+/**
+ * Synthesize sine waves for a particular subband.
+ *
+ * @param[in]   ch_unit   pointer to the channel unit context
+ * @param[in]   fdsp      pointer to float DSP context
+ * @param[in]   ch_num    which channel to process
+ * @param[in]   sb        which subband to process
+ * @param[out]  out       receives processed data
+ */
+void ff_atrac3p_generate_tones(Atrac3pChanUnitCtx *ch_unit,
+                               int ch_num, int sb, float *out);
+
+#ifdef __cplusplus
+}
+#endif /* __cplusplus */
+
+#endif /* AVCODEC_ATRAC3PLUS_H */
diff --git a/src/atrac/at3p/ff/atrac3plusdsp.c b/src/atrac/at3p/ff/atrac3plusdsp.c
new file mode 100644
index 0000000..566f160
--- /dev/null
+++ b/src/atrac/at3p/ff/atrac3plusdsp.c
@@ -0,0 +1,203 @@
+/*
+ * ATRAC3+ compatible decoder
+ *
+ * Copyright (c) 2010-2013 Maxim Poliakovski
+ *
+ * This file is part of FFmpeg.
+ *
+ * FFmpeg is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public
+ * License as published by the Free Software Foundation; either
+ * version 2.1 of the License, or (at your option) any later version.
+ *
+ * FFmpeg 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
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with FFmpeg; if not, write to the Free Software
+ * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
+ */
+
+/**
+ *  @file
+ *  DSP functions for ATRAC3+ decoder.
+ */
+
+#include <math.h>
+#include <string.h>
+
+#include "atrac3plus.h"
+
+#define TWOPI (2 * M_PI)
+#define DEQUANT_PHASE(ph) (((ph) & 0x1F) << 6)
+
+static float sine_table[2048]; ///< wave table
+static float hann_window[256]; ///< Hann windowing function
+static float amp_sf_tab[64];   ///< scalefactors for quantized amplitudes
+
+static void vector_fmul(float *dst, const float *src0, const float *src1,
+                          int len)
+{
+    int i;
+    for (i = 0; i < len; i++)
+        dst[i] = src0[i] * src1[i];
+}
+
+void ff_atrac3p_init_dsp_static(void)
+{
+    int i;
+
+    /* generate sine wave table */
+    for (i = 0; i < 2048; i++)
+        sine_table[i] = sin(TWOPI * i / 2048);
+
+    /* generate Hann window */
+    for (i = 0; i < 256; i++)
+        hann_window[i] = (1.0f - cos(TWOPI * i / 256.0f)) * 0.5f;
+
+    /* generate amplitude scalefactors table */
+    for (i = 0; i < 64; i++)
+        amp_sf_tab[i] = exp2f((i - 3) / 4.0f);
+}
+
+/**
+ *  Synthesize sine waves according to given parameters.
+ *
+ *  @param[in]    synth_param   ptr to common synthesis parameters
+ *  @param[in]    waves_info    parameters for each sine wave
+ *  @param[in]    envelope      envelope data for all waves in a group
+ *  @param[in]    invert_phase  flag indicating 180° phase shift
+ *  @param[in]    reg_offset    region offset for trimming envelope data
+ *  @param[out]   out           receives sythesized data
+ */
+
+static void waves_synth(Atrac3pWaveSynthParams *synth_param,
+                        Atrac3pWavesData *waves_info,
+                        Atrac3pWaveEnvelope *envelope,
+                        int invert_phase, int reg_offset, float *out)
+{
+    int i, wn, inc, pos;
+    double amp;
+    Atrac3pWaveParam *wave_param = &synth_param->waves[waves_info->start_index];
+
+    for (wn = 0; wn < waves_info->num_wavs; wn++, wave_param++) {
+        /* amplitude dequantization */
+        amp = amp_sf_tab[wave_param->amp_sf] *
+              (!synth_param->amplitude_mode
+               ? (wave_param->amp_index + 1) / 15.13f
+               : 1.0f);
+
+        inc = wave_param->freq_index;
+        pos = DEQUANT_PHASE(wave_param->phase_index) - (reg_offset ^ 128) * inc & 2047;
+
+        /* waveform generation */
+        for (i = 0; i < 128; i++) {
+            out[i] += sine_table[pos] * amp;
+            pos     = (pos + inc) & 2047;
+        }
+    }
+
+    /* invert phase if requested */
+    if (invert_phase)
+        for (i = 0; i < 128; i++)
+            out[i] *= -1.0f;
+
+    /* fade in with steep Hann window if requested */
+    if (envelope->has_start_point) {
+        pos = (envelope->start_pos << 2) - reg_offset;
+        if (pos > 0 && pos <= 128) {
+            memset(out, 0, pos * sizeof(*out));
+            if (!envelope->has_stop_point ||
+                envelope->start_pos != envelope->stop_pos) {
+                out[pos + 0] *= hann_window[0];
+                out[pos + 1] *= hann_window[32];
+                out[pos + 2] *= hann_window[64];
+                out[pos + 3] *= hann_window[96];
+            }
+        }
+    }
+
+    /* fade out with steep Hann window if requested */
+    if (envelope->has_stop_point) {
+        pos = (envelope->stop_pos + 1 << 2) - reg_offset;
+        if (pos > 0 && pos <= 128) {
+            out[pos - 4] *= hann_window[96];
+            out[pos - 3] *= hann_window[64];
+            out[pos - 2] *= hann_window[32];
+            out[pos - 1] *= hann_window[0];
+            memset(&out[pos], 0, (128 - pos) * sizeof(out[pos]));
+        }
+    }
+}
+
+void ff_atrac3p_generate_tones(Atrac3pChanUnitCtx *ch_unit,
+                               int ch_num, int sb, float *out)
+{
+    float wavreg1[128] = { 0 };
+    float wavreg2[128] = { 0 };
+    int i, reg1_env_nonzero, reg2_env_nonzero;
+    Atrac3pWavesData *tones_now  = &ch_unit->channels[ch_num].tones_info_prev[sb];
+    Atrac3pWavesData *tones_next = &ch_unit->channels[ch_num].tones_info[sb];
+
+    /* reconstruct full envelopes for both overlapping regions
+     * from truncated bitstream data */
+    if (tones_next->pend_env.has_start_point &&
+        tones_next->pend_env.start_pos < tones_next->pend_env.stop_pos) {
+        tones_next->curr_env.has_start_point = 1;
+        tones_next->curr_env.start_pos       = tones_next->pend_env.start_pos + 32;
+    } else if (tones_now->pend_env.has_start_point) {
+        tones_next->curr_env.has_start_point = 1;
+        tones_next->curr_env.start_pos       = tones_now->pend_env.start_pos;
+    } else {
+        tones_next->curr_env.has_start_point = 0;
+        tones_next->curr_env.start_pos       = 0;
+    }
+
+    if (tones_now->pend_env.has_stop_point &&
+        tones_now->pend_env.stop_pos >= tones_next->curr_env.start_pos) {
+        tones_next->curr_env.has_stop_point = 1;
+        tones_next->curr_env.stop_pos       = tones_now->pend_env.stop_pos;
+    } else if (tones_next->pend_env.has_stop_point) {
+        tones_next->curr_env.has_stop_point = 1;
+        tones_next->curr_env.stop_pos       = tones_next->pend_env.stop_pos + 32;
+    } else {
+        tones_next->curr_env.has_stop_point = 0;
+        tones_next->curr_env.stop_pos       = 64;
+    }
+
+    /* is the visible part of the envelope non-zero? */
+    reg1_env_nonzero = (tones_now->curr_env.stop_pos    < 32) ? 0 : 1;
+    reg2_env_nonzero = (tones_next->curr_env.start_pos >= 32) ? 0 : 1;
+
+    /* synthesize waves for both overlapping regions */
+    if (tones_now->num_wavs && reg1_env_nonzero) {
+        waves_synth(ch_unit->waves_info_prev, tones_now, &tones_now->curr_env,
+                    ch_unit->waves_info_prev->invert_phase[sb] & ch_num,
+                    128, wavreg1);
+    }
+
+    if (tones_next->num_wavs && reg2_env_nonzero) {
+        waves_synth(ch_unit->waves_info, tones_next, &tones_next->curr_env,
+                    ch_unit->waves_info->invert_phase[sb] & ch_num, 0, wavreg2);
+    }
+
+    /* Hann windowing for non-faded wave signals */
+    if (tones_now->num_wavs && tones_next->num_wavs &&
+        reg1_env_nonzero && reg2_env_nonzero) {
+        vector_fmul(wavreg1, wavreg1, &hann_window[128], 128);
+        vector_fmul(wavreg2, wavreg2,  hann_window,      128);
+    } else {
+        if (tones_now->num_wavs && !tones_now->curr_env.has_stop_point)
+            vector_fmul(wavreg1, wavreg1, &hann_window[128], 128);
+
+        if (tones_next->num_wavs && !tones_next->curr_env.has_start_point)
+            vector_fmul(wavreg2, wavreg2, hann_window, 128);
+    }
+
+    /* Overlap and subtract to get residual */
+    for (i = 0; i < 128; i++) {
+        out[i] -= wavreg1[i] + wavreg2[i];
+    }
+}