1
2 #ifdef HAVE_CONFIG_H
3 #include "config.h"
4 #endif
5 #include <schroedinger/schro.h>
6 #include <schroedinger/schrohistogram.h>
7 #include <string.h>
8 #include <math.h>
9
10 //#define DUMP_SUBBAND_CURVES
11
12 void schro_encoder_choose_quantisers_simple (SchroEncoderFrame * frame);
13 void schro_encoder_choose_quantisers_rdo_bit_allocation (SchroEncoderFrame *
14 frame);
15 void schro_encoder_choose_quantisers_lossless (SchroEncoderFrame * frame);
16 void schro_encoder_choose_quantisers_lowdelay (SchroEncoderFrame * frame);
17 void schro_encoder_choose_quantisers_rdo_lambda (SchroEncoderFrame * frame);
18 void schro_encoder_choose_quantisers_rdo_cbr (SchroEncoderFrame * frame);
19 void schro_encoder_choose_quantisers_constant_error (SchroEncoderFrame * frame);
20
21 double schro_encoder_entropy_to_lambda (SchroEncoderFrame * frame,
22 double entropy);
23 static double schro_encoder_lambda_to_entropy (SchroEncoderFrame * frame,
24 double lambda);
25 static void schro_encoder_generate_subband_histograms (SchroEncoderFrame *
26 frame);
27
28
29 static int schro_subband_pick_quant (SchroEncoderFrame * frame,
30 int component, int i, double lambda);
31
32 #define CURVE_SIZE 128
33
34 double
schro_encoder_perceptual_weight_moo(double cpd)35 schro_encoder_perceptual_weight_moo (double cpd)
36 {
37 /* I pretty much pulled this out of my ass (ppd = pixels per degree) */
38 if (cpd < 4)
39 return 1;
40 return 0.68 * cpd * exp (-0.25 * cpd);
41 }
42
43 double
schro_encoder_perceptual_weight_constant(double cpd)44 schro_encoder_perceptual_weight_constant (double cpd)
45 {
46 return 1;
47 }
48
49 double
schro_encoder_perceptual_weight_ccir959(double cpd)50 schro_encoder_perceptual_weight_ccir959 (double cpd)
51 {
52 double w;
53 w = 0.255 * pow (1 + 0.2561 * cpd * cpd, -0.75);
54
55 /* return normalized value */
56 return w / 0.255;
57 }
58
59 double
schro_encoder_perceptual_weight_manos_sakrison(double cpd)60 schro_encoder_perceptual_weight_manos_sakrison (double cpd)
61 {
62 double w;
63
64 if (cpd < 4)
65 return 1;
66 w = 2.6 * (0.0192 + 0.114 * cpd) * exp (-pow (0.114 * cpd, 1.1));
67
68 /* return normalized value */
69 return w / 0.980779694777866;
70 }
71
72 static double
weighted_sum(const float * h1,const float * v1,double * weight)73 weighted_sum (const float *h1, const float *v1, double *weight)
74 {
75 int i, j;
76 double sum;
77 double rowsum;
78
79 sum = 0;
80 for (j = 0; j < CURVE_SIZE; j++) {
81 rowsum = 0;
82 for (i = 0; i < CURVE_SIZE; i++) {
83 rowsum += h1[i] * v1[j] * weight[CURVE_SIZE * j + i];
84 }
85 sum += rowsum;
86 }
87 return sum;
88 }
89
90 static double
dot_product(const float * h1,const float * v1,const float * h2,const float * v2,double * weight)91 dot_product (const float *h1, const float *v1, const float *h2,
92 const float *v2, double *weight)
93 {
94 int i, j;
95 double sum;
96 double rowsum;
97
98 sum = 0;
99 for (j = 0; j < CURVE_SIZE; j++) {
100 rowsum = 0;
101 for (i = 0; i < CURVE_SIZE; i++) {
102 rowsum +=
103 h1[i] * v1[j] * h2[i] * v2[j] * weight[CURVE_SIZE * j +
104 i] * weight[CURVE_SIZE * j + i];
105 }
106 sum += rowsum;
107 }
108 return sum;
109 }
110
111 static void
solve(double * matrix,double * column,int n)112 solve (double *matrix, double *column, int n)
113 {
114 int i;
115 int j;
116 int k;
117 double x;
118
119 for (i = 0; i < n; i++) {
120 x = 1 / matrix[i * n + i];
121 for (k = i; k < n; k++) {
122 matrix[i * n + k] *= x;
123 }
124 column[i] *= x;
125
126 for (j = i + 1; j < n; j++) {
127 x = matrix[j * n + i];
128 for (k = i; k < n; k++) {
129 matrix[j * n + k] -= matrix[i * n + k] * x;
130 }
131 column[j] -= column[i] * x;
132 }
133 }
134
135 for (i = n - 1; i > 0; i--) {
136 for (j = i - 1; j >= 0; j--) {
137 column[j] -= matrix[j * n + i] * column[i];
138 matrix[j * n + i] = 0;
139 }
140 }
141 }
142
143 #ifdef unused
144 void
schro_encoder_set_default_subband_weights(SchroEncoder * encoder)145 schro_encoder_set_default_subband_weights (SchroEncoder * encoder)
146 {
147 schro_encoder_calculate_subband_weights (encoder,
148 schro_encoder_perceptual_weight_constant);
149 }
150 #endif
151
152 void
schro_encoder_calculate_subband_weights(SchroEncoder * encoder,double (* perceptual_weight)(double))153 schro_encoder_calculate_subband_weights (SchroEncoder * encoder,
154 double (*perceptual_weight) (double))
155 {
156 int wavelet;
157 int n_levels;
158 double *matrix_intra;
159 double *matrix_inter;
160 int n;
161 int i, j;
162 double column_intra[SCHRO_LIMIT_SUBBANDS];
163 double column_inter[SCHRO_LIMIT_SUBBANDS];
164 double *weight_intra;
165 double *weight_inter;
166
167 matrix_intra =
168 schro_malloc (sizeof (double) * SCHRO_LIMIT_SUBBANDS *
169 SCHRO_LIMIT_SUBBANDS);
170 matrix_inter =
171 schro_malloc (sizeof (double) * SCHRO_LIMIT_SUBBANDS *
172 SCHRO_LIMIT_SUBBANDS);
173 weight_intra = schro_malloc (sizeof (double) * CURVE_SIZE * CURVE_SIZE);
174 weight_inter = schro_malloc (sizeof (double) * CURVE_SIZE * CURVE_SIZE);
175
176 for (j = 0; j < CURVE_SIZE; j++) {
177 for (i = 0; i < CURVE_SIZE; i++) {
178 double fv_intra =
179 j * encoder->cycles_per_degree_vert * (1.0 / CURVE_SIZE);
180 double fh_intra =
181 i * encoder->cycles_per_degree_horiz * (1.0 / CURVE_SIZE);
182
183 double fv_inter =
184 j * encoder->cycles_per_degree_vert * (1.0 / CURVE_SIZE) *
185 encoder->magic_inter_cpd_scale;
186 double fh_inter =
187 i * encoder->cycles_per_degree_horiz * (1.0 / CURVE_SIZE) *
188 encoder->magic_inter_cpd_scale;
189
190 weight_intra[j * CURVE_SIZE + i] =
191 perceptual_weight (sqrt (fv_intra * fv_intra + fh_intra * fh_intra));
192 weight_inter[j * CURVE_SIZE + i] =
193 perceptual_weight (sqrt (fv_intra * fv_inter + fh_inter * fh_inter));
194 }
195 }
196
197 for (wavelet = 0; wavelet < SCHRO_N_WAVELETS; wavelet++) {
198 for (n_levels = 1; n_levels <= SCHRO_LIMIT_ENCODER_TRANSFORM_DEPTH;
199 n_levels++) {
200 const float *h_curve[SCHRO_LIMIT_SUBBANDS];
201 const float *v_curve[SCHRO_LIMIT_SUBBANDS];
202 int hi[SCHRO_LIMIT_SUBBANDS];
203 int vi[SCHRO_LIMIT_SUBBANDS];
204
205 n = 3 * n_levels + 1;
206
207 for (i = 0; i < n; i++) {
208 int position = schro_subband_get_position (i);
209 int n_transforms;
210
211 n_transforms = n_levels - SCHRO_SUBBAND_SHIFT (position);
212 if (position & 1) {
213 hi[i] = (n_transforms - 1) * 2;
214 } else {
215 hi[i] = (n_transforms - 1) * 2 + 1;
216 }
217 if (position & 2) {
218 vi[i] = (n_transforms - 1) * 2;
219 } else {
220 vi[i] = (n_transforms - 1) * 2 + 1;
221 }
222 h_curve[i] = schro_tables_wavelet_noise_curve[wavelet][hi[i]];
223 v_curve[i] = schro_tables_wavelet_noise_curve[wavelet][vi[i]];
224 }
225
226 if (0) {
227 for (i = 0; i < n; i++) {
228 column_intra[i] = weighted_sum (h_curve[i], v_curve[i], weight_intra);
229 column_inter[i] = weighted_sum (h_curve[i], v_curve[i], weight_inter);
230 matrix_intra[i * n + i] = dot_product (h_curve[i], v_curve[i],
231 h_curve[i], v_curve[i], weight_intra);
232 matrix_inter[i * n + i] = dot_product (h_curve[i], v_curve[i],
233 h_curve[i], v_curve[i], weight_inter);
234 for (j = i + 1; j < n; j++) {
235 matrix_intra[i * n + j] = dot_product (h_curve[i], v_curve[i],
236 h_curve[j], v_curve[j], weight_intra);
237 matrix_intra[j * n + i] = matrix_intra[i * n + j];
238 matrix_inter[i * n + j] = dot_product (h_curve[i], v_curve[i],
239 h_curve[j], v_curve[j], weight_inter);
240 matrix_inter[j * n + i] = matrix_inter[i * n + j];
241 }
242 }
243
244 solve (matrix_intra, column_intra, n);
245 solve (matrix_inter, column_inter, n);
246
247 for (i = 0; i < n; i++) {
248 if (column_intra[i] < 0 || column_inter[i] < 0) {
249 SCHRO_ERROR ("BROKEN wavelet %d n_levels %d", wavelet, n_levels);
250 break;
251 }
252 }
253
254 SCHRO_DEBUG ("wavelet %d n_levels %d", wavelet, n_levels);
255 for (i = 0; i < n; i++) {
256 SCHRO_DEBUG ("%g", 1.0 / sqrt (column_intra[i]));
257 SCHRO_DEBUG ("%g", 1.0 / sqrt (column_inter[i]));
258 encoder->intra_subband_weights[wavelet][n_levels - 1][i] =
259 sqrt (column_intra[i]);
260 encoder->inter_subband_weights[wavelet][n_levels - 1][i] =
261 sqrt (column_inter[i]);
262 }
263 } else {
264 for (i = 0; i < n; i++) {
265 int position = schro_subband_get_position (i);
266 int n_transforms;
267 double size;
268
269 n_transforms = n_levels - SCHRO_SUBBAND_SHIFT (position);
270 size = (1.0 / CURVE_SIZE) * (1 << n_transforms);
271 encoder->intra_subband_weights[wavelet][n_levels - 1][i] =
272 1.0 / (size * sqrt (weighted_sum (h_curve[i], v_curve[i],
273 weight_intra)));
274 encoder->inter_subband_weights[wavelet][n_levels - 1][i] =
275 1.0 / (size * sqrt (weighted_sum (h_curve[i], v_curve[i],
276 weight_inter)));
277 }
278 }
279 }
280 }
281
282 schro_free (weight_intra);
283 schro_free (matrix_intra);
284 schro_free (weight_inter);
285 schro_free (matrix_inter);
286 }
287
288 void
schro_encoder_choose_quantisers(SchroEncoderFrame * frame)289 schro_encoder_choose_quantisers (SchroEncoderFrame * frame)
290 {
291 switch (frame->encoder->quantiser_engine) {
292 case SCHRO_QUANTISER_ENGINE_SIMPLE:
293 schro_encoder_choose_quantisers_simple (frame);
294 break;
295 case SCHRO_QUANTISER_ENGINE_RDO_BIT_ALLOCATION:
296 schro_encoder_choose_quantisers_rdo_bit_allocation (frame);
297 break;
298 case SCHRO_QUANTISER_ENGINE_CBR:
299 schro_encoder_choose_quantisers_rdo_cbr (frame);
300 break;
301 case SCHRO_QUANTISER_ENGINE_LOSSLESS:
302 schro_encoder_choose_quantisers_lossless (frame);
303 break;
304 case SCHRO_QUANTISER_ENGINE_LOWDELAY:
305 schro_encoder_choose_quantisers_lowdelay (frame);
306 break;
307 case SCHRO_QUANTISER_ENGINE_RDO_LAMBDA:
308 schro_encoder_choose_quantisers_rdo_lambda (frame);
309 break;
310 case SCHRO_QUANTISER_ENGINE_CONSTANT_ERROR:
311 schro_encoder_choose_quantisers_constant_error (frame);
312 break;
313 default:
314 SCHRO_ASSERT (0);
315 }
316 }
317
318 void
schro_encoder_choose_quantisers_lossless(SchroEncoderFrame * frame)319 schro_encoder_choose_quantisers_lossless (SchroEncoderFrame * frame)
320 {
321 int i;
322 int component;
323
324 for (component = 0; component < 3; component++) {
325 for (i = 0; i < 1 + 3 * frame->params.transform_depth; i++) {
326 schro_encoder_frame_set_quant_index (frame, component, i, -1, -1, 0);
327 }
328 }
329 }
330
331 void
schro_encoder_choose_quantisers_simple(SchroEncoderFrame * frame)332 schro_encoder_choose_quantisers_simple (SchroEncoderFrame * frame)
333 {
334 SchroParams *params = &frame->params;
335 int i;
336 int component;
337 double noise_amplitude;
338 double a;
339 double max;
340 double *table;
341 double range;
342
343 /* FIXME this should really be based directly on excursion */
344 range = (1<<frame->encoder->bit_depth) - 1.0;
345
346 noise_amplitude = range * pow (0.1, frame->encoder->noise_threshold * 0.05);
347 SCHRO_DEBUG ("noise %g", noise_amplitude);
348
349 if (frame->num_refs == 0) {
350 table = frame->encoder->intra_subband_weights[params->wavelet_filter_index]
351 [MAX (0, params->transform_depth - 1)];
352 } else {
353 table = frame->encoder->inter_subband_weights[params->wavelet_filter_index]
354 [MAX (0, params->transform_depth - 1)];
355 }
356
357 for (component = 0; component < 3; component++) {
358 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
359 a = noise_amplitude * table[i];
360
361 schro_encoder_frame_set_quant_index (frame, component, i, -1, -1,
362 schro_utils_multiplier_to_quant_index (a));
363 }
364 }
365
366 max = 1.0;
367
368 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
369 params->quant_matrix[i] =
370 schro_utils_multiplier_to_quant_index (max / table[i]);
371 SCHRO_DEBUG ("%g %g %d", table[i], max / table[i], params->quant_matrix[i]);
372 }
373 }
374
375 void
schro_encoder_choose_quantisers_lowdelay(SchroEncoderFrame * frame)376 schro_encoder_choose_quantisers_lowdelay (SchroEncoderFrame * frame)
377 {
378 SchroParams *params = &frame->params;
379 int i;
380 int component;
381 int base;
382 const int *table;
383
384 /* completely made up */
385 base = 12 + (30 - frame->encoder->noise_threshold) / 2;
386
387 table = schro_tables_lowdelay_quants[params->wavelet_filter_index]
388 [MAX (0, params->transform_depth - 1)];
389
390 for (component = 0; component < 3; component++) {
391 schro_encoder_frame_set_quant_index (frame, component, 0, -1, -1,
392 base - table[0]);
393
394 for (i = 0; i < params->transform_depth; i++) {
395 schro_encoder_frame_set_quant_index (frame, component, 1 + 3 * i + 0, -1,
396 -1, base - table[1 + 2 * i + 0]);
397 schro_encoder_frame_set_quant_index (frame, component, 1 + 3 * i + 1, -1,
398 -1, base - table[1 + 2 * i + 0]);
399 schro_encoder_frame_set_quant_index (frame, component, 1 + 3 * i + 2, -1,
400 -1, base - table[1 + 2 * i + 1]);
401 }
402 }
403
404 }
405
406 #ifdef DUMP_SUBBAND_CURVES
407 static double
pow2(double x)408 pow2 (double x)
409 {
410 return x * x;
411 }
412 #endif
413
414 #ifdef DUMP_SUBBAND_CURVES
415 static double
measure_error_subband(SchroEncoderFrame * frame,int component,int index,int quant_index)416 measure_error_subband (SchroEncoderFrame * frame, int component, int index,
417 int quant_index)
418 {
419 int i;
420 int j;
421 int16_t *line;
422 int skip = 1;
423 double error = 0;
424 int q;
425 int quant_factor;
426 int quant_offset;
427 int value;
428 int position;
429 SchroFrameData fd;
430
431 position = schro_subband_get_position (index);
432 schro_subband_get_frame_data (&fd, frame->iwt_frame, component, position,
433 &frame->params);
434
435 quant_factor = schro_table_quant[quant_index];
436 if (frame->params.num_refs > 0) {
437 quant_offset = schro_table_offset_3_8[quant_index];
438 } else {
439 quant_offset = schro_table_offset_1_2[quant_index];
440 }
441
442 error = 0;
443 if (index == 0) {
444 for (j = 0; j < fd.height; j += skip) {
445 line = SCHRO_FRAME_DATA_GET_LINE (&fd, j);
446 for (i = 1; i < fd.width; i += skip) {
447 q = schro_quantise (abs (line[i] - line[i - 1]), quant_factor,
448 quant_offset);
449 value = schro_dequantise (q, quant_factor, quant_offset);
450 error += pow2 (value - abs (line[i] - line[i - 1]));
451 }
452 }
453 } else {
454 for (j = 0; j < fd.height; j += skip) {
455 line = SCHRO_FRAME_DATA_GET_LINE (&fd, j);
456 for (i = 0; i < fd.width; i += skip) {
457 q = schro_quantise (line[i], quant_factor, quant_offset);
458 value = schro_dequantise (q, quant_factor, quant_offset);
459 error += pow2 (value - line[i]);
460 }
461 }
462 }
463 error *= skip * skip;
464
465 return error;
466 }
467 #endif
468
469 typedef struct _ErrorFuncInfo ErrorFuncInfo;
470 struct _ErrorFuncInfo
471 {
472 int quant_factor;
473 int quant_offset;
474 double power;
475 };
476
477 static double
error_pow(int x,void * priv)478 error_pow (int x, void *priv)
479 {
480 ErrorFuncInfo *efi = priv;
481 int q;
482 int value;
483 int y;
484
485 q = schro_quantise (x, efi->quant_factor, efi->quant_offset);
486 value = schro_dequantise (q, efi->quant_factor, efi->quant_offset);
487
488 y = abs (value - x);
489
490 return pow (y, efi->power);
491 }
492
493 void
schro_encoder_init_error_tables(SchroEncoder * encoder)494 schro_encoder_init_error_tables (SchroEncoder * encoder)
495 {
496 int i;
497
498 for (i = 0; i < 60; i++) {
499 ErrorFuncInfo efi;
500
501 efi.quant_factor = schro_table_quant[i];
502 efi.quant_offset = schro_table_offset_1_2[i];
503 efi.power = encoder->magic_error_power;
504
505 schro_histogram_table_generate (encoder->intra_hist_tables + i,
506 error_pow, &efi);
507 }
508
509 }
510
511 #ifdef unused
512 static double
schro_histogram_estimate_error(SchroHistogram * hist,int quant_index,int num_refs)513 schro_histogram_estimate_error (SchroHistogram * hist, int quant_index,
514 int num_refs)
515 {
516 SchroHistogramTable *table;
517
518 if (num_refs == 0) {
519 table =
520 (SchroHistogramTable
521 *) (schro_table_error_hist_shift3_1_2[quant_index]);
522 } else {
523 /* FIXME the 3/8 table doesn't exist yet */
524 //table = (SchroHistogramTable *)(schro_table_error_hist_shift3_3_8[quant_index]);
525 table =
526 (SchroHistogramTable
527 *) (schro_table_error_hist_shift3_1_2[quant_index]);
528 }
529 return schro_histogram_apply_table (hist, table);
530 }
531 #endif
532
533 #ifdef DUMP_SUBBAND_CURVES
534 static double
schro_encoder_estimate_subband_arith(SchroEncoderFrame * frame,int component,int index,int quant_index)535 schro_encoder_estimate_subband_arith (SchroEncoderFrame * frame, int component,
536 int index, int quant_index)
537 {
538 int i;
539 int j;
540 int16_t *line;
541 int q;
542 int quant_factor;
543 int quant_offset;
544 int estimated_entropy;
545 SchroArith *arith;
546 int position;
547 SchroFrameData fd;
548
549 arith = schro_arith_new ();
550 schro_arith_estimate_init (arith);
551
552 position = schro_subband_get_position (index);
553 schro_subband_get_frame_data (&fd, frame->iwt_frame, component, position,
554 &frame->params);
555
556 quant_factor = schro_table_quant[quant_index];
557 quant_offset = schro_table_offset_1_2[quant_index];
558
559 if (index == 0) {
560 for (j = 0; j < fd.height; j++) {
561 line = SCHRO_FRAME_DATA_GET_LINE (&fd, j);
562 schro_arith_estimate_sint (arith,
563 SCHRO_CTX_ZPZN_F1, SCHRO_CTX_COEFF_DATA, SCHRO_CTX_SIGN_ZERO, 0);
564 for (i = 1; i < fd.width; i++) {
565 q = schro_quantise (line[i] - line[i - 1], quant_factor, quant_offset);
566 schro_arith_estimate_sint (arith,
567 SCHRO_CTX_ZPZN_F1, SCHRO_CTX_COEFF_DATA, SCHRO_CTX_SIGN_ZERO, q);
568 }
569 }
570 } else {
571 for (j = 0; j < fd.height; j++) {
572 line = SCHRO_FRAME_DATA_GET_LINE (&fd, j);
573 for (i = 0; i < fd.width; i++) {
574 q = schro_quantise (line[i], quant_factor, quant_offset);
575 schro_arith_estimate_sint (arith,
576 SCHRO_CTX_ZPZN_F1, SCHRO_CTX_COEFF_DATA, SCHRO_CTX_SIGN_ZERO, q);
577 }
578 }
579 }
580
581 estimated_entropy = 0;
582
583 estimated_entropy += arith->contexts[SCHRO_CTX_ZPZN_F1].n_bits;
584 estimated_entropy += arith->contexts[SCHRO_CTX_ZP_F2].n_bits;
585 estimated_entropy += arith->contexts[SCHRO_CTX_ZP_F3].n_bits;
586 estimated_entropy += arith->contexts[SCHRO_CTX_ZP_F4].n_bits;
587 estimated_entropy += arith->contexts[SCHRO_CTX_ZP_F5].n_bits;
588 estimated_entropy += arith->contexts[SCHRO_CTX_ZP_F6p].n_bits;
589
590 estimated_entropy += arith->contexts[SCHRO_CTX_COEFF_DATA].n_bits;
591
592 estimated_entropy += arith->contexts[SCHRO_CTX_SIGN_ZERO].n_bits;
593
594 schro_arith_free (arith);
595
596 return estimated_entropy;
597 }
598 #endif
599
600 static void
schro_encoder_generate_subband_histogram(SchroEncoderFrame * frame,int component,int index,SchroHistogram * hist,int skip)601 schro_encoder_generate_subband_histogram (SchroEncoderFrame * frame,
602 int component, int index, SchroHistogram * hist, int skip)
603 {
604 int position;
605 SchroFrameData fd;
606
607 position = schro_subband_get_position (index);
608 schro_subband_get_frame_data (&fd, frame->iwt_frame, component, position,
609 &frame->params);
610
611 if (index == 0 && frame->num_refs == 0) {
612 schro_frame_data_generate_histogram_dc_predict (&fd, hist, skip, 0, 0);
613 } else {
614 schro_frame_data_generate_histogram (&fd, hist, skip);
615 }
616 }
617
618 static void
schro_encoder_generate_subband_histograms(SchroEncoderFrame * frame)619 schro_encoder_generate_subband_histograms (SchroEncoderFrame * frame)
620 {
621 SchroParams *params = &frame->params;
622 int i;
623 int component;
624 int pos;
625 int skip;
626
627 for (component = 0; component < 3; component++) {
628 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
629 pos = schro_subband_get_position (i);
630 skip = 1 << MAX (0, SCHRO_SUBBAND_SHIFT (pos) - 1);
631 schro_encoder_generate_subband_histogram (frame, component, i,
632 &frame->subband_hists[component][i], skip);
633 }
634 }
635
636 frame->have_histograms = TRUE;
637 }
638
639 #ifdef DUMP_SUBBAND_CURVES
640 static void
schro_encoder_dump_subband_curves(SchroEncoderFrame * frame)641 schro_encoder_dump_subband_curves (SchroEncoderFrame * frame)
642 {
643 SchroParams *params = &frame->params;
644 int i;
645 int component;
646 int position;
647
648 SCHRO_ASSERT (frame->have_histograms);
649
650 for (component = 0; component < 3; component++) {
651 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
652 int vol;
653 SchroFrameData fd;
654 int j;
655 SchroHistogram *hist;
656
657 position = schro_subband_get_position (i);
658 schro_subband_get_frame_data (&fd, frame->iwt_frame, component,
659 position, &frame->params);
660 vol = fd.width * fd.height;
661 hist = &frame->subband_hists[component][i];
662
663 for (j = 0; j < 60; j++) {
664 double est_entropy;
665 double error;
666 double est_error;
667 double arith_entropy;
668
669 error = measure_error_subband (frame, component, i, j);
670 est_entropy =
671 schro_histogram_estimate_entropy (&frame->subband_hists[component]
672 [i], j, params->is_noarith);
673 est_error =
674 schro_histogram_apply_table (hist,
675 &frame->encoder->intra_hist_tables[j]);
676 arith_entropy =
677 schro_encoder_estimate_subband_arith (frame, component, i, j);
678
679 schro_dump (SCHRO_DUMP_SUBBAND_CURVE, "%d %d %d %g %g %g %g\n",
680 component, i, j, est_entropy / vol, arith_entropy / vol,
681 est_error / vol, error / vol);
682 }
683 }
684 }
685 }
686 #endif
687
688 static void
schro_encoder_calc_estimates(SchroEncoderFrame * frame)689 schro_encoder_calc_estimates (SchroEncoderFrame * frame)
690 {
691 SchroParams *params = &frame->params;
692 int i;
693 int j;
694 int component;
695 double *arith_context_ratios;
696
697 SCHRO_ASSERT (frame->have_histograms);
698
699 #ifdef DUMP_SUBBAND_CURVES
700 schro_encoder_dump_subband_curves (frame);
701 #endif
702
703 for (component = 0; component < 3; component++) {
704 if (frame->num_refs == 0) {
705 arith_context_ratios =
706 frame->encoder->average_arith_context_ratios_intra[component];
707 } else {
708 arith_context_ratios =
709 frame->encoder->average_arith_context_ratios_inter[component];
710 }
711
712 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
713 for (j = 0; j < 60; j++) {
714 int position;
715 SchroHistogram *hist;
716 SchroFrameData fd;
717
718 position = schro_subband_get_position (i);
719 schro_subband_get_frame_data (&fd, frame->iwt_frame, component,
720 position, &frame->params);
721
722 hist = &frame->subband_hists[component][i];
723 frame->est_entropy[component][i][j] =
724 schro_histogram_estimate_entropy (hist, j, params->is_noarith);
725 frame->est_entropy[component][i][j] *= arith_context_ratios[i];
726 frame->est_error[component][i][j] = schro_histogram_apply_table (hist,
727 &frame->encoder->intra_hist_tables[j]);
728 }
729 }
730 }
731 frame->have_estimate_tables = TRUE;
732 }
733
734 /*
735 * Quantiser engine which picks the best RDO quantisers to fit in
736 * a frame bit allocation.
737 */
738 void
schro_encoder_choose_quantisers_rdo_bit_allocation(SchroEncoderFrame * frame)739 schro_encoder_choose_quantisers_rdo_bit_allocation (SchroEncoderFrame * frame)
740 {
741 double frame_lambda;
742 int bits;
743
744 schro_encoder_generate_subband_histograms (frame);
745 schro_encoder_calc_estimates (frame);
746
747 SCHRO_ASSERT (frame->have_estimate_tables);
748
749 bits = frame->allocated_residual_bits;
750
751 frame_lambda = schro_encoder_entropy_to_lambda (frame, bits);
752 frame->frame_lambda = frame_lambda;
753 SCHRO_DEBUG ("LAMBDA: %d %g %d", frame->frame_number, frame_lambda, bits);
754
755 schro_encoder_lambda_to_entropy (frame, frame_lambda);
756 }
757
758 void
schro_encoder_choose_quantisers_rdo_lambda(SchroEncoderFrame * frame)759 schro_encoder_choose_quantisers_rdo_lambda (SchroEncoderFrame * frame)
760 {
761 SCHRO_DEBUG ("Using rdo_lambda quant selection on frame %d with lambda %g",
762 frame->frame_number, frame->frame_lambda);
763
764 schro_encoder_generate_subband_histograms (frame);
765 schro_encoder_calc_estimates (frame);
766
767 SCHRO_ASSERT (frame->have_estimate_tables);
768
769 schro_encoder_lambda_to_entropy (frame, frame->frame_lambda);
770 }
771
772
773 void
schro_encoder_choose_quantisers_rdo_cbr(SchroEncoderFrame * frame)774 schro_encoder_choose_quantisers_rdo_cbr (SchroEncoderFrame * frame)
775 {
776 schro_encoder_generate_subband_histograms (frame);
777 schro_encoder_calc_estimates (frame);
778
779 SCHRO_ASSERT (frame->have_estimate_tables);
780
781 schro_encoder_lambda_to_entropy (frame, frame->frame_lambda);
782 }
783
784 void
schro_encoder_estimate_entropy(SchroEncoderFrame * frame)785 schro_encoder_estimate_entropy (SchroEncoderFrame * frame)
786 {
787 SchroParams *params = &frame->params;
788 int i;
789 int component;
790 int n = 0;
791
792 for (component = 0; component < 3; component++) {
793 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
794 n += frame->
795 est_entropy[component][i][frame->quant_indices[component][i][0]];
796 }
797 }
798 frame->estimated_residual_bits = n;
799
800 if (frame->allocated_residual_bits > 0 &&
801 frame->estimated_residual_bits >
802 2 * frame->encoder->bits_per_picture + frame->allocated_residual_bits) {
803 SCHRO_WARNING ("%d: estimated entropy too big (%d vs %d)",
804 frame->frame_number,
805 frame->estimated_residual_bits, frame->allocated_residual_bits);
806 }
807 }
808
809 static int
schro_subband_pick_quant(SchroEncoderFrame * frame,int component,int i,double lambda)810 schro_subband_pick_quant (SchroEncoderFrame * frame, int component, int i,
811 double lambda)
812 {
813 double x;
814 double min;
815 int j;
816 int j_min;
817 double entropy;
818 double error;
819
820 SCHRO_ASSERT (frame->have_estimate_tables);
821
822 j_min = -1;
823 min = 0;
824 for (j = 0; j < 60; j++) {
825 entropy = frame->est_entropy[component][i][j];
826 error = frame->est_error[component][i][j];
827
828 x = entropy + lambda * error;
829 if (j == 0 || x < min) {
830 j_min = j;
831 min = x;
832 }
833 }
834
835 return j_min;
836 }
837
838 static double
schro_encoder_lambda_to_entropy(SchroEncoderFrame * frame,double frame_lambda)839 schro_encoder_lambda_to_entropy (SchroEncoderFrame * frame, double frame_lambda)
840 {
841 SchroParams *params = &frame->params;
842 int i;
843 int component;
844 double entropy = 0;
845 double *table;
846
847 if (frame->num_refs == 0) {
848 table = frame->encoder->intra_subband_weights[params->wavelet_filter_index]
849 [MAX (0, params->transform_depth - 1)];
850 } else {
851 table = frame->encoder->inter_subband_weights[params->wavelet_filter_index]
852 [MAX (0, params->transform_depth - 1)];
853 }
854
855 for (component = 0; component < 3; component++) {
856 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
857 double lambda;
858 double weight;
859 int quant_index;
860 int position = schro_subband_get_position (i);
861
862 lambda = frame_lambda;
863
864 if (i == 0) {
865 lambda *= frame->encoder->magic_subband0_lambda_scale;
866 }
867 if (component > 0) {
868 lambda *= frame->encoder->magic_chroma_lambda_scale;
869 }
870 if (SCHRO_SUBBAND_IS_DIAGONALLY_ORIENTED (position)) {
871 lambda *= frame->encoder->magic_diagonal_lambda_scale;
872 }
873
874 weight = table[i];
875 lambda /= weight * weight;
876
877 quant_index = schro_subband_pick_quant (frame, component, i, lambda);
878 entropy += frame->est_entropy[component][i][quant_index];
879 schro_encoder_frame_set_quant_index (frame, component, i, -1, -1,
880 quant_index);
881 }
882 }
883
884 return entropy;
885 }
886
887 double
schro_encoder_entropy_to_lambda(SchroEncoderFrame * frame,double entropy)888 schro_encoder_entropy_to_lambda (SchroEncoderFrame * frame, double entropy)
889 {
890 int j;
891 double lambda_hi, lambda_lo, lambda_mid;
892 double entropy_hi, entropy_lo, entropy_mid;
893
894 lambda_hi = 1;
895 entropy_hi = schro_encoder_lambda_to_entropy (frame, lambda_hi);
896 SCHRO_DEBUG ("start target=%g lambda=%g entropy=%g",
897 entropy, lambda_hi, entropy_hi);
898
899 if (entropy_hi < entropy) {
900 entropy_lo = entropy_hi;
901 lambda_lo = lambda_hi;
902
903 for (j = 0; j < 5; j++) {
904 lambda_hi = lambda_lo * 100;
905 entropy_hi = schro_encoder_lambda_to_entropy (frame, lambda_hi);
906
907 SCHRO_DEBUG ("have: lambda=[%g,%g] entropy=[%g,%g] target=%g",
908 lambda_lo, lambda_hi, entropy_lo, entropy_hi, entropy);
909 if (entropy_hi > entropy)
910 break;
911
912 SCHRO_DEBUG ("--> step up");
913
914 entropy_lo = entropy_hi;
915 lambda_lo = lambda_hi;
916 }
917 SCHRO_DEBUG ("--> stopping");
918 } else {
919 for (j = 0; j < 5; j++) {
920 lambda_lo = lambda_hi * 0.01;
921 entropy_lo = schro_encoder_lambda_to_entropy (frame, lambda_lo);
922
923 SCHRO_DEBUG ("have: lambda=[%g,%g] entropy=[%g,%g] target=%g",
924 lambda_lo, lambda_hi, entropy_lo, entropy_hi, entropy);
925
926 SCHRO_DEBUG ("--> step down");
927 if (entropy_lo < entropy)
928 break;
929
930 entropy_hi = entropy_lo;
931 lambda_hi = lambda_lo;
932 }
933 SCHRO_DEBUG ("--> stopping");
934 }
935 if (entropy_lo == entropy_hi) {
936 return sqrt (lambda_lo * lambda_hi);
937 }
938
939 if (entropy_lo > entropy || entropy_hi < entropy) {
940 SCHRO_ERROR ("entropy not bracketed");
941 }
942
943 for (j = 0; j < 7; j++) {
944 if (entropy_hi == entropy_lo)
945 break;
946
947 SCHRO_DEBUG ("have: lambda=[%g,%g] entropy=[%g,%g] target=%g",
948 lambda_lo, lambda_hi, entropy_lo, entropy_hi, entropy);
949
950 lambda_mid = sqrt (lambda_lo * lambda_hi);
951 entropy_mid = schro_encoder_lambda_to_entropy (frame, lambda_mid);
952
953 SCHRO_DEBUG ("picking lambda_mid=%g entropy=%g", lambda_mid, entropy_mid);
954
955 if (entropy_mid > entropy) {
956 lambda_hi = lambda_mid;
957 entropy_hi = entropy_mid;
958 SCHRO_DEBUG ("--> focus up");
959 } else {
960 lambda_lo = lambda_mid;
961 entropy_lo = entropy_mid;
962 SCHRO_DEBUG ("--> focus down");
963 }
964 }
965
966 lambda_mid = sqrt (lambda_hi * lambda_lo);
967 SCHRO_DEBUG ("done %g", lambda_mid);
968 return lambda_mid;
969 }
970
971 static double
schro_encoder_lambda_to_error(SchroEncoderFrame * frame,double frame_lambda)972 schro_encoder_lambda_to_error (SchroEncoderFrame * frame, double frame_lambda)
973 {
974 SchroParams *params = &frame->params;
975 int i;
976 int component;
977 double error = 0;
978 double *table;
979
980 if (frame->num_refs == 0) {
981 table = frame->encoder->intra_subband_weights[params->wavelet_filter_index]
982 [MAX (0, params->transform_depth - 1)];
983 } else {
984 table = frame->encoder->inter_subband_weights[params->wavelet_filter_index]
985 [MAX (0, params->transform_depth - 1)];
986 }
987
988 for (component = 0; component < 3; component++) {
989 for (i = 0; i < 1 + 3 * params->transform_depth; i++) {
990 double lambda;
991 double weight;
992 int quant_index;
993
994 lambda = frame_lambda;
995
996 if (i == 0) {
997 lambda *= frame->encoder->magic_subband0_lambda_scale;
998 }
999 if (component > 0) {
1000 lambda *= frame->encoder->magic_chroma_lambda_scale;
1001 }
1002
1003 weight = table[i];
1004 lambda /= weight * weight;
1005
1006 quant_index = schro_subband_pick_quant (frame, component, i, lambda);
1007 error += frame->est_error[component][i][quant_index];
1008 schro_encoder_frame_set_quant_index (frame, component, i, -1, -1,
1009 quant_index);
1010 }
1011 }
1012
1013 return error;
1014 }
1015
1016 static double
schro_encoder_error_to_lambda(SchroEncoderFrame * frame,double error)1017 schro_encoder_error_to_lambda (SchroEncoderFrame * frame, double error)
1018 {
1019 int j;
1020 double lambda_hi, lambda_lo, lambda_mid;
1021 double error_hi, error_lo, error_mid;
1022
1023 lambda_lo = 1;
1024 error_lo = schro_encoder_lambda_to_error (frame, lambda_lo);
1025 SCHRO_DEBUG ("start target=%g lambda=%g error=%g",
1026 error, lambda_lo, error_lo, lambda_lo, error);
1027
1028 if (error < error_lo) {
1029 error_hi = error_lo;
1030 lambda_hi = lambda_lo;
1031
1032 for (j = 0; j < 5; j++) {
1033 lambda_lo = lambda_hi * 100;
1034 error_lo = schro_encoder_lambda_to_error (frame, lambda_lo);
1035
1036 SCHRO_DEBUG ("have: lambda=[%g,%g] error=[%g,%g] target=%g",
1037 lambda_lo, lambda_hi, error_lo, error_hi, error);
1038 if (error > error_lo)
1039 break;
1040
1041 SCHRO_DEBUG ("--> step up");
1042
1043 error_hi = error_lo;
1044 lambda_hi = lambda_lo;
1045 }
1046 SCHRO_DEBUG ("--> stopping");
1047 } else {
1048 for (j = 0; j < 5; j++) {
1049 lambda_hi = lambda_lo * 0.01;
1050 error_hi = schro_encoder_lambda_to_error (frame, lambda_hi);
1051
1052 SCHRO_DEBUG ("have: lambda=[%g,%g] error=[%g,%g] target=%g",
1053 lambda_lo, lambda_hi, error_lo, error_hi, error);
1054
1055 SCHRO_DEBUG ("--> step down");
1056 if (error < error_hi)
1057 break;
1058
1059 error_lo = error_hi;
1060 lambda_lo = lambda_hi;
1061 }
1062 SCHRO_DEBUG ("--> stopping");
1063 }
1064 if (error_lo == error_hi) {
1065 return sqrt (lambda_lo * lambda_hi);
1066 }
1067
1068 if (error_lo > error || error_hi < error) {
1069 SCHRO_DEBUG ("error not bracketed");
1070 }
1071
1072 for (j = 0; j < 14; j++) {
1073 if (error_hi == error_lo)
1074 break;
1075
1076 SCHRO_DEBUG ("have: lambda=[%g,%g] error=[%g,%g] target=%g",
1077 lambda_lo, lambda_hi, error_lo, error_hi, error);
1078
1079 lambda_mid = sqrt (lambda_lo * lambda_hi);
1080 error_mid = schro_encoder_lambda_to_error (frame, lambda_mid);
1081
1082 SCHRO_DEBUG ("picking lambda_mid=%g error=%g", lambda_mid, error_mid);
1083
1084 if (error_mid > error) {
1085 lambda_hi = lambda_mid;
1086 error_hi = error_mid;
1087 SCHRO_DEBUG ("--> focus up");
1088 } else {
1089 lambda_lo = lambda_mid;
1090 error_lo = error_mid;
1091 SCHRO_DEBUG ("--> focus down");
1092 }
1093 }
1094
1095 lambda_mid = sqrt (lambda_hi * lambda_lo);
1096 SCHRO_DEBUG ("done %g", lambda_mid);
1097 return lambda_mid;
1098 }
1099
1100 void
schro_encoder_choose_quantisers_constant_error(SchroEncoderFrame * frame)1101 schro_encoder_choose_quantisers_constant_error (SchroEncoderFrame * frame)
1102 {
1103 double frame_lambda;
1104 double error;
1105
1106 schro_encoder_generate_subband_histograms (frame);
1107 schro_encoder_calc_estimates (frame);
1108
1109 SCHRO_ASSERT (frame->have_estimate_tables);
1110
1111 error = 255.0 * pow (0.1, frame->encoder->noise_threshold * 0.05);
1112 error *= frame->params.video_format->width *
1113 frame->params.video_format->height;
1114
1115 frame_lambda = schro_encoder_error_to_lambda (frame, error);
1116
1117 frame->frame_lambda = frame_lambda;
1118 SCHRO_DEBUG ("LAMBDA: %d %g", frame->frame_number, frame_lambda);
1119 }
1120