Ticket #2686: aac-improvements-wip-v2-rclookahead.2.patch

libavcodec/aaccoder.c

@@ -57 +57 @@
     run_value_bits_long, run_value_bits_short
 };
 
-
 /**
  * Quantize one coefficient.
  * @return absolute value of the quantized coefficient
@@ -66 +65 @@
 static av_always_inline int quant(float coef, const float Q)
 {
     float a = coef * Q;
-    return sqrtf(a * sqrtf(a)) + 0.4054;
+    return sqrtf(a * sqrtf(a)) + 0.4054f;
 }
 
 static void quantize_bands(int *out, const float *in, const float *scaled,
@@ -76 +75 @@
     double qc;
     for (i = 0; i < size; i++) {
         qc = scaled[i] * Q34;
-        out[i] = (int)FFMIN(qc + 0.4054, (double)maxval);
+        out[i] = (int)FFMIN(qc + 0.4054f, (double)maxval);
         if (is_signed && in[i] < 0.0f) {
             out[i] = -out[i];
         }
@@ -282 +281 @@
     return maxval;
 }
 
+static float find_max_absval(int group_len, int swb_size, const float *scaled) {
+    float maxval = 0.0f;
+    int w2, i;
+    for (w2 = 0; w2 < group_len; w2++) {
+        for (i = 0; i < swb_size; i++) {
+            maxval = FFMAX(maxval, fabs(scaled[w2*128+i]));
+        }
+    }
+    return maxval;
+}
+
 static int find_min_book(float maxval, int sf) {
     float Q = ff_aac_pow2sf_tab[POW_SF2_ZERO - sf + SCALE_ONE_POS - SCALE_DIV_512];
     float Q34 = sqrtf(Q * sqrtf(Q));
@@ -701 +711 @@
                 sce->sf_idx[(w+w2)*16+g] = sce->sf_idx[w*16+g];
 }
 
+#define sclip(x) av_clip(x,60,218)
+
 /**
  * two-loop quantizers search taken from ISO 13818-7 Appendix C
  */
 static void search_for_quantizers_twoloop(AVCodecContext *avctx,
                                           AACEncContext *s,
                                           SingleChannelElement *sce,
-                                          const float lambda)
+                                          float lambda)
 {
     int start = 0, i, w, w2, g;
-    int destbits = avctx->bit_rate * 1024.0 / avctx->sample_rate / avctx->channels * (lambda / 120.f);
-    float dists[128] = { 0 }, uplims[128];
+    int destbits = avctx->bit_rate * 1024.0 / avctx->sample_rate 
+        / ((avctx->flags & CODEC_FLAG_QSCALE) ? 2.0f : avctx->channels) 
+        * (lambda / 120.f);
+    int toomanybits, toofewbits;
+    float dists[128] = { 0 }, uplims[128], energies[128];
     float maxvals[128];
-    int fflag, minscaler;
+    
+    /* rdlambda controls the maximum tolerated distortion. Twoloop
+     * will keep iterating until it fails to lower it or it reaches
+     * ulimit * rdlambda. Keeping it low increases quality on difficult
+     * signals, but lower it too much, and bits will be taken from weak
+     * signals, creating "holes". A balance is necesary. 
+     */
+    float rdlambda = av_clipf(2 * 120.f / lambda, 0.0625f, 16.0f);
+    float rdmax = (avctx->flags & CODEC_FLAG_QSCALE) ? 4.0f : 1.0f;
+    
+    /* sfoffs controls an offset of optmium allocation that will be
+     * applied based on lambda. Keep it real and modest, the loop
+     * will take care of the rest, this just accelerates convergence
+     */
+    float sfoffs = av_clipf((120.0f - lambda) * 0.1f, -5, 10);
+    
+    int fflag, minscaler, nminscaler, minrdsf;
     int its  = 0;
     int allz = 0;
-    float minthr = INFINITY;
+    int tbits;
+    
+    /* zeroscale controls a multiplier of the threshold, if band energy
+     * is below this, a zero is forced. Keep it lower than 1, unless
+     * low lambda is used, because energy < threshold doesn't mean there's
+     * no audible signal outright, it's just energy. Also make it rise
+     * slower than rdlambda, as rdscale has due compensation with
+     * noisy band depriorization below, whereas zeroing logic is rather dumb
+     */
+    float zeroscale;
+    if (lambda > 120.f)
+        zeroscale = av_clipf(powf(120.f / lambda, 0.25f), 0.0625f, 1.0f);
+    else if (lambda < 80.f)
+        zeroscale = av_clipf(powf(80.f / lambda, 0.25f), 1.0f, 2.0f);
+    else
+        zeroscale = 1.f;
+    
+    if (s->psy.bitres.alloc >= 0) {
+        // Psy granted us extra bits to use, from the reservoire
+        destbits = s->psy.bitres.alloc * (lambda / 120.f);
+    }
+    
+    if (avctx->flags & CODEC_FLAG_QSCALE) {
+        // When using a constant Q-scale, be lenient on bit over/underflow
+        toomanybits = destbits * 3;
+        toofewbits = destbits / 2;
+    } else {
+        // When using ABR, be strict
+        toomanybits = destbits + destbits/16;
+        toofewbits = destbits - destbits/16;
+    }
 
     // for values above this the decoder might end up in an endless loop
     // due to always having more bits than what can be encoded.
     destbits = FFMIN(destbits, 5800);
+    toomanybits = FFMIN(toomanybits, 5800);
+    toofewbits = FFMIN(toofewbits, 5800);
     //XXX: some heuristic to determine initial quantizers will reduce search time
     //determine zero bands and upper limits
     for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) {
         for (g = 0;  g < sce->ics.num_swb; g++) {
             int nz = 0;
-            float uplim = 0.0f;
+            float uplim = INFINITY;
+            float energy = 0.0f;
             for (w2 = 0; w2 < sce->ics.group_len[w]; w2++) {
                 FFPsyBand *band = &s->psy.ch[s->cur_channel].psy_bands[(w+w2)*16+g];
-                uplim += band->threshold;
-                if (band->energy <= band->threshold || band->threshold == 0.0f) {
+                if (band->energy <= band->threshold * zeroscale || band->threshold == 0.0f) {
                     sce->zeroes[(w+w2)*16+g] = 1;
                     continue;
                 }
+                if (uplim > band->threshold) {
+                    uplim = band->threshold;
+                    energy = band->energy;
+                }
                 nz = 1;
             }
-            uplims[w*16+g] = uplim *512;
+            if (!nz)
+                uplim = 0.0f;
+            uplims[w*16+g] = uplim;
+            energies[w*16+g] = energy;
             sce->zeroes[w*16+g] = !nz;
-            if (nz)
-                minthr = FFMIN(minthr, uplim);
             allz |= nz;
         }
     }
+    
+    /* Compute initial scalers */
+    minscaler = 65535;
     for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) {
         for (g = 0;  g < sce->ics.num_swb; g++) {
             if (sce->zeroes[w*16+g]) {
                 sce->sf_idx[w*16+g] = SCALE_ONE_POS;
                 continue;
             }
-            sce->sf_idx[w*16+g] = SCALE_ONE_POS + FFMIN(log2f(uplims[w*16+g]/minthr)*4,59);
+            /* log2f-to-distortion ratio is, technically, 2 (1.5db = 4, but it's power vs level so it's 2). 
+             * But, as offsets are applied, low-frequency signals are too sensitive to the induced distortion,
+             * so we make scaling more conservative by choosing a lower log2f-to-distortion ratio, and thus
+             * more robust.
+             */
+            sce->sf_idx[w*16+g] = SCALE_ONE_POS + 1.75*log2f(FFMAX(0.125f,uplims[w*16+g])) + sfoffs;
+            minscaler = FFMIN(minscaler, sce->sf_idx[w*16+g]);
         }
     }
-
+    
+    /* Clip */
+    minscaler = av_clip(minscaler, SCALE_ONE_POS - SCALE_DIV_512, SCALE_MAX_POS - SCALE_DIV_512);
+    for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) 
+        for (g = 0;  g < sce->ics.num_swb; g++) 
+            if (!sce->zeroes[w*16+g]) 
+                sce->sf_idx[w*16+g] = av_clip(sce->sf_idx[w*16+g], minscaler, minscaler + SCALE_MAX_DIFF - 1);
+    
     if (!allz)
         return;
     abs_pow34_v(s->scoefs, sce->coeffs, 1024);
@@ -766 +850 @@
         }
     }
 
+    /* Scale uplims to match rate distortion to quality 
+     * and apply noisy band depriorization.
+     * Maxval-energy ratio gives us an idea of how noisy the band is. 
+     * If maxval^2 ~ energy, then that band is mostly noise, and we can relax
+     * rate distortion requirements.
+     */
+    for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) {
+        start = w*128;
+        for (g = 0;  g < sce->ics.num_swb; g++) {
+            float max = find_max_absval(sce->ics.group_len[w], sce->ics.swb_sizes[g], sce->coeffs + start);
+            if (max > 0) {
+                float energy2uplim = energies[w*16+g] / (max*max*sce->ics.swb_sizes[g]);
+                energy2uplim = FFMAX(0.0625f, FFMIN(1.0f,energy2uplim));
+                uplims[w*16+g] *= av_clipf(rdlambda * rdlambda * energy2uplim, 0.0625f, rdmax);
+                start += sce->ics.swb_sizes[g];
+            }
+        }
+    }
+
     //perform two-loop search
     //outer loop - improve quality
     do {
-        int tbits, qstep;
-        minscaler = sce->sf_idx[0];
+        int qstep;
         //inner loop - quantize spectrum to fit into given number of bits
         qstep = its ? 1 : 32;
         do {
@@ -790 +892 @@
                         start += sce->ics.swb_sizes[g];
                         continue;
                     }
-                    minscaler = FFMIN(minscaler, sce->sf_idx[w*16+g]);
                     cb = find_min_book(maxvals[w*16+g], sce->sf_idx[w*16+g]);
                     for (w2 = 0; w2 < sce->ics.group_len[w]; w2++) {
                         int b;
@@ -813 +914 @@
                     prev = sce->sf_idx[w*16+g];
                 }
             }
-            if (tbits > destbits) {
+            if (tbits > toomanybits) {
                 for (i = 0; i < 128; i++)
-                    if (sce->sf_idx[i] < 218 - qstep)
-                        sce->sf_idx[i] += qstep;
-            } else {
+                    if (sce->sf_idx[i] < (SCALE_MAX_POS - SCALE_DIV_512))
+                        sce->sf_idx[i] = FFMIN(SCALE_MAX_POS, sce->sf_idx[i] + qstep);
+            } else if (tbits < toofewbits) {
                 for (i = 0; i < 128; i++)
-                    if (sce->sf_idx[i] > 60 - qstep)
-                        sce->sf_idx[i] -= qstep;
+                    if (sce->sf_idx[i] > SCALE_ONE_POS)
+                        sce->sf_idx[i] = FFMAX(SCALE_ONE_POS, sce->sf_idx[i] - qstep);
             }
             qstep >>= 1;
-            if (!qstep && tbits > destbits*1.02 && sce->sf_idx[0] < 217)
+            if (!qstep && tbits > toomanybits && sce->sf_idx[0] < 217)
                 qstep = 1;
         } while (qstep);
 
+        minscaler = SCALE_MAX_POS;
+        for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w])
+            for (g = 0;  g < sce->ics.num_swb; g++)
+                if (!sce->zeroes[w*16+g])
+                    minscaler = FFMIN(minscaler, sce->sf_idx[w*16+g]);
+
         fflag = 0;
-        minscaler = av_clip(minscaler, 60, 255 - SCALE_MAX_DIFF);
+        minscaler = nminscaler = av_clip(minscaler, SCALE_ONE_POS - SCALE_DIV_512, SCALE_MAX_POS - SCALE_DIV_512);
+        minrdsf = (avctx->flags & CODEC_FLAG_QSCALE) ? 60 : minscaler;
         for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) {
             for (g = 0; g < sce->ics.num_swb; g++) {
                 int prevsc = sce->sf_idx[w*16+g];
-                if (dists[w*16+g] > uplims[w*16+g] && sce->sf_idx[w*16+g] > 60) {
+                if (dists[w*16+g] > uplims[w*16+g] && sce->sf_idx[w*16+g] > minrdsf) {
                     if (find_min_book(maxvals[w*16+g], sce->sf_idx[w*16+g]-1))
                         sce->sf_idx[w*16+g]--;
                     else //Try to make sure there is some energy in every band
                         sce->sf_idx[w*16+g]-=2;
                 }
-                sce->sf_idx[w*16+g] = av_clip(sce->sf_idx[w*16+g], minscaler, minscaler + SCALE_MAX_DIFF);
-                sce->sf_idx[w*16+g] = FFMIN(sce->sf_idx[w*16+g], 219);
+                sce->sf_idx[w*16+g] = av_clip(sce->sf_idx[w*16+g], minrdsf, minscaler + SCALE_MAX_DIFF);
+                sce->sf_idx[w*16+g] = FFMIN(sce->sf_idx[w*16+g], SCALE_MAX_POS - SCALE_DIV_512);
                 if (sce->sf_idx[w*16+g] != prevsc)
                     fflag = 1;
+                nminscaler = FFMIN(nminscaler, sce->sf_idx[w*16+g]);
                 sce->band_type[w*16+g] = find_min_book(maxvals[w*16+g], sce->sf_idx[w*16+g]);
             }
         }
+        if (nminscaler < minscaler) {
+            // Drecreased some scalers below minscaler. Must re-clamp.
+            minscaler = nminscaler;
+            for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) 
+                for (g = 0; g < sce->ics.num_swb; g++) 
+                    sce->sf_idx[w*16+g] = av_clip(sce->sf_idx[w*16+g], minscaler, minscaler + SCALE_MAX_DIFF);
+        }
         its++;
-    } while (fflag && its < 10);
+    } while (fflag && its < 20);
+    
+    /* Fill implicit zeroes */
+    for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) 
+        for (g = 0; g < sce->ics.num_swb; g++) 
+            if (sce->band_type[w*16+g] == 0)
+                sce->zeroes[w*16+g] = 1;
 }
 
 static void search_for_quantizers_faac(AVCodecContext *avctx, AACEncContext *s,
@@ -1021 +1143 @@
                                        SingleChannelElement *sce,
                                        const float lambda)
 {
-    int i, w, w2, g;
-    int minq = 255;
+    int w, w2, g;
+    float lowlambda = av_clipf(120.f / lambda, 0.85f, 1.f);
+    float rlambda = av_clipf(120.f / lambda, 0.75f, 10.f);
+    const int minq = av_clip(2 * log2f(120.f / lambda) + 150, 100, 218 - SCALE_MAX_DIFF);
+    const int maxq = minq + SCALE_MAX_DIFF - 1;
 
     memset(sce->sf_idx, 0, sizeof(sce->sf_idx));
     for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) {
         for (g = 0; g < sce->ics.num_swb; g++) {
             for (w2 = 0; w2 < sce->ics.group_len[w]; w2++) {
                 FFPsyBand *band = &s->psy.ch[s->cur_channel].psy_bands[(w+w2)*16+g];
-                if (band->energy <= band->threshold) {
-                    sce->sf_idx[(w+w2)*16+g] = 218;
+                if (band->energy <= 0.05 * lowlambda * band->threshold) {
+                    sce->sf_idx[(w+w2)*16+g] = maxq;
                     sce->zeroes[(w+w2)*16+g] = 1;
                 } else {
-                    sce->sf_idx[(w+w2)*16+g] = av_clip(SCALE_ONE_POS - SCALE_DIV_512 + log2f(band->threshold), 80, 218);
+                    sce->sf_idx[(w+w2)*16+g] = av_clip(SCALE_ONE_POS - SCALE_DIV_512 + 1.414*log2f(band->threshold * rlambda), minq, maxq);
                     sce->zeroes[(w+w2)*16+g] = 0;
                 }
-                minq = FFMIN(minq, sce->sf_idx[(w+w2)*16+g]);
             }
         }
     }
-    for (i = 0; i < 128; i++) {
-        sce->sf_idx[i] = 140;
-        //av_clip(sce->sf_idx[i], minq, minq + SCALE_MAX_DIFF - 1);
-    }
     //set the same quantizers inside window groups
-    for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w])
-        for (g = 0;  g < sce->ics.num_swb; g++)
-            for (w2 = 1; w2 < sce->ics.group_len[w]; w2++)
-                sce->sf_idx[(w+w2)*16+g] = sce->sf_idx[w*16+g];
+    if (sce->ics.num_windows > 1) {
+        for (w = 0; w < sce->ics.num_windows; w += sce->ics.group_len[w]) {
+            for (g = 0;  g < sce->ics.num_swb; g++) {
+                if (sce->ics.group_len[w] > 1) {
+                    int avg_sf_idx = 0;
+                    for (w2 = 0; w2 < sce->ics.group_len[w]; w2++)
+                        avg_sf_idx += sce->sf_idx[w*16+g];
+                    avg_sf_idx /= sce->ics.group_len[w];
+                    for (w2 = 1; w2 < sce->ics.group_len[w]; w2++)
+                        sce->sf_idx[(w+w2)*16+g] = avg_sf_idx;
+                }
+            }
+        }
+    }
+}
+
+static float bval2bmax(float b)
+{
+    /* approximates exp10f(-3.0f*(0.5f + 0.5f * cosf(FFMIN(b,15.5f) / 15.5f))) */
+    return 0.001f + 0.0035f * (b*b*b) / (15.5f*15.5f*15.5f);
 }
 
 static void search_for_ms(AACEncContext *s, ChannelElement *cpe,
@@ -1059 +1195 @@
     float *L34 = s->scoefs, *R34 = s->scoefs + 128, *M34 = s->scoefs + 128*2, *S34 = s->scoefs + 128*3;
     SingleChannelElement *sce0 = &cpe->ch[0];
     SingleChannelElement *sce1 = &cpe->ch[1];
+    
     if (!cpe->common_window)
         return;
+    
     for (w = 0; w < sce0->ics.num_windows; w += sce0->ics.group_len[w]) {
+        int min_sf_idx_mid = SCALE_MAX_POS;
+        int min_sf_idx_side = SCALE_MAX_POS;
+        for (g = 0; g < sce0->ics.num_swb; g++) {
+            if (!sce0->zeroes[w*16+g])
+                min_sf_idx_mid = FFMIN(min_sf_idx_mid, sce0->sf_idx[w*16+g]);
+            if (!sce1->zeroes[w*16+g])
+                min_sf_idx_side = FFMIN(min_sf_idx_side, sce1->sf_idx[w*16+g]);
+        }
+        
         for (g = 0;  g < sce0->ics.num_swb; g++) {
+            float bmax = bval2bmax(g * 17.0f / sce0->ics.num_swb) / 0.0045f;
             if (!cpe->ch[0].zeroes[w*16+g] && !cpe->ch[1].zeroes[w*16+g]) {
                 float dist1 = 0.0f, dist2 = 0.0f;
+                int B0 = 0, B1 = 0;
+                int minidx;
+                int mididx, sididx;
+                float Mmax = 0.0f, Smax = 0.0f;
+                int midcb, sidcb;
+                
+                /* Must compute mid/side SF and book for the whole window group */
+                minidx = FFMIN(sce0->sf_idx[w*16+g], sce1->sf_idx[w*16+g]);
+                for (w2 = 0; w2 < sce0->ics.group_len[w]; w2++) {
+                    for (i = 0; i < sce0->ics.swb_sizes[g]; i++) {
+                        M[i] = (sce0->coeffs[start+w2*128+i]
+                              + sce1->coeffs[start+w2*128+i]) * 0.5;
+                        S[i] =  M[i]
+                              - sce1->coeffs[start+w2*128+i];
+                    }
+                    abs_pow34_v(M34, M, sce0->ics.swb_sizes[g]);
+                    abs_pow34_v(S34, S, sce0->ics.swb_sizes[g]);
+                    for (i = 0; i < sce0->ics.swb_sizes[g]; i++ ) {
+                        Mmax = FFMAX(Mmax, M34[i]);
+                        Smax = FFMAX(Smax, S34[i]);
+                    }
+                }
+                mididx = av_clip(minidx, min_sf_idx_mid, min_sf_idx_mid + SCALE_MAX_DIFF);
+                sididx = av_clip(minidx, min_sf_idx_side, min_sf_idx_side + SCALE_MAX_DIFF);
+                midcb = find_min_book(Mmax, mididx);
+                sidcb = find_min_book(Smax, sididx);
+                
                 for (w2 = 0; w2 < sce0->ics.group_len[w]; w2++) {
                     FFPsyBand *band0 = &s->psy.ch[s->cur_channel+0].psy_bands[(w+w2)*16+g];
                     FFPsyBand *band1 = &s->psy.ch[s->cur_channel+1].psy_bands[(w+w2)*16+g];
                     float minthr = FFMIN(band0->threshold, band1->threshold);
-                    float maxthr = FFMAX(band0->threshold, band1->threshold);
+                    int b1,b2,b3,b4;
                     for (i = 0; i < sce0->ics.swb_sizes[g]; i++) {
                         M[i] = (sce0->coeffs[start+w2*128+i]
                               + sce1->coeffs[start+w2*128+i]) * 0.5;
@@ -1085 +1260 @@
                                                 sce0->ics.swb_sizes[g],
                                                 sce0->sf_idx[(w+w2)*16+g],
                                                 sce0->band_type[(w+w2)*16+g],
-                                                lambda / band0->threshold, INFINITY, NULL);
+                                                lambda / band0->threshold, INFINITY, &b1);
                     dist1 += quantize_band_cost(s, sce1->coeffs + start + w2*128,
                                                 R34,
                                                 sce1->ics.swb_sizes[g],
                                                 sce1->sf_idx[(w+w2)*16+g],
                                                 sce1->band_type[(w+w2)*16+g],
-                                                lambda / band1->threshold, INFINITY, NULL);
+                                                lambda / band1->threshold, INFINITY, &b2);
                     dist2 += quantize_band_cost(s, M,
                                                 M34,
                                                 sce0->ics.swb_sizes[g],
-                                                sce0->sf_idx[(w+w2)*16+g],
-                                                sce0->band_type[(w+w2)*16+g],
-                                                lambda / maxthr, INFINITY, NULL);
+                                                mididx,
+                                                midcb,
+                                                lambda / minthr, INFINITY, &b3);
                     dist2 += quantize_band_cost(s, S,
                                                 S34,
                                                 sce1->ics.swb_sizes[g],
-                                                sce1->sf_idx[(w+w2)*16+g],
-                                                sce1->band_type[(w+w2)*16+g],
-                                                lambda / minthr, INFINITY, NULL);
+                                                sididx,
+                                                sidcb,
+                                                lambda * (lambda / 120.0f) / (minthr * bmax), INFINITY, &b4);
+                    B0 += b1+b2;
+                    B1 += b3+b4;
+                    dist1 -= B0;
+                    dist2 -= B1;
                 }
-                cpe->ms_mask[w*16+g] = dist2 < dist1;
+                cpe->ms_mask[w*16+g] = dist2 <= dist1 && B1 < B0;
+                if (cpe->ms_mask[w*16+g]) {
+                    for (w2 = 0; w2 < sce0->ics.group_len[w]; w2++) {
+                        sce0->sf_idx[(w+w2)*16+g] = mididx;
+                        sce0->band_type[(w+w2)*16+g] = midcb;
+                        sce1->sf_idx[(w+w2)*16+g] = sididx;
+                        sce1->band_type[(w+w2)*16+g] = sidcb;
+                    }
+                }
+            } else {
+                cpe->ms_mask[w*16+g] = 0;
             }
             start += sce0->ics.swb_sizes[g];
         }

libavcodec/aacenc.c

@@ -317 +317 @@
                 if (cpe->common_window && !ch && cpe->ms_mask[w + g]) {
                     for (i = 0; i < ics->swb_sizes[g]; i++) {
                         cpe->ch[0].coeffs[start+i] = (cpe->ch[0].coeffs[start+i] + cpe->ch[1].coeffs[start+i]) / 2.0;
-                        cpe->ch[1].coeffs[start+i] =  cpe->ch[0].coeffs[start+i] - cpe->ch[1].coeffs[start+i];
+                        cpe->ch[1].coeffs[start+i] = cpe->ch[0].coeffs[start+i] - cpe->ch[1].coeffs[start+i];
                     }
                 }
                 start += ics->swb_sizes[g];
@@ -507 +507 @@
     AACEncContext *s = avctx->priv_data;
     float **samples = s->planar_samples, *samples2, *la, *overlap;
     ChannelElement *cpe;
-    int i, ch, w, g, chans, tag, start_ch, ret;
+    int i, ch, w, g, chans, tag, start_ch, ret, frame_bits, its;
     int chan_el_counter[4];
     FFPsyWindowInfo windows[AAC_MAX_CHANNELS];
 
@@ -572 +572 @@
     }
     if ((ret = ff_alloc_packet2(avctx, avpkt, 8192 * s->channels)) < 0)
         return ret;
+    
+    frame_bits = 0;
+    its = 0;
     do {
-        int frame_bits;
-
+        int target_bits, too_many_bits, too_few_bits;
+        
         init_put_bits(&s->pb, avpkt->data, avpkt->size);
 
         if ((avctx->frame_number & 0xFF)==1 && !(avctx->flags & CODEC_FLAG_BITEXACT))
             put_bitstream_info(s, LIBAVCODEC_IDENT);
         start_ch = 0;
+        target_bits = 0;
         memset(chan_el_counter, 0, sizeof(chan_el_counter));
         for (i = 0; i < s->chan_map[0]; i++) {
             FFPsyWindowInfo* wi = windows + start_ch;
@@ -591 +595 @@
             put_bits(&s->pb, 4, chan_el_counter[tag]++);
             for (ch = 0; ch < chans; ch++)
                 coeffs[ch] = cpe->ch[ch].coeffs;
+            s->psy.bitres.alloc = -1;
             s->psy.model->analyze(&s->psy, start_ch, coeffs, wi);
+            target_bits += s->psy.bitres.alloc * (120.f / s->lambda);
             for (ch = 0; ch < chans; ch++) {
                 s->cur_channel = start_ch + ch;
                 s->coder->search_for_quantizers(avctx, s, &cpe->ch[ch], s->lambda);
@@ -635 +641 @@
             start_ch += chans;
         }
 
+        if (avctx->flags & CODEC_FLAG_QSCALE) {
+            // When using a constant Q-scale, don't mess with lambda
+            break;
+        }
+        
+        // rate control stuff
+        // target either the nominal bitrate, or what psy's bit reservoir says to target
+        // whichever is greatest
         frame_bits = put_bits_count(&s->pb);
-        if (frame_bits <= 6144 * s->channels - 3) {
-            s->psy.bitres.bits = frame_bits / s->channels;
+        target_bits = FFMAX(target_bits, avctx->bit_rate * 1024 / avctx->sample_rate);
+        target_bits = FFMIN(target_bits, 6144 * s->channels - 3);
+        
+        // When using ABR, be strict (but only for increasing)
+        too_many_bits = target_bits + target_bits/2;
+        too_few_bits = target_bits - target_bits/8;
+        
+        if (   its == 0 /* for steady-state Q-scale tracking */
+            || (its < 5 && (frame_bits < too_few_bits || frame_bits > too_many_bits))
+            || frame_bits >= 6144 * s->channels - 3  ) 
+        {
+            float prev_lambda = s->lambda;
+            float ratio = ((float)target_bits) / frame_bits;
+            s->lambda = FFMIN(s->lambda * ratio, 65536.f);
+            
+            if (prev_lambda == s->lambda)
+                break;
+            
+            // Keep iterating if we must reduce and lambda is in the sky
+            if (ratio > 0.9f || s->lambda <= 300.f)
+                its++;
+            
+            if (frame_bits >= too_few_bits && frame_bits <= too_many_bits) {
+                /*
+                This path is for steady-state Q-scale tracking
+                When frame bits fall within the stable range, we still need to adjust
+                lambda to maintain it like so in a stable fashion (large jumps in lambda
+                create artifacts and shoulda be avoided)
+                */
+                break;
+            }
+        } else {
             break;
         }
-
-        s->lambda *= avctx->bit_rate * 1024.0f / avctx->sample_rate / frame_bits;
-
     } while (1);
 
     put_bits(&s->pb, 3, TYPE_END);
     flush_put_bits(&s->pb);
     avctx->frame_bits = put_bits_count(&s->pb);
-
-    // rate control stuff
-    if (!(avctx->flags & CODEC_FLAG_QSCALE)) {
-        float ratio = avctx->bit_rate * 1024.0f / avctx->sample_rate / avctx->frame_bits;
-        s->lambda *= ratio;
-        s->lambda = FFMIN(s->lambda, 65536.f);
-    }
+    s->psy.bitres.bits = avctx->frame_bits / s->channels;
 
     if (!frame)
         s->last_frame++;

libavcodec/aacpsy.c

@@ -24 +24 @@
  * AAC encoder psychoacoustic model
  */
 
-#include "libavutil/attributes.h"
 #include "libavutil/libm.h"
 
 #include "avcodec.h"
@@ -255 +254 @@
 /**
  * LAME psy model specific initialization
  */
-static av_cold void lame_window_init(AacPsyContext *ctx, AVCodecContext *avctx)
-{
+static void lame_window_init(AacPsyContext *ctx, AVCodecContext *avctx) {
     int i, j;
 
     for (i = 0; i < avctx->channels; i++) {
@@ -299 +297 @@
     float bark;
     int i, j, g, start;
     float prev, minscale, minath, minsnr, pe_min;
-    const int chan_bitrate = ctx->avctx->bit_rate / ctx->avctx->channels;
+    const int chan_bitrate = ctx->avctx->bit_rate / ((ctx->avctx->flags & CODEC_FLAG_QSCALE) ? 2.0f : ctx->avctx->channels);
     const int bandwidth    = ctx->avctx->cutoff ? ctx->avctx->cutoff : AAC_CUTOFF(ctx->avctx);
     const float num_bark   = calc_bark((float)bandwidth);
 
@@ -391 +389 @@
                                                  int channel, int prev_type)
 {
     int i, j;
-    int br               = ctx->avctx->bit_rate / ctx->avctx->channels;
+    int br               = ((AacPsyContext*)ctx->model_priv_data)->chan_bitrate;
     int attack_ratio     = br <= 16000 ? 18 : 10;
     AacPsyContext *pctx = (AacPsyContext*) ctx->model_priv_data;
     AacPsyChannel *pch  = &pctx->ch[channel];
@@ -681 +679 @@
         desired_pe *= av_clipf(pctx->pe.previous / PSY_3GPP_BITS_TO_PE(ctx->bitres.bits),
                                0.85f, 1.15f);
     pctx->pe.previous = PSY_3GPP_BITS_TO_PE(desired_bits);
+    ctx->bitres.alloc = desired_bits;
 
     if (desired_pe < pe) {
         /* 5.6.1.3.4 "First Estimation of the reduction value" */

libavcodec/psymodel.h

@@ -88 +88 @@
     struct {
         int size;                     ///< size of the bitresevoir in bits
         int bits;                     ///< number of bits used in the bitresevoir
+        int alloc;                    ///< number of bits allocated by the psy, or -1 if no allocation was done
     } bitres;
 
     void* model_priv_data;            ///< psychoacoustic model implementation private data