/* * matrices.c: Output of transitions matrices * * Written by: Ullrich Hafner * * This file is part of FIASCO («F»ractal «I»mage «A»nd «S»equence «CO»dec) * Copyright (C) 1994-2000 Ullrich Hafner */ /* NETPBM: When you call delta_encoding() with last_domain < 4, it crashes. And we have seen it happen. */ /* * $Date: 2000/06/14 20:50:31 $ * $Author: hafner $ * $Revision: 5.1 $ * $State: Exp $ */ #include "config.h" #include #include "pm_c_util.h" #include "types.h" #include "macros.h" #include "error.h" #include "wfa.h" #include "bit-io.h" #include "arith.h" #include "misc.h" #include "wfalib.h" #include "matrices.h" /***************************************************************************** prototypes *****************************************************************************/ static unsigned delta_encoding (bool_t use_normal_domains, bool_t use_delta_domains, const wfa_t *wfa, unsigned last_domain, bitfile_t *output); static unsigned column_0_encoding (const wfa_t *wfa, unsigned last_row, bitfile_t *output); static unsigned chroma_encoding (const wfa_t *wfa, bitfile_t *output); /***************************************************************************** public code *****************************************************************************/ unsigned write_matrices (bool_t use_normal_domains, bool_t use_delta_domains, const wfa_t *wfa, bitfile_t *output) /* * Write transition matrices of 'wfa' to stream 'output'. * * Return value: * number of transitions encoded */ { unsigned root_state; /* root of luminance */ unsigned total = 0; /* number of transitions */ root_state = wfa->wfainfo->color ? wfa->tree [wfa->tree [wfa->root_state][0]][0] : wfa->root_state; total = column_0_encoding (wfa, root_state, output); total += delta_encoding (use_normal_domains, use_delta_domains, wfa, root_state, output); if (wfa->wfainfo->color) total += chroma_encoding (wfa, output); return total; } /***************************************************************************** private code *****************************************************************************/ static unsigned delta_encoding (bool_t use_normal_domains, bool_t use_delta_domains, const wfa_t *wfa, unsigned last_domain, bitfile_t *output) /* * Write transition matrices with delta coding to stream 'input'. * 'last_domain' is the maximum state number used as domain image. * * Return value: * number of non-zero matrix elements (WFA edges) */ { range_sort_t rs; /* ranges are sorted as in the coder */ unsigned max_domain; /* dummy used for recursion */ unsigned total = 0; /* * Generate a list of range blocks. * The order is the same as in the coder. */ rs.range_state = Calloc ((last_domain + 1) * MAXLABELS, sizeof (u_word_t)); rs.range_label = Calloc ((last_domain + 1) * MAXLABELS, sizeof (byte_t)); rs.range_max_domain = Calloc ((last_domain + 1) * MAXLABELS, sizeof (u_word_t)); rs.range_subdivided = Calloc ((last_domain + 1) * MAXLABELS, sizeof (bool_t)); rs.range_no = 0; max_domain = wfa->basis_states - 1; sort_ranges (last_domain, &max_domain, &rs, wfa); /* * Compute and write distribution of #edges */ { unsigned state, label; unsigned edge; unsigned count [MAXEDGES + 1]; unsigned n; unsigned edges = 0; unsigned M = 0; unsigned bits = bits_processed (output); for (n = 0; n < MAXEDGES + 1; n++) count [n] = 0; for (state = wfa->basis_states; state <= last_domain; state++) for (label = 0; label < MAXLABELS; label++) if (isrange (wfa->tree [state][label])) { for (edge = 0; isedge (wfa->into [state][label][edge]); edge++) ; count [edge]++; edges++; M = MAX(edge, M); } write_rice_code (M, 3, output); for (n = 0; n <= M; n++) /* NETPBM: The following causes a crash when last_domain < 4, because it requests writing of a negative number of bits. And we have seen last_domain = 3. But we have no clue what last_domain means, or even what a rice code is, so we don't know where the error lies. -Bryan 2001.02.09 */ write_rice_code (count [n], (int) log2 (last_domain) - 2, output); /* * Arithmetic coding of values */ { unsigned range; model_t *elements = alloc_model (M + 1, 0, 0, count); arith_t *encoder = alloc_encoder (output); for (range = 0; range < rs.range_no; range++) if (!rs.range_subdivided [range]) { state = rs.range_state [range]; label = rs.range_label [range]; for (edge = 0; isedge (wfa->into [state][label][edge]); edge++) ; encode_symbol (edge, encoder, elements); } free_encoder (encoder); free_model (elements); } debug_message ("delta-#edges: %5d bits. (%5d symbols => %5.2f bps)", bits_processed (output) - bits, edges, edges > 0 ? ((bits_processed (output) - bits) / (double) edges) : 0); } /* * Write matrix elements */ { unsigned bits = bits_processed (output); u_word_t *mapping1 = Calloc (wfa->states, sizeof (u_word_t)); u_word_t *mapping2 = Calloc (wfa->states, sizeof (u_word_t)); unsigned range; put_bit (output, use_normal_domains); put_bit (output, use_delta_domains); /* * Generate array of states which are admitted domains. * When coding intra frames 'mapping1' == 'mapping2' otherwise * 'mapping1' is a list of 'normal' domains which are admitted for * coding intra blocks * 'mapping2' is a list of 'delta' domains which are admitted for * coding the motion compensated prediction error */ { unsigned n1, n2, state; for (n1 = n2 = state = 0; state < wfa->states; state++) { mapping1 [state] = n1; if (usedomain (state, wfa) && (state < wfa->basis_states || use_delta_domains || !wfa->delta_state [state])) n1++; mapping2 [state] = n2; if (usedomain (state, wfa) && (state < wfa->basis_states || use_normal_domains || wfa->delta_state [state])) n2++; } debug_message ("# normal states = %d, # delta states = %d," " # WFA states = %d", n1, n2, wfa->states); } for (range = 0; range < rs.range_no; range++) if (!rs.range_subdivided [range]) { unsigned state = rs.range_state [range]; unsigned label = rs.range_label [range]; unsigned last = 1; u_word_t *mapping; unsigned max_value; unsigned edge; word_t domain; if (wfa->delta_state [state] || wfa->mv_tree [state][label].type != NONE) mapping = mapping2; else mapping = mapping1; max_value = mapping [rs.range_max_domain [range]]; for (edge = 0; isedge (domain = wfa->into [state][label][edge]); edge++) if (domain > 0) { total++; if (max_value - last) { write_bin_code (mapping [domain] - last, max_value - last, output); last = mapping [domain] + 1; } } } debug_message ("delta-index: %5d bits. (%5d symbols => %5.2f bps)", bits_processed (output) - bits, total, total > 0 ? ((bits_processed (output) - bits) / (double) total) : 0); Free (mapping1); Free (mapping2); } Free (rs.range_state); Free (rs.range_label); Free (rs.range_max_domain); Free (rs.range_subdivided); return total; } static unsigned column_0_encoding (const wfa_t *wfa, unsigned last_row, bitfile_t *output) /* * Write column 0 of the transition matrices of the 'wfa' to stream 'output' * with quasi arithmetic coding. * All rows from 'wfa->basis_states' up to 'last_row' are decoded. * * Return value: * number of non-zero matrix elements (WFA edges) */ { u_word_t high; /* Start of the current code range */ u_word_t low; /* End of the current code range */ unsigned *prob; /* probability array */ unsigned row; /* current matrix row */ unsigned label; /* current matrix label */ unsigned underflow; /* Underflow bits */ unsigned index; /* probability index */ unsigned total = 0; /* Number of '1' elements */ unsigned bits = bits_processed (output); /* * Compute the probability array: * prob[] = { 1/2, 1/2, 1/4, 1/4, 1/4, 1/4, * 1/8, ... , 1/16, ..., 1/(MAXPROB+1)} */ { unsigned n; unsigned exp; /* current exponent */ prob = Calloc (1 << (MAX_PROB + 1), sizeof (unsigned)); for (index = 0, n = MIN_PROB; n <= MAX_PROB; n++) for (exp = 0; exp < 1U << n; exp++, index++) prob [index] = n; } high = HIGH; /* 1.0 */ low = LOW; /* 0.0 */ underflow = 0; /* no underflow bits */ index = 0; /* * Encode column 0 with a quasi arithmetic coder (QAC). * Advantage of this QAC with respect to a binary AC: * Instead of using time consuming multiplications and divisions * to compute the probability of the most probable symbol (MPS) and * the range of the interval, a table look up procedure linked * with a shift operation is used for both computations. */ for (row = wfa->basis_states; row <= last_row; row++) for (label = 0; label < MAXLABELS; label++) if (isrange (wfa->tree [row][label])) { if (wfa->into [row][label][0] != 0) { /* * encode the MPS '0' */ high = high - ((high - low) >> prob [index]) - 1; RESCALE_OUTPUT_INTERVAL; if (index < 1020) index++; } else { /* * encode the LPS '1' */ low = high - ((high - low) >> prob [index]); RESCALE_OUTPUT_INTERVAL; total++; index >>= 1; } } /* * Flush the quasi-arithmetic encoder */ low = high; RESCALE_OUTPUT_INTERVAL; OUTPUT_BYTE_ALIGN (output); Free (prob); debug_message ("delta-state0: %5d bits. (%5d symbols => %5.2f bps)", bits_processed (output) - bits, total, total > 0 ? ((bits_processed (output) - bits) / (double) total) : 0); return total; } static unsigned chroma_encoding (const wfa_t *wfa, bitfile_t *output) /* * Write transition matrices of 'wfa' states which are part of the * chroma channels Cb and Cr to stream 'output'. * * Return value: * number of non-zero matrix elements (WFA edges) */ { unsigned domain; /* current domain, counter */ unsigned label; /* current label */ unsigned total = 0; /* number of '1' elements */ u_word_t high; /* Start of the current code range */ u_word_t low; /* End of the current code range */ unsigned underflow; /* underflow bits */ unsigned *prob; /* probability array */ unsigned index; /* probability index, counter */ unsigned next_index; /* probability of last domain */ unsigned row; /* current matrix row */ word_t *y_domains; unsigned count = 0; /* number of transitions for part 1 */ unsigned bits = bits_processed (output); /* * Compute the asymmetric probability array * prob[] = { 1/2, 1/2, 1/4, 1/4, 1/4, 1/4, * 1/8, ... , 1/16, ..., 1/(MAXPROB+1)} */ { unsigned n; unsigned exp; /* current exponent */ prob = Calloc (1 << (MAX_PROB + 1), sizeof (unsigned)); for (index = 0, n = MIN_PROB; n <= MAX_PROB; n++) for (exp = 0; exp < 1U << n; exp++, index++) prob [index] = n; } high = HIGH; /* 1.0 */ low = LOW; /* 0.0 */ underflow = 0; /* no underflow bits */ next_index = index = 0; y_domains = compute_hits (wfa->basis_states, wfa->tree [wfa->tree [wfa->root_state][0]][0], wfa->wfainfo->chroma_max_states, wfa); /* * First of all, read all matrix columns given in the list 'y_domains' * which note all admitted domains. * These matrix elements are stored with QAC (see column_0_encoding ()). */ for (domain = 0; y_domains [domain] != -1; domain++) { bool_t save_index = YES; /* YES: store current prob. index */ row = wfa->tree [wfa->tree [wfa->root_state][0]][0] + 1; index = next_index; for (; row < wfa->states; row++) { for (label = 0; label < MAXLABELS; label++) if (isrange (wfa->tree [row][label])) { unsigned edge; int into; bool_t match; /* approx with current domain found */ for (match = NO, edge = 0; isedge (into = wfa->into [row][label][edge]) && (unsigned) into < row; edge++) if (into == y_domains [domain] && into != wfa->y_state [row][label]) match = YES; if (!match) { /* * encode the MPS '0' */ high = high - ((high - low) >> prob [index]) - 1; RESCALE_OUTPUT_INTERVAL; if (index < 1020) index++; } else { /* * encode the LPS '1' */ low = high - ((high - low) >> prob [index]); RESCALE_OUTPUT_INTERVAL; total++; index >>= 1; } } if (save_index) { next_index = index; save_index = NO; } } } debug_message ("CbCr_matrix: %5d bits. (%5d symbols => %5.2f bps)", bits_processed (output) - bits, total, total > 0 ? ((bits_processed (output) - bits) / (double) total) : 0); count = total; bits = bits_processed (output); /* * Encode the additional column which indicates whether there * are transitions to a state with same spatial coordinates * in the Y component. * * Again, quasi arithmetic coding is used for this task. */ next_index = index = 0; for (row = wfa->tree [wfa->tree [wfa->root_state][0]][0] + 1; row < wfa->states; row++) for (label = 0; label < MAXLABELS; label++) if (!wfa->y_column [row][label]) { /* * encode the MPS '0' */ high = high - ((high - low) >> prob [index]) - 1; RESCALE_OUTPUT_INTERVAL; if (index < 1020) index++; } else { /* * encode the LPS '1' */ low = high - ((high - low) >> prob [index]); RESCALE_OUTPUT_INTERVAL; index >>= 1; total++; } /* * Flush the quasi-arithmetic encoder */ low = high; RESCALE_OUTPUT_INTERVAL; OUTPUT_BYTE_ALIGN (output); debug_message ("Yreferences: %5d bits. (%5d symbols => %5.2f bps)", bits_processed (output) - bits, total - count, total - count > 0 ? ((bits_processed (output) - bits) / (double) (total - count)) : 0); Free (prob); Free (y_domains); return total; }