aee.c 17.1 KB
Newer Older
1 2 3 4 5
/* Adaptive Entropy Encoder            */
/* CCSDS 121.0-B-1 and CCSDS 120.0-G-2 */

#include <stdio.h>
#include <stdlib.h>
6
#include <unistd.h>
7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30
#include <inttypes.h>
#include <string.h>

#include "libae.h"

#define ROS 5

#define MIN(a, b) (((a) < (b))? (a): (b))

enum
{
    M_NEW_BLOCK,
    M_GET_BLOCK,
    M_CHECK_ZERO_BLOCK,
    M_SELECT_CODE_OPTION,
    M_ENCODE_SPLIT,
    M_FLUSH_BLOCK,
    M_FLUSH_BLOCK_LOOP,
    M_ENCODE_UNCOMP,
    M_ENCODE_SE,
    M_ENCODE_ZERO,
};

typedef struct internal_state {
31 32 33
    int id_len;             /* bit length of code option identification key */
    int64_t last_in;        /* previous input for preprocessing */
    int64_t (*get_sample)(ae_streamp);
34 35 36
    int64_t xmin;           /* minimum integer for preprocessing */
    int64_t xmax;           /* maximum integer for preprocessing */
    int mode;               /* current mode of FSM */
37 38
    int i;                  /* counter for samples */
    int64_t *block_in;      /* input block buffer */
39 40
    uint8_t *block_out;     /* output block buffer */
    uint8_t *bp_out;        /* pointer to current output */
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
41
    size_t total_blocks;
42 43
    int bitp;               /* bit pointer to the next unused bit in accumulator */
    int block_deferred;     /* there is a block in the input buffer
44
                               but we first have to emit a zero block */
45 46 47 48 49 50
    int ref;                /* length of reference sample in current block
                               i.e. 0 or 1 depending on whether the block has
                               a reference sample or not */
    int zero_ref;           /* current zero block has a reference sample */
    int64_t zero_ref_sample;/* reference sample of zero block */
    int zero_blocks;        /* number of contiguous zero blocks */
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
51 52 53
#ifdef PROFILE
    int *prof;
#endif
54 55
} encode_state;

Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
56 57 58
static int64_t get_lsb_32(ae_streamp strm)
{
    int64_t data;
59

60 61 62 63
    data = (strm->next_in[3] << 24)
        | (strm->next_in[2] << 16)
        | (strm->next_in[1] << 8)
        | strm->next_in[0];
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
64

65
    strm->next_in += 4;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
66
    strm->total_in += 4;
67
    strm->avail_in -= 4;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
68 69 70 71 72 73 74
    return data;
}

static int64_t get_lsb_16(ae_streamp strm)
{
    int64_t data;

75
    data = (strm->next_in[1] << 8) | strm->next_in[0];
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
76

77
    strm->next_in += 2;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
78
    strm->total_in += 2;
79
    strm->avail_in -= 2;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
80 81 82 83 84 85
    return data;
}

static int64_t get_msb_32(ae_streamp strm)
{
    int64_t data;
86

87 88 89 90
    data = (strm->next_in[0] << 24)
        | (strm->next_in[1] << 16)
        | (strm->next_in[2] << 8)
        | strm->next_in[3];
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
91

92
    strm->next_in += 4;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
93
    strm->total_in += 4;
94
    strm->avail_in -= 4;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
95 96 97 98 99 100 101
    return data;
}

static int64_t get_msb_16(ae_streamp strm)
{
    int64_t data;

102
    data = (strm->next_in[0] << 8) | strm->next_in[1];
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
103

104
    strm->next_in += 2;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
105
    strm->total_in += 2;
106
    strm->avail_in -= 2;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
107 108 109 110 111 112 113 114 115
    return data;
}

static int64_t get_8(ae_streamp strm)
{
    strm->avail_in--;
    strm->total_in++;
    return *strm->next_in++;
}
116

117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135
int ae_encode_init(ae_streamp strm)
{
    int blklen;
    encode_state *state;

    /* Some sanity checks */
    if (strm->bit_per_sample > 32 || strm->bit_per_sample == 0)
    {
        return AE_ERRNO;
    }

    /* Internal state for encoder */
    state = (encode_state *) malloc(sizeof(encode_state));
    if (state == NULL)
    {
        return AE_MEM_ERROR;
    }
    strm->state = state;

136 137
    if (strm->bit_per_sample > 16)
    {
138
        state->id_len = 5;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
139 140
        if (strm->flags & AE_DATA_MSB)
            state->get_sample = get_msb_32;
141
        else
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
142
            state->get_sample = get_lsb_32;
143 144 145
    }
    else if (strm->bit_per_sample > 8)
    {
146
        state->id_len = 4;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
147 148
        if (strm->flags & AE_DATA_MSB)
            state->get_sample = get_msb_16;
149
        else
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
150
            state->get_sample = get_lsb_16;
151
    }
152
    else
153
    {
154
        state->id_len = 3;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
155
        state->get_sample = get_8;
156 157 158 159 160 161 162 163 164 165 166 167
    }

    if (strm->flags & AE_DATA_SIGNED)
    {
        state->xmin = -(1ULL << (strm->bit_per_sample - 1));
        state->xmax = (1ULL << (strm->bit_per_sample - 1)) - 1;
    }
    else
    {
        state->xmin = 0;
        state->xmax = (1ULL << strm->bit_per_sample) - 1;
    }
168

Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
169 170 171 172 173 174 175 176 177
#ifdef PROFILE
    state->prof = (int *)malloc((strm->bit_per_sample + 2) * sizeof(int));
    if (state->prof == NULL)
    {
        return AE_MEM_ERROR;
    }
    memset(state->prof, 0, (strm->bit_per_sample + 2) * sizeof(int));
#endif

178
    state->block_in = (int64_t *)malloc(strm->block_size * sizeof(int64_t));
179 180 181 182 183 184
    if (state->block_in == NULL)
    {
        return AE_MEM_ERROR;
    }

    blklen = (strm->block_size * strm->bit_per_sample
185
              + state->id_len) / 8 + 16;
186 187 188 189 190 191 192 193 194 195 196 197 198 199

    state->block_out = (uint8_t *)malloc(blklen);
    if (state->block_out == NULL)
    {
        return AE_MEM_ERROR;
    }
    state->bp_out = state->block_out;
    state->bitp = 8;

    strm->total_in = 0;
    strm->total_out = 0;

    state->mode = M_NEW_BLOCK;

Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
200
    state->total_blocks = 0;
201
    state->block_deferred = 0;
202
    state->zero_blocks = 0;
203 204 205 206 207 208
    state->zero_ref = 0;
    state->ref = 0;

    return AE_OK;
}

209 210 211 212 213 214 215 216 217 218 219 220 221 222
int ae_encode_end(ae_streamp strm)
{
    encode_state *state;

    state = strm->state;
#ifdef PROFILE
    free(state->prof);
#endif
    free(state->block_in);
    free(state->block_out);
    free(state);
    return AE_OK;
}

223
static inline void emit(encode_state *state, int64_t data, int bits)
224 225 226
{
    while(bits)
    {
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
227
        data &= ((1ULL << bits) - 1);
228 229 230
        if (bits <= state->bitp)
        {
            state->bitp -= bits;
231
            *state->bp_out += data << state->bitp;
232 233 234 235 236
            bits = 0;
        }
        else
        {
            bits -= state->bitp;
237
            *state->bp_out += data >> bits;
238 239 240 241 242 243
            *++state->bp_out = 0;
            state->bitp = 8;
        }
    }
}

244
static inline void emitfs(encode_state *state, int fs)
245
{
246
    emit(state, 1, fs + 1);
247
}
248

Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273
#ifdef PROFILE
static inline void profile_print(ae_streamp strm)
{
    int i, total;
    encode_state *state;

    state = strm->state;
    fprintf(stderr, "Blocks encoded by each coding option\n");
    fprintf(stderr, "Zero blocks:  %i\n", state->prof[0]);
    total = state->prof[0];
    fprintf(stderr, "Second Ext.:  %i\n", state->prof[strm->bit_per_sample+1]);
    total += state->prof[strm->bit_per_sample+1];
    fprintf(stderr, "FS:           %i\n", state->prof[1]);
    total += state->prof[1];
    for (i = 2; i < strm->bit_per_sample - 1; i++)
    {
        fprintf(stderr, "k = %02i:       %i\n", i-1, state->prof[i]);
        total += state->prof[i];
    }
    fprintf(stderr, "Uncompressed: %i\n", state->prof[strm->bit_per_sample]);
    total += state->prof[strm->bit_per_sample];
    fprintf(stderr, "Total blocks: %i\n", total);
}
#endif

274 275 276 277 278 279 280
int ae_encode(ae_streamp strm, int flush)
{
    /**
       Finite-state machine implementation of the adaptive entropy
       encoder.
    */

281 282 283
    int i, j, k, zb, this_bs;
    int64_t split_len;
    int64_t split_len_min, se_len, fs_len;
284
    int64_t d;
285
    int64_t theta, Delta;
286
    size_t avail_out, total_out;
287 288 289 290

    encode_state *state;

    state = strm->state;
291 292
    total_out = strm->total_out;
    avail_out = strm->avail_out;
293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324

    for (;;)
    {
        switch(state->mode)
        {
        case M_NEW_BLOCK:
            if (state->zero_blocks == 0)
            {
                /* copy leftover from last block */
                *state->block_out = *state->bp_out;
                state->bp_out = state->block_out;
            }

            if(state->block_deferred)
            {
                state->block_deferred = 0;
                state->mode = M_SELECT_CODE_OPTION;
                break;
            }

            state->i = 0;
            state->mode = M_GET_BLOCK;

        case M_GET_BLOCK:
            do
            {
                if (strm->avail_in == 0)
                {
                    if (flush == AE_FLUSH)
                    {
                        if (state->i > 0)
                        {
325
                            /* pad block with last sample if we have
326
                               a partial block */
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
327
                            state->block_in[state->i] = state->block_in[state->i - 1];
328 329 330 331 332 333 334
                        }
                        else
                        {
                            /* Pad last output byte with 1 bits
                               if user wants to flush, i.e. we got
                               all input there is.
                            */
335
                            emit(state, 0xff, state->bitp);
336
                            *strm->next_out++ = *state->bp_out;
337 338
                            avail_out--;
                            total_out++;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
339 340 341 342
#ifdef PROFILE
                            profile_print(strm);
#endif
                            goto req_buffer;
343 344
                        }
                    }
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
345 346 347 348
                    else
                    {
                        goto req_buffer;
                    }
349 350 351
                }
                else
                {
352
                    state->block_in[state->i] = state->get_sample(strm);
353 354 355 356
                }
            }
            while (++state->i < strm->block_size);

357 358
            state->total_blocks++;

359
            /* preprocess block if needed */
360
            if (strm->flags & AE_DATA_PREPROCESS)
361 362 363 364
            {
                /* If this is the first block in a segment
                   then we need to insert a reference sample.
                */
365
                if(state->total_blocks % strm->segment_size == 1)
366 367 368 369 370 371 372 373 374 375 376 377 378
                {
                    state->ref = 1;
                    state->last_in = state->block_in[0];
                }
                else
                {
                    state->ref = 0;
                }

                for (i = state->ref; i < strm->block_size; i++)
                {
                    theta = MIN(state->last_in - state->xmin,
                                state->xmax - state->last_in);
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
379
                    Delta = state->block_in[i] - state->last_in;
380
                    state->last_in = state->block_in[i];
381
                    if (0 <= Delta && Delta <= theta)
382 383 384
                    {
                        state->block_in[i] = 2 * Delta;
                    }
385
                    else if (-theta <= Delta && Delta < 0)
386 387 388 389
                    {
                        d = Delta < 0 ? -(uint64_t)Delta : Delta;
                        state->block_in[i] = 2 * d - 1;
                    }
390
                    else
391 392 393 394
                    {
                        state->block_in[i] = theta +
                            (Delta < 0 ? -(uint64_t)Delta : Delta);
                    }
395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414
                }
            }
            state->mode = M_CHECK_ZERO_BLOCK;

        case M_CHECK_ZERO_BLOCK:
            zb = 1;
            for (i = state->ref; i < strm->block_size && zb; i++)
                if (state->block_in[i] != 0) zb = 0;

            if (zb)
            {
                /* remember ref on first zero block */
                if (state->zero_blocks == 0)
                {
                    state->zero_ref = state->ref;
                    state->zero_ref_sample = state->block_in[0];
                }

                state->zero_blocks++;

415
                if (state->total_blocks % strm->segment_size == 0)
416 417 418
                {
                    if (state->zero_blocks > ROS)
                        state->zero_blocks = ROS;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
419 420 421
#ifdef PROFILE
                    state->prof[0] += state->zero_blocks;
#endif
422 423 424 425 426 427 428 429
                    state->mode = M_ENCODE_ZERO;
                    break;
                }
                state->mode = M_NEW_BLOCK;
                break;
            }
            else if (state->zero_blocks)
            {
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
430 431 432
#ifdef PROFILE
                state->prof[0] += state->zero_blocks;
#endif
433 434 435 436 437 438 439 440 441 442 443 444
                state->mode = M_ENCODE_ZERO;
                /* The current block isn't zero but we have to
                   emit a previous zero block first. The
                   current block has to be handled later.
                */
                state->block_deferred = 1;
                break;
            }

            state->mode = M_SELECT_CODE_OPTION;

        case M_SELECT_CODE_OPTION:
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
445 446 447
            /* If zero block isn't an option then count length of
               sample splitting options */

448
            /* Baseline is the size of an uncompressed block */
449
            split_len_min = (strm->block_size - state->ref) * strm->bit_per_sample;
450 451
            k = strm->bit_per_sample;

452 453 454 455
            /* Length of this block minus reference sample if present */
            this_bs = strm->block_size - state->ref;

            /* Add FS encoded to unencoded parts */
456 457
            for (j = 0; j < strm->bit_per_sample - 2; j++)
            {
458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474
#ifdef UNROLL_BLOCK_8
                fs_len = (state->block_in[1] >> j)
                    + (state->block_in[2] >> j)
                    + (state->block_in[3] >> j)
                    + (state->block_in[4] >> j)
                    + (state->block_in[5] >> j)
                    + (state->block_in[6] >> j)
                    + (state->block_in[7] >> j);
                if (state->ref == 0)
                    fs_len += (state->block_in[0] >> j);
#else
                fs_len = 0;
                for (i = state->ref; i < strm->block_size; i++)
                    fs_len += state->block_in[i] >> j;
#endif
                split_len = fs_len + this_bs * (j + 1);
                if (split_len < split_len_min)
475
                {
476
                    split_len_min = split_len;
477
                    k = j;
478 479 480 481 482 483 484 485

                    if (fs_len < this_bs)
                    {
                        /* Next can't get better because what we lose
                           by additional uncompressed bits isn't
                           compensated by a smaller FS part. */
                        break;
                    }
486
                }
487 488
                else
                    break;
489 490 491
            }

            /* Count bits for 2nd extension */
492 493
            se_len = 1;
            for (i = 0; i < strm->block_size && split_len_min > se_len; i+= 2)
494 495
            {
                d = state->block_in[i] + state->block_in[i + 1];
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
496
                /* we have to worry about overflow here */
497
                if (d > split_len_min)
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
498 499 500
                    se_len = d;
                else
                    se_len += d * (d + 1) / 2 + state->block_in[i + 1];
501 502 503
            }

            /* Decide which option to use */
504
            if (split_len_min <= se_len)
505 506 507
            {
                if (k == strm->bit_per_sample)
                {
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
508 509 510
#ifdef PROFILE
                    state->prof[k]++;
#endif
511 512 513 514 515
                    state->mode = M_ENCODE_UNCOMP;
                    break;
                }
                else
                {
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
516 517 518
#ifdef PROFILE
                    state->prof[k + 1]++;
#endif
519 520 521 522 523
                    state->mode = M_ENCODE_SPLIT;
                }
            }
            else
            {
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
524 525 526
#ifdef PROFILE
                state->prof[strm->bit_per_sample + 1]++;
#endif
527 528 529 530
                state->mode = M_ENCODE_SE;
                break;
            }

531
            emit(state, k + 1, state->id_len);
532
            if (state->ref)
533
                emit(state, state->block_in[0], strm->bit_per_sample);
534 535

            for (i = state->ref; i < strm->block_size; i++)
536
                emitfs(state, state->block_in[i] >> k);
537 538

            for (i = state->ref; i < strm->block_size; i++)
539
                emit(state, state->block_in[i], k);
540 541 542 543 544 545 546

            state->mode = M_FLUSH_BLOCK;

        case M_FLUSH_BLOCK:
            if (strm->avail_in == 0 && flush == AE_FLUSH)
            {
                /* pad last byte with 1 bits */
547
                emit(state, 0xff, state->bitp);
548 549 550 551 552 553 554
            }
            state->i = 0;
            state->mode = M_FLUSH_BLOCK_LOOP;

        case M_FLUSH_BLOCK_LOOP:
            while(state->block_out + state->i < state->bp_out)
            {
555
                if (avail_out == 0)
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
556 557 558 559
                {
#ifdef PROFILE
                    profile_print(strm);
#endif
560
                    goto req_buffer;
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
561
                }
562

563
                *strm->next_out++ = state->block_out[state->i];
564 565
                avail_out--;
                total_out++;
566 567 568 569 570 571
                state->i++;
            }
            state->mode = M_NEW_BLOCK;
            break;

        case M_ENCODE_UNCOMP:
572
            emit(state, 0x1f, state->id_len);
573
            for (i = 0; i < strm->block_size; i++)
574
                emit(state, state->block_in[i], strm->bit_per_sample);
575 576 577 578 579

            state->mode = M_FLUSH_BLOCK;
            break;

        case M_ENCODE_SE:
Mathis Rosenhauer's avatar
Mathis Rosenhauer committed
580
            emit(state, 1, state->id_len + 1);
581
            if (state->ref)
582
                emit(state, state->block_in[0], strm->bit_per_sample);
583 584 585 586

            for (i = 0; i < strm->block_size; i+= 2)
            {
                d = state->block_in[i] + state->block_in[i + 1];
587
                emitfs(state, d * (d + 1) / 2 + state->block_in[i + 1]);
588 589 590 591 592 593
            }

            state->mode = M_FLUSH_BLOCK;
            break;

        case M_ENCODE_ZERO:
594
            emit(state, 0, state->id_len + 1);
595 596
            if (state->zero_ref)
            {
597
                emit(state, state->zero_ref_sample, strm->bit_per_sample);
598
            }
599
            emitfs(state, state->zero_blocks - 1);
600 601 602 603 604 605 606 607 608 609
            state->zero_blocks = 0;
            state->mode = M_FLUSH_BLOCK;
            break;

        default:
            return AE_STREAM_ERROR;
        }
    }

req_buffer:
610 611
    strm->total_out = total_out;
    strm->avail_out = avail_out;
612 613
    return AE_OK;
}