atrac3 analyze and synthesis MLT and possibility to apply gain info

12 files changed, 519 insertions, 49 deletions

diff --git a/src/atrac/atrac3.cpp b/src/atrac/atrac3.cpp
index 54bf497..cdf6b3f 100644
--- a/src/atrac/atrac3.cpp
+++ b/src/atrac/atrac3.cpp
@@ -15,7 +15,10 @@ constexpr TAtrac3Data::THuffTablePair TAtrac3Data::HuffTables[7];
 
 constexpr TContainerParams TAtrac3Data::ContainerParams[8];
 double TAtrac3Data::EncodeWindow[256] = {0};
+double TAtrac3Data::DecodeWindow[256] = {0};
 double TAtrac3Data::ScaleTable[64] = {0};
+double TAtrac3Data::GainLevel[16];
+double TAtrac3Data::GainInterpolation[31];
 
 const TContainerParams* TAtrac3Data::GetContainerParamsForBitrate(uint32_t bitrate) {
     std::cout << bitrate << std::endl;
diff --git a/src/atrac/atrac3.h b/src/atrac/atrac3.h
index 476d0f0..02e694b 100644
--- a/src/atrac/atrac3.h
+++ b/src/atrac/atrac3.h
@@ -31,8 +31,16 @@ inline bool operator> (const TContainerParams& x, const unsigned int y)
 
 class TAtrac3Data {
 protected:
+    static const uint32_t MDCTSz = 512;
     static double ScaleTable[64];
     static double EncodeWindow[256];
+    static double DecodeWindow[256];
+    static double GainLevel[16];
+    static double GainInterpolation[31];
+    static constexpr int32_t ExponentOffset = 4;
+    static constexpr int32_t LocScale = 3;
+    static constexpr int32_t LocSz = 1 << LocScale;
+    static constexpr int32_t GainInterpolationPosShift = 15;
 
     static const uint32_t NumSamples=1024;
     static const uint32_t frameSz = 152;
@@ -140,6 +148,17 @@ public:
         for (int i = 0; i < 256; i++) {
             EncodeWindow[i] = (sin(((i + 0.5) / 256.0 - 0.5) * M_PI) + 1.0) * 0.5;
         }
+        for (int i = 0; i < 256; i++) {
+            const double a = EncodeWindow[i];
+            const double b = EncodeWindow[255-i];
+            DecodeWindow[i] = 2.0 * a / (a*a + b*b);
+        }
+        for (int i = 0; i < 16; i++) {
+            GainLevel[i] = pow(2.0, ExponentOffset - i);
+        }
+        for (int i = 0; i < 31; i++) {
+            GainInterpolation[i] = pow(2.0, -1.0 / LocSz * (i - 15));
+        }
     }
     static uint32_t MantissaToCLcIdx(int32_t mantissa) {
         assert(mantissa > -3 && mantissa < 2);
@@ -164,27 +183,28 @@ public:
         { 352800, 1024, false }
     };
     static const TContainerParams* GetContainerParamsForBitrate(uint32_t bitrate);
-};
 
-class TAtrac3SubbandInfo {
-public:
-    struct TGainPoint {
-        const uint32_t Level;
-        const uint32_t Location;
+    class SubbandInfo {
+    public:
+        struct TGainPoint {
+            const uint32_t Level;
+            const uint32_t Location;
+        };
+    private:
+        std::vector<std::vector<TGainPoint>> Info;
+    public:
+        SubbandInfo()
+        {
+            Info.resize(4);
+        }
+        void AddSubbandCurve(uint16_t n, std::vector<TGainPoint>&& curve) {
+            Info[n] = std::move(curve);
+        }
+        uint32_t GetQmfNum() const {
+            return Info.size();
+        }
+        const std::vector<TGainPoint>& GetGainPoints(uint32_t i) const {
+            return Info[i];
+        }
     };
-private:
-    std::vector<std::vector<TGainPoint>> Info;
-public:
-    TAtrac3SubbandInfo()
-    {
-        Info.resize(4);
-    }
-    uint32_t GetQmfNum() const {
-        return Info.size();
-    }
-    const std::vector<TGainPoint>& GetGainPoints(uint32_t i) const {
-        assert(Info[i].size() == 0);
-        return Info[i];
-    }
 };
-
diff --git a/src/atrac/atrac3_bitstream.cpp b/src/atrac/atrac3_bitstream.cpp
index 395bb1c..ff28132 100644
--- a/src/atrac/atrac3_bitstream.cpp
+++ b/src/atrac/atrac3_bitstream.cpp
@@ -122,7 +122,7 @@ void TAtrac3BitStreamWriter::EncodeSpecs(const vector<TScaledBlock>& scaledBlock
     }
 }
 
-void TAtrac3BitStreamWriter::WriteSoundUnit(const TAtrac3SubbandInfo& subbandInfo, const vector<TScaledBlock>& scaledBlocks) {
+void TAtrac3BitStreamWriter::WriteSoundUnit(const TAtrac3Data::SubbandInfo& subbandInfo, const vector<TScaledBlock>& scaledBlocks) {
     NBitStream::TBitStream bitStream;
     if (Params.Js) {
         //TODO
@@ -134,11 +134,9 @@ void TAtrac3BitStreamWriter::WriteSoundUnit(const TAtrac3SubbandInfo& subbandInf
 
     //write gain info
     for (uint32_t band = 0; band < numQmfBand; ++band) {
-        const vector<TAtrac3SubbandInfo::TGainPoint>& GainPoints = subbandInfo.GetGainPoints(band);
+        const vector<TAtrac3Data::SubbandInfo::TGainPoint>& GainPoints = subbandInfo.GetGainPoints(band);
         bitStream.Write(GainPoints.size(), 3);
-        assert(GainPoints.size() == 0);
-        for (const TAtrac3SubbandInfo::TGainPoint& point : GainPoints) {
-            abort();
+        for (const TAtrac3Data::SubbandInfo::TGainPoint& point : GainPoints) {
             bitStream.Write(point.Level, 4);
             bitStream.Write(point.Location, 5);
         }
diff --git a/src/atrac/atrac3_bitstream.h b/src/atrac/atrac3_bitstream.h
index af9763c..60f1e4a 100644
--- a/src/atrac/atrac3_bitstream.h
+++ b/src/atrac/atrac3_bitstream.h
@@ -25,5 +25,5 @@ public:
     {
 
     }
-    void WriteSoundUnit(const TAtrac3SubbandInfo& subbandInfo, const std::vector<NAtracDEnc::TScaledBlock>& scaledBlocks);
+    void WriteSoundUnit(const TAtrac3Data::SubbandInfo& subbandInfo, const std::vector<NAtracDEnc::TScaledBlock>& scaledBlocks);
 };
diff --git a/src/atrac3denc.cpp b/src/atrac3denc.cpp
index 9b5b9e6..7557d5e 100644
--- a/src/atrac3denc.cpp
+++ b/src/atrac3denc.cpp
@@ -1,5 +1,6 @@
 #include "atrac3denc.h"
 #include "atrac/atrac3_bitstream.h"
+#include "util.h"
 #include <assert.h>
 
 #include <iostream>
@@ -9,11 +10,16 @@ namespace NAtracDEnc {
 using namespace NMDCT;
 using std::vector;
 
-void TAtrac3MDCT::Mdct(double specs[1024], double* bands[4]) {
+void TAtrac3MDCT::Mdct(double specs[1024], double* bands[4], TGainModulatorArray gainModulators) {
     for (int band = 0; band < 4; ++band) {
         double* srcBuff = bands[band];
+        double* const curSpec = &specs[band*256];
+        TGainModulator modFn = gainModulators[band];
         vector<double> tmp(512);
         memcpy(&tmp[0], &srcBuff[256], 256 * sizeof(double));
+        if (modFn) {
+            modFn(tmp.data(), srcBuff);
+        }
         for (int i = 0; i < 256; i++) {
             srcBuff[256+i] = TAtrac3Data::EncodeWindow[i] * srcBuff[i];
             srcBuff[i] = TAtrac3Data::EncodeWindow[255-i] * srcBuff[i];
@@ -21,14 +27,36 @@ void TAtrac3MDCT::Mdct(double specs[1024], double* bands[4]) {
         memcpy(&tmp[256], &srcBuff[0], 256 * sizeof(double));
         const vector<double>& sp = Mdct512(&tmp[0]);
         assert(sp.size() == 256);
-        memcpy(&specs[band*256], sp.data(), 256*sizeof(double));
+        memcpy(curSpec, sp.data(), 256 * sizeof(double));
+        if (band & 1) {
+            SwapArray(curSpec, 256);
+        }
+    }
+}
+
+void TAtrac3MDCT::Midct(double specs[1024], double* bands[4], TGainDemodulatorArray gainDemodulators) {
+    for (int band = 0; band < 4; ++band) {
+        double* dstBuff = bands[band];
+        double* curSpec = &specs[band*256];
+        double* prevBuff = dstBuff + 256;
+        TAtrac3GainProcessor::TGainDemodulator demodFn = gainDemodulators[band];
         if (band & 1) {
-            for (uint32_t j = 0; j < sp.size() / 2; j++) {
-                double tmp = specs[band*256 +j];
-                specs[band*256 + j] = specs[band*256 + sp.size() - 1 -j];
-                specs[band*256 + sp.size() - 1 -j] = tmp;
+            SwapArray(curSpec, 256);
+        }
+        vector<double> inv  = Midct512(curSpec);
+        assert(inv.size()/2 == 256);
+        for (int j = 0; j < 256; ++j) {
+            inv[j] *= 2 * DecodeWindow[j];
+            inv[511 - j] *= 2 * DecodeWindow[j];
+        }
+        if (demodFn) {
+            demodFn(dstBuff, inv.data(), prevBuff);
+        } else {
+            for (uint32_t j = 0; j < 256; ++j) {
+                dstBuff[j] = inv[j] + prevBuff[j];
             }
         }
+        memcpy(prevBuff, &inv[256], sizeof(double)*256);
     }
 }
 
@@ -40,6 +68,35 @@ TAtrac3Processor::TAtrac3Processor(TAeaPtr&& oma, const TContainerParams& params
 TAtrac3Processor::~TAtrac3Processor()
 {}
 
+TAtrac3MDCT::TGainModulatorArray TAtrac3MDCT::MakeGainModulatorArray(const TAtrac3Data::SubbandInfo& si) {
+    switch (si.GetQmfNum()) {
+        case 1:
+        {
+            return {{ GainProcessor.Modulate(si.GetGainPoints(0)), TAtrac3MDCT::TGainModulator(),
+                TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator() }};
+        }
+        case 2:
+        {
+            return {{ GainProcessor.Modulate(si.GetGainPoints(0)), GainProcessor.Modulate(si.GetGainPoints(1)),
+                TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator() }};
+        }
+        case 3:
+        {
+            return {{ GainProcessor.Modulate(si.GetGainPoints(0)), GainProcessor.Modulate(si.GetGainPoints(1)),
+                GainProcessor.Modulate(si.GetGainPoints(2)), TAtrac3MDCT::TGainModulator() }};
+        }
+        case 4:
+        {
+            return {{ GainProcessor.Modulate(si.GetGainPoints(0)), GainProcessor.Modulate(si.GetGainPoints(1)),
+                GainProcessor.Modulate(si.GetGainPoints(2)), GainProcessor.Modulate(si.GetGainPoints(3)) }};
+        }
+        default:
+            assert(false);
+            return {};
+
+    }
+}
+
 TPCMEngine<double>::TProcessLambda TAtrac3Processor::GetEncodeLambda() {
     TOma* omaptr = dynamic_cast<TOma*>(Oma.get());
     if (!omaptr) {
@@ -53,14 +110,17 @@ TPCMEngine<double>::TProcessLambda TAtrac3Processor::GetEncodeLambda() {
             vector<double> specs(1024);
             double src[NumSamples];
             for (int i = 0; i < NumSamples; ++i) {
-                src[i] = data[meta.Channels == 1 ? i : (i * 2 + channel)]; //no mono mode in atrac1. //TODO we can double frame after encoding
+                src[i] = data[meta.Channels == 1 ? i : (i * 2 + channel)]; //no mono mode in atrac3. //TODO we can double frame after encoding
             }
 
             double* p[4] = {&PcmBuffer[channel][0][0], &PcmBuffer[channel][1][0], &PcmBuffer[channel][2][0], &PcmBuffer[channel][3][0]};
             SplitFilterBank[channel].Split(&src[0], p);
-            Mdct(specs.data(), p);
+
+            TAtrac3Data::SubbandInfo siCur;
+
+            Mdct(specs.data(), p, MakeGainModulatorArray(siCur));
             const TBlockSize blockSize(false, false, false);
-            bitStreamWriter->WriteSoundUnit(TAtrac3SubbandInfo(), Scaler.Scale(specs, blockSize)); 
+            bitStreamWriter->WriteSoundUnit(siCur, Scaler.Scale(specs, blockSize)); 
         }
     };
 }
diff --git a/src/atrac3denc.h b/src/atrac3denc.h
index 0b20061..149f1c6 100644
--- a/src/atrac3denc.h
+++ b/src/atrac3denc.h
@@ -7,13 +7,27 @@
 
 #include "atrac/atrac_scale.h"
 #include "mdct/mdct.h"
+#include "gain_processor.h"
 
+#include <functional>
+#include <array>
 namespace NAtracDEnc {
 
 class TAtrac3MDCT : public virtual TAtrac3Data {
     NMDCT::TMDCT<512> Mdct512;
+    NMDCT::TMIDCT<512> Midct512;
 public:
-    void Mdct(double specs[1024], double* bands[4]);
+    typedef TGainProcessor<TAtrac3Data> TAtrac3GainProcessor;
+    TAtrac3GainProcessor GainProcessor;
+public:
+    using TGainModulator = TAtrac3GainProcessor::TGainModulator;
+    using TGainDemodulator = TAtrac3GainProcessor::TGainDemodulator;
+    typedef std::array<TGainDemodulator, 4> TGainDemodulatorArray;
+    typedef std::array<TGainModulator, 4> TGainModulatorArray;
+    void Mdct(double specs[1024], double* bands[4], TGainModulatorArray gainModulators = TGainModulatorArray());
+    void Midct(double specs[1024], double* bands[4], TGainDemodulatorArray gainDemodulators = TGainDemodulatorArray());
+protected:
+    TAtrac3MDCT::TGainModulatorArray MakeGainModulatorArray(const TAtrac3Data::SubbandInfo& si);
 };
 
 class TAtrac3Processor : public IProcessor<double>, public TAtrac3MDCT, public virtual TAtrac3Data {
diff --git a/src/atrac3denc_ut.cpp b/src/atrac3denc_ut.cpp
new file mode 100644
index 0000000..f15c799
--- /dev/null
+++ b/src/atrac3denc_ut.cpp
@@ -0,0 +1,253 @@
+#include "atrac3denc.h"
+#include <gtest/gtest.h>
+
+#include <vector>
+#include <cmath>
+using std::vector;
+using namespace NAtracDEnc;
+
+static void GenerateSignal(double* buf, size_t n, double f, double a) {
+    for (size_t i = 0; i < n; ++i) {
+        buf[i] = a * sin((M_PI/2) * i * f);
+    }
+}
+
+class TWindowTest : public TAtrac3Data {
+public:
+    void RunTest() {
+        for (size_t i = 0; i < 256; i++) {
+            const double ha1 = EncodeWindow[i] / 2.0; //compensation
+            const double hs1 = DecodeWindow[i];
+            const double hs2 = DecodeWindow[255-i];
+            const double res = hs1 / (hs1 * hs1 + hs2 * hs2);
+            EXPECT_NEAR(ha1, res, 0.000000001);
+        }
+    }
+};
+
+template<class T>
+class TAtrac3MDCTWorkBuff {
+    T* Buffer;
+public:
+    static const size_t BandBuffSz = 256;
+    static const size_t BandBuffAndOverlapSz = BandBuffSz * 2;
+    static const size_t BuffSz = BandBuffAndOverlapSz * (4 + 4); 
+    T* const Band0;
+    T* const Band1;
+    T* const Band2;
+    T* const Band3;
+    T* const Band0Res;
+    T* const Band1Res;
+    T* const Band2Res;
+    T* const Band3Res;
+    TAtrac3MDCTWorkBuff()
+        : Buffer(new T[BuffSz])
+        , Band0(Buffer)
+        , Band1(Band0 + BandBuffAndOverlapSz)
+        , Band2(Band1 + BandBuffAndOverlapSz)
+        , Band3(Band2 + BandBuffAndOverlapSz)
+        , Band0Res(Band3 + BandBuffAndOverlapSz)
+        , Band1Res(Band0Res + BandBuffAndOverlapSz)
+        , Band2Res(Band1Res + BandBuffAndOverlapSz)
+        , Band3Res(Band2Res + BandBuffAndOverlapSz)
+    {
+        memset(Buffer, 0, sizeof(T)*BuffSz);
+    }
+    ~TAtrac3MDCTWorkBuff()
+    {
+        delete[] Buffer;
+    }
+};
+
+
+TEST(TAtrac3MDCT, TAtrac3MDCTZeroOneBlock) {
+    TAtrac3MDCT mdct;
+    TAtrac3MDCTWorkBuff<double> buff;
+    size_t workSz = TAtrac3MDCTWorkBuff<double>::BandBuffSz;
+
+    vector<double> specs(1024);
+
+    double* p[4] = { buff.Band0, buff.Band1, buff.Band2, buff.Band3 };
+
+    mdct.Mdct(specs.data(), p);
+    for(auto s: specs)
+        EXPECT_NEAR(s, 0.0, 0.0000000001);
+
+    double* t[4] = { buff.Band0Res, buff.Band1Res, buff.Band2Res, buff.Band3Res };
+    mdct.Midct(specs.data(), p);
+
+    for(size_t i = 0; i < workSz; ++i)
+        EXPECT_NEAR(buff.Band0Res[i], 0.0, 0.0000000001);
+
+   for(size_t i = 0; i < workSz; ++i)
+        EXPECT_NEAR(buff.Band1Res[i], 0.0, 0.0000000001);
+
+   for(size_t i = 0; i < workSz; ++i)
+        EXPECT_NEAR(buff.Band2Res[i], 0.0, 0.0000000001);
+
+   for(size_t i = 0; i < workSz; ++i)
+        EXPECT_NEAR(buff.Band3Res[i], 0.0, 0.0000000001);
+
+
+}
+
+TEST(TAtrac3MDCT, TAtrac3MDCTSignal) {
+    TAtrac3MDCT mdct;
+    TAtrac3MDCTWorkBuff<double> buff;
+    size_t workSz = TAtrac3MDCTWorkBuff<double>::BandBuffSz;
+
+    const size_t len = 1024;
+    vector<double> signal(len);
+    vector<double> signalRes(len);
+    GenerateSignal(signal.data(), signal.size(), 0.25, 32768);
+    
+    for (size_t pos = 0; pos < len; pos += workSz) {
+        vector<double> specs(1024);
+        memcpy(buff.Band0, signal.data() + pos, workSz * sizeof(double));
+
+        double* p[4] = { buff.Band0, buff.Band1, buff.Band2, buff.Band3 };
+        mdct.Mdct(specs.data(), p);
+
+        double* t[4] = { buff.Band0Res, buff.Band1Res, buff.Band2Res, buff.Band3Res };
+        mdct.Midct(specs.data(), t);
+
+        memcpy(signalRes.data() + pos, buff.Band0Res, workSz * sizeof(double));
+    }
+
+    for (int i = workSz; i < len; ++i)
+        EXPECT_NEAR(signal[i - workSz], signalRes[i], 0.00000001);
+}
+
+TEST(TAtrac3MDCT, TAtrac3MDCTSignalWithGainCompensation) {
+    TAtrac3MDCT mdct;
+    TAtrac3MDCTWorkBuff<double> buff;
+    size_t workSz = TAtrac3MDCTWorkBuff<double>::BandBuffSz;
+
+    const size_t len = 4096;
+    vector<double> signal(len, 8000);
+    vector<double> signalRes(len);
+    GenerateSignal(signal.data() + 1024, signal.size()-1024, 0.25, 32768);
+    
+    for (size_t pos = 0; pos < len; pos += workSz) {
+        vector<double> specs(1024);
+        memcpy(buff.Band0, signal.data() + pos, workSz * sizeof(double));
+
+        double* p[4] = { buff.Band0, buff.Band1, buff.Band2, buff.Band3 };
+
+        if (pos == 256) { //apply gain modulation
+            TAtrac3Data::SubbandInfo siCur;
+            siCur.AddSubbandCurve(0, {{3, 2}});
+
+            mdct.Mdct(specs.data(), p, { mdct.GainProcessor.Modulate(siCur.GetGainPoints(0)),
+                    TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator()});
+        } else if (pos == 1024) {
+            TAtrac3Data::SubbandInfo siCur;
+            siCur.AddSubbandCurve(0, {{3, 2}});
+            siCur.AddSubbandCurve(0, {{2, 5}});
+
+            mdct.Mdct(specs.data(), p, { mdct.GainProcessor.Modulate(siCur.GetGainPoints(0)),
+                    TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator()});
+        } else if (pos == 1024 + 256) {
+            TAtrac3Data::SubbandInfo siCur;
+            siCur.AddSubbandCurve(0, {{1, 0}});
+
+            mdct.Mdct(specs.data(), p, { mdct.GainProcessor.Modulate(siCur.GetGainPoints(0)),
+                    TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator()});
+        } else if (pos == 2048) {
+            TAtrac3Data::SubbandInfo siCur;
+            siCur.AddSubbandCurve(0, {{4, 2}});
+            siCur.AddSubbandCurve(0, {{1, 5}});
+
+            mdct.Mdct(specs.data(), p, { mdct.GainProcessor.Modulate(siCur.GetGainPoints(0)),
+                    TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator(), TAtrac3MDCT::TGainModulator()});
+        } else {
+            mdct.Mdct(specs.data(), p);
+        }
+
+        double* t[4] = { buff.Band0Res, buff.Band1Res, buff.Band2Res, buff.Band3Res };
+
+        if (pos == 256) { //restore gain modulation
+            TAtrac3Data::SubbandInfo siCur;
+            TAtrac3Data::SubbandInfo siNext;
+            siNext.AddSubbandCurve(0, {{3, 2}});
+ 
+            mdct.Midct(specs.data(), t, {mdct.GainProcessor.Demodulate(siCur.GetGainPoints(0), siNext.GetGainPoints(0)),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator()});
+        } else if (pos == 512) {
+            TAtrac3Data::SubbandInfo siNext;
+            TAtrac3Data::SubbandInfo siCur;
+            siCur.AddSubbandCurve(0, {{3, 2}});
+ 
+            mdct.Midct(specs.data(), t, {mdct.GainProcessor.Demodulate(siCur.GetGainPoints(0), siNext.GetGainPoints(0)),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator()});
+        } else if (pos == 1024) {
+            TAtrac3Data::SubbandInfo siCur;
+            TAtrac3Data::SubbandInfo siNext;
+            siNext.AddSubbandCurve(0, {{3, 2}});
+            siNext.AddSubbandCurve(0, {{2, 5}});
+
+            mdct.Midct(specs.data(), t, {mdct.GainProcessor.Demodulate(siCur.GetGainPoints(0), siNext.GetGainPoints(0)),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator()});
+        } else if (pos == 1024 + 256) {
+            TAtrac3Data::SubbandInfo siNext;
+            TAtrac3Data::SubbandInfo siCur;
+            siNext.AddSubbandCurve(0, {{1, 0}});
+            siCur.AddSubbandCurve(0, {{3, 2}});
+            siCur.AddSubbandCurve(0, {{2, 5}});
+
+            mdct.Midct(specs.data(), t, {mdct.GainProcessor.Demodulate(siCur.GetGainPoints(0), siNext.GetGainPoints(0)),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator()});
+        } else if (pos == 1024 + 256 + 256) {
+            TAtrac3Data::SubbandInfo siNext;
+            TAtrac3Data::SubbandInfo siCur;
+            siCur.AddSubbandCurve(0, {{1, 0}});
+
+            mdct.Midct(specs.data(), t, {mdct.GainProcessor.Demodulate(siCur.GetGainPoints(0), siNext.GetGainPoints(0)),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator()});
+        } else if (pos == 2048) {
+            TAtrac3Data::SubbandInfo siCur;
+            TAtrac3Data::SubbandInfo siNext;
+            siNext.AddSubbandCurve(0, {{4, 2}});
+            siNext.AddSubbandCurve(0, {{1, 5}});
+
+            mdct.Midct(specs.data(), t, {mdct.GainProcessor.Demodulate(siCur.GetGainPoints(0), siNext.GetGainPoints(0)),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator()});
+        } else if (pos == 2048 + 256) {
+            TAtrac3Data::SubbandInfo siNext;
+            TAtrac3Data::SubbandInfo siCur;
+            siCur.AddSubbandCurve(0, {{4, 2}});
+            siCur.AddSubbandCurve(0, {{1, 5}});
+
+            mdct.Midct(specs.data(), t, {mdct.GainProcessor.Demodulate(siCur.GetGainPoints(0), siNext.GetGainPoints(0)),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator(),
+                    TAtrac3MDCT::TAtrac3GainProcessor::TGainDemodulator()});
+        } else {
+            mdct.Midct(specs.data(), t);
+        }
+        memcpy(signalRes.data() + pos, buff.Band0Res,  workSz * sizeof(double));
+    }
+    for (int i = workSz; i < len; ++i) {
+        //std::cout << "res: " << i << " " << signalRes[i] << std::endl;
+        EXPECT_NEAR(signal[i - workSz], signalRes[i], 0.00000001);
+    }
+}
+
+
+
+TEST(TAtrac3MDCT, TAtrac3MDCTWindow) {
+    TWindowTest test;
+    test.RunTest();
+}
diff --git a/src/atracdenc.cpp b/src/atracdenc.cpp
index 94cf394..4754256 100644
--- a/src/atracdenc.cpp
+++ b/src/atracdenc.cpp
@@ -6,6 +6,7 @@
 #include "atrac/atrac1_dequantiser.h"
 #include "atrac/atrac1_qmf.h"
 #include "atrac/atrac1_bitalloc.h"
+#include "util.h"
 
 namespace NAtracDEnc {
 using namespace std;
@@ -85,11 +86,7 @@ void TAtrac1MDCT::Mdct(double Specs[512], double* low, double* mid, double* hi,
                 Specs[blockPos + pos + i] = sp[i] * multiple;
             }
             if (band) {
-                for (uint32_t j = 0; j < sp.size() / 2; j++) {
-                    double tmp = Specs[blockPos + pos +j];
-                    Specs[blockPos + pos + j] = Specs[blockPos + pos + sp.size() - 1 -j];
-                    Specs[blockPos + pos + sp.size() - 1 -j] = tmp;
-                }
+                SwapArray(&Specs[blockPos + pos], sp.size());
             }
 
             blockPos += 32;
@@ -110,15 +107,9 @@ void TAtrac1MDCT::IMdct(double Specs[512], const TBlockSize& mode, double* low,
         vector<double> invBuf(512);
         double* prevBuf = &dstBuf[bufSz * 2  - 16];
         for (uint32_t block = 0; block < numMdctBlocks; block++) {
-
             if (band) {
-                for (uint32_t j = 0; j < blockSz/2; j++) {
-                    double tmp = Specs[pos+j];
-                    Specs[pos+j] = Specs[pos + blockSz - 1 -j];
-                    Specs[pos + blockSz - 1 -j] = tmp;
-                }
+                SwapArray(&Specs[pos], blockSz);
             }
-
             vector<double> inv = (numMdctBlocks != 1) ? midct(&Specs[pos], blockSz) : (bufSz == 128) ? Midct256(&Specs[pos]) : Midct512(&Specs[pos]);
             for (int i = 0; i < (inv.size()/2); i++) {
                 invBuf[start+i] = inv[i + inv.size()/4];
diff --git a/src/gain_processor.h b/src/gain_processor.h
new file mode 100644
index 0000000..8369ccf
--- /dev/null
+++ b/src/gain_processor.h
@@ -0,0 +1,82 @@
+#include <functional>
+
+template<class T>
+class TGainProcessor : public virtual T {
+
+public:
+    typedef std::function<void(double* out, double* cur, double* prev)> TGainDemodulator;
+    /*
+     * example GainModulation:
+     * PCMinput:
+     *     N   b    N        N
+     * --------|--------|--------|--------
+     * |       | - mdct #1
+     *     |        | - mdct #2
+     *     a
+     *         |        | - mdct #3
+     *         ^^^^^ - modulated by previous step
+     * lets consider a case we want to modulate mdct #2.
+     *     bufCur - is a buffer of first half of mdct transformation (a)
+     *     bufNext - is a buffer of second half of mdct transformation and overlaping (i.e the input buffer started at b point)
+     * so next transformation (mdct #3) gets modulated first part
+     */
+    typedef std::function<void(double* bufCur, double* bufNext)> TGainModulator;
+    TGainDemodulator Demodulate(const std::vector<typename T::SubbandInfo::TGainPoint>& giNow, const std::vector<typename T::SubbandInfo::TGainPoint>& giNext) {
+        return [=](double* out, double* cur, double* prev) {
+            uint32_t pos = 0;
+            const double scale = giNext.size() ? T::GainLevel[giNext[0].Level] : 1;
+            for (uint32_t i = 0; i < giNow.size(); ++i) {
+                uint32_t lastPos = giNow[i].Location << T::LocScale;
+                const uint32_t levelPos = giNow[i].Level;
+                assert(levelPos < sizeof(T::GainLevel)/sizeof(T::GainLevel[0]));
+                double level = T::GainLevel[levelPos];
+                const int incPos = ((i + 1) < giNow.size() ? giNow[i + 1].Level : T::ExponentOffset) - giNow[i].Level + T::GainInterpolationPosShift;
+                double gainInc = T::GainInterpolation[incPos];
+                for (; pos < lastPos; pos++) {
+                    //std::cout << "pos: " << pos << " scale: " << scale << " level: " << level << std::endl;
+                    out[pos] = (cur[pos] * scale + prev[pos]) * level;
+                }
+                for (; pos < lastPos + T::LocSz; pos++) {
+                    //std::cout << "pos: " << pos << " scale: " << scale << " level: " << level << " gainInc: " << gainInc << std::endl;
+                    out[pos] = (cur[pos] * scale + prev[pos]) * level;
+                    level *= gainInc;
+                }
+            }
+            for (; pos < T::MDCTSz/2; pos++) {
+                //std::cout << "pos: " << pos << " scale: " << scale << std::endl;
+                out[pos] = cur[pos] * scale + prev[pos];
+            }
+        };
+    }
+    TGainModulator Modulate(const std::vector<typename T::SubbandInfo::TGainPoint>& giCur) {
+        if (giCur.empty())
+            return {};
+        return [=](double* bufCur, double* bufNext) {
+            uint32_t pos = 0;
+            const double scale = T::GainLevel[giCur[0].Level];
+            for (uint32_t i = 0; i < giCur.size(); ++i) {
+                uint32_t lastPos = giCur[i].Location << T::LocScale;
+                const uint32_t levelPos = giCur[i].Level;
+                assert(levelPos < sizeof(T::GainLevel)/sizeof(T::GainLevel[0]));
+                double level = T::GainLevel[levelPos];
+                const int incPos = ((i + 1) < giCur.size() ? giCur[i + 1].Level : T::ExponentOffset) - giCur[i].Level + T::GainInterpolationPosShift;
+                double gainInc = T::GainInterpolation[incPos];
+                for (; pos < lastPos; pos++) {
+                    bufCur[pos] /= scale;
+                    bufNext[pos] /= level;
+                    //std::cout << "mod pos: " << pos << " scale: " << scale << " level: " << level << std::endl;
+                }
+                for (; pos < lastPos + T::LocSz; pos++) {
+                    bufCur[pos] /= scale;
+                    bufNext[pos] /= level;
+                    //std::cout << "mod pos: " << pos << " scale: " << scale << " level: " << level << " gainInc: " << gainInc << std::endl;
+                    level *= gainInc;
+                }
+            }
+            for (; pos < T::MDCTSz/2; pos++) {
+                bufCur[pos] /= scale;
+                //std::cout << "mod pos: " << pos << " scale: " << scale << std::endl;
+            }
+        };
+    }
+};
diff --git a/src/util.h b/src/util.h
new file mode 100644
index 0000000..e167655
--- /dev/null
+++ b/src/util.h
@@ -0,0 +1,12 @@
+#pragma once
+#include <cstdint>
+
+
+template<class T>
+inline void SwapArray(T* p, const size_t len) {
+    for (size_t i = 0, j = len - 1; i < len / 2; ++i, --j) {
+        T tmp = p[i];
+        p[i] = p[j];
+        p[j] = tmp;
+    }
+}
diff --git a/src/util_ut.cpp b/src/util_ut.cpp
new file mode 100644
index 0000000..8c5ea29
--- /dev/null
+++ b/src/util_ut.cpp
@@ -0,0 +1,15 @@
+#include "util.h"
+#include <gtest/gtest.h>
+
+#include <vector>
+
+
+TEST(Util, SwapArrayTest) {
+
+    double arr[8] = { 0, 1, 2, 3, 4, 5, 6, 7 };
+    SwapArray(arr, 8);
+    for (size_t i = 0; i < 8; ++i) {
+        EXPECT_NEAR((double)i, arr[7-i], 0.000000000001);
+    }
+
+}
diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt
index 6e8517b..82adb57 100644
--- a/test/CMakeLists.txt
+++ b/test/CMakeLists.txt
@@ -1,3 +1,5 @@
+set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -g -fsanitize=address -fno-omit-frame-pointer")
+
 include_directories(${gtest_SOURCE_DIR}/include)
 
 set(mdct_test_sources
@@ -28,3 +30,23 @@ set(atrac1mdct_test_sources
 )
 add_executable(atrac1mdct_test ${atrac1mdct_test_sources})
 target_link_libraries(atrac1mdct_test mdct_impl gtest_main)
+
+set(atrac3mdct_test_sources
+    ../src/atrac3denc.cpp
+    ../src/transient_detector.cpp
+    ../src/bitstream/bitstream.cpp
+    ../src/atrac/atrac_scale.cpp
+    ../src/atrac/atrac3_bitstream.cpp
+    ../src/atrac/atrac1.cpp #atrac_scale has explicit instantiation
+    ../src/atrac/atrac3.cpp
+    ../src/atrac3denc_ut.cpp
+    ../src/oma.cpp
+)
+add_executable(atrac3mdct_test ${atrac3mdct_test_sources})
+target_link_libraries(atrac3mdct_test mdct_impl oma gtest_main)
+
+set(util_test_sources
+    ../src/util_ut.cpp
+)
+add_executable(util_test ${util_test_sources})
+target_link_libraries(util_test gtest_main)