jpc_qmfb.c (30572B)
1 /* 2 * Copyright (c) 1999-2000 Image Power, Inc. and the University of 3 * British Columbia. 4 * Copyright (c) 2001-2003 Michael David Adams. 5 * All rights reserved. 6 */ 7 8 /* __START_OF_JASPER_LICENSE__ 9 * 10 * JasPer License Version 2.0 11 * 12 * Copyright (c) 1999-2000 Image Power, Inc. 13 * Copyright (c) 1999-2000 The University of British Columbia 14 * Copyright (c) 2001-2003 Michael David Adams 15 * 16 * All rights reserved. 17 * 18 * Permission is hereby granted, free of charge, to any person (the 19 * "User") obtaining a copy of this software and associated documentation 20 * files (the "Software"), to deal in the Software without restriction, 21 * including without limitation the rights to use, copy, modify, merge, 22 * publish, distribute, and/or sell copies of the Software, and to permit 23 * persons to whom the Software is furnished to do so, subject to the 24 * following conditions: 25 * 26 * 1. The above copyright notices and this permission notice (which 27 * includes the disclaimer below) shall be included in all copies or 28 * substantial portions of the Software. 29 * 30 * 2. The name of a copyright holder shall not be used to endorse or 31 * promote products derived from the Software without specific prior 32 * written permission. 33 * 34 * THIS DISCLAIMER OF WARRANTY CONSTITUTES AN ESSENTIAL PART OF THIS 35 * LICENSE. NO USE OF THE SOFTWARE IS AUTHORIZED HEREUNDER EXCEPT UNDER 36 * THIS DISCLAIMER. THE SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS 37 * "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING 38 * BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A 39 * PARTICULAR PURPOSE AND NONINFRINGEMENT OF THIRD PARTY RIGHTS. IN NO 40 * EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, OR ANY SPECIAL 41 * INDIRECT OR CONSEQUENTIAL DAMAGES, OR ANY DAMAGES WHATSOEVER RESULTING 42 * FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, 43 * NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION 44 * WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. NO ASSURANCES ARE 45 * PROVIDED BY THE COPYRIGHT HOLDERS THAT THE SOFTWARE DOES NOT INFRINGE 46 * THE PATENT OR OTHER INTELLECTUAL PROPERTY RIGHTS OF ANY OTHER ENTITY. 47 * EACH COPYRIGHT HOLDER DISCLAIMS ANY LIABILITY TO THE USER FOR CLAIMS 48 * BROUGHT BY ANY OTHER ENTITY BASED ON INFRINGEMENT OF INTELLECTUAL 49 * PROPERTY RIGHTS OR OTHERWISE. AS A CONDITION TO EXERCISING THE RIGHTS 50 * GRANTED HEREUNDER, EACH USER HEREBY ASSUMES SOLE RESPONSIBILITY TO SECURE 51 * ANY OTHER INTELLECTUAL PROPERTY RIGHTS NEEDED, IF ANY. THE SOFTWARE 52 * IS NOT FAULT-TOLERANT AND IS NOT INTENDED FOR USE IN MISSION-CRITICAL 53 * SYSTEMS, SUCH AS THOSE USED IN THE OPERATION OF NUCLEAR FACILITIES, 54 * AIRCRAFT NAVIGATION OR COMMUNICATION SYSTEMS, AIR TRAFFIC CONTROL 55 * SYSTEMS, DIRECT LIFE SUPPORT MACHINES, OR WEAPONS SYSTEMS, IN WHICH 56 * THE FAILURE OF THE SOFTWARE OR SYSTEM COULD LEAD DIRECTLY TO DEATH, 57 * PERSONAL INJURY, OR SEVERE PHYSICAL OR ENVIRONMENTAL DAMAGE ("HIGH 58 * RISK ACTIVITIES"). THE COPYRIGHT HOLDERS SPECIFICALLY DISCLAIM ANY 59 * EXPRESS OR IMPLIED WARRANTY OF FITNESS FOR HIGH RISK ACTIVITIES. 60 * 61 * __END_OF_JASPER_LICENSE__ 62 */ 63 64 /* 65 * Quadrature Mirror-Image Filter Bank (QMFB) Library 66 * 67 * $Id: jpc_qmfb.c 1918 2005-07-24 14:12:08Z baford $ 68 */ 69 70 /******************************************************************************\ 71 * Includes. 72 \******************************************************************************/ 73 74 #include <assert.h> 75 76 #include "jasper/jas_fix.h" 77 #include "jasper/jas_malloc.h" 78 #include "jasper/jas_math.h" 79 80 #include "jpc_qmfb.h" 81 #include "jpc_tsfb.h" 82 #include "jpc_math.h" 83 84 /******************************************************************************\ 85 * 86 \******************************************************************************/ 87 88 static jpc_qmfb1d_t *jpc_qmfb1d_create(void); 89 90 static int jpc_ft_getnumchans(jpc_qmfb1d_t *qmfb); 91 static int jpc_ft_getanalfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters); 92 static int jpc_ft_getsynfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters); 93 static void jpc_ft_analyze(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x); 94 static void jpc_ft_synthesize(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x); 95 96 static int jpc_ns_getnumchans(jpc_qmfb1d_t *qmfb); 97 static int jpc_ns_getanalfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters); 98 static int jpc_ns_getsynfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters); 99 static void jpc_ns_analyze(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x); 100 static void jpc_ns_synthesize(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x); 101 102 /******************************************************************************\ 103 * 104 \******************************************************************************/ 105 106 jpc_qmfb1dops_t jpc_ft_ops = { 107 jpc_ft_getnumchans, 108 jpc_ft_getanalfilters, 109 jpc_ft_getsynfilters, 110 jpc_ft_analyze, 111 jpc_ft_synthesize 112 }; 113 114 jpc_qmfb1dops_t jpc_ns_ops = { 115 jpc_ns_getnumchans, 116 jpc_ns_getanalfilters, 117 jpc_ns_getsynfilters, 118 jpc_ns_analyze, 119 jpc_ns_synthesize 120 }; 121 122 /******************************************************************************\ 123 * 124 \******************************************************************************/ 125 126 static void jpc_qmfb1d_setup(jpc_fix_t *startptr, int startind, int endind, 127 int intrastep, jpc_fix_t **lstartptr, int *lstartind, int *lendind, 128 jpc_fix_t **hstartptr, int *hstartind, int *hendind) 129 { 130 *lstartind = JPC_CEILDIVPOW2(startind, 1); 131 *lendind = JPC_CEILDIVPOW2(endind, 1); 132 *hstartind = JPC_FLOORDIVPOW2(startind, 1); 133 *hendind = JPC_FLOORDIVPOW2(endind, 1); 134 *lstartptr = startptr; 135 *hstartptr = &startptr[(*lendind - *lstartind) * intrastep]; 136 } 137 138 static void jpc_qmfb1d_split(jpc_fix_t *startptr, int startind, int endind, 139 register int step, jpc_fix_t *lstartptr, int lstartind, int lendind, 140 jpc_fix_t *hstartptr, int hstartind, int hendind) 141 { 142 int bufsize = JPC_CEILDIVPOW2(endind - startind, 2); 143 #if !defined(HAVE_VLA) 144 #define QMFB_SPLITBUFSIZE 4096 145 jpc_fix_t splitbuf[QMFB_SPLITBUFSIZE]; 146 #else 147 jpc_fix_t splitbuf[bufsize]; 148 #endif 149 jpc_fix_t *buf = splitbuf; 150 int llen; 151 int hlen; 152 int twostep; 153 jpc_fix_t *tmpptr; 154 register jpc_fix_t *ptr; 155 register jpc_fix_t *hptr; 156 register jpc_fix_t *lptr; 157 register int n; 158 int state; 159 160 twostep = step << 1; 161 llen = lendind - lstartind; 162 hlen = hendind - hstartind; 163 164 #if !defined(HAVE_VLA) 165 /* Get a buffer. */ 166 if (bufsize > QMFB_SPLITBUFSIZE) { 167 if (!(buf = jas_malloc(bufsize * sizeof(jpc_fix_t)))) { 168 /* We have no choice but to commit suicide in this case. */ 169 abort(); 170 } 171 } 172 #endif 173 174 if (hstartind < lstartind) { 175 /* The first sample in the input signal is to appear 176 in the highpass subband signal. */ 177 /* Copy the appropriate samples into the lowpass subband 178 signal, saving any samples destined for the highpass subband 179 signal as they are overwritten. */ 180 tmpptr = buf; 181 ptr = &startptr[step]; 182 lptr = lstartptr; 183 n = llen; 184 state = 1; 185 while (n-- > 0) { 186 if (state) { 187 *tmpptr = *lptr; 188 ++tmpptr; 189 } 190 *lptr = *ptr; 191 ptr += twostep; 192 lptr += step; 193 state ^= 1; 194 } 195 /* Copy the appropriate samples into the highpass subband 196 signal. */ 197 /* Handle the nonoverwritten samples. */ 198 hptr = &hstartptr[(hlen - 1) * step]; 199 ptr = &startptr[(((llen + hlen - 1) >> 1) << 1) * step]; 200 n = hlen - (tmpptr - buf); 201 while (n-- > 0) { 202 *hptr = *ptr; 203 hptr -= step; 204 ptr -= twostep; 205 } 206 /* Handle the overwritten samples. */ 207 n = tmpptr - buf; 208 while (n-- > 0) { 209 --tmpptr; 210 *hptr = *tmpptr; 211 hptr -= step; 212 } 213 } else { 214 /* The first sample in the input signal is to appear 215 in the lowpass subband signal. */ 216 /* Copy the appropriate samples into the lowpass subband 217 signal, saving any samples for the highpass subband 218 signal as they are overwritten. */ 219 state = 0; 220 ptr = startptr; 221 lptr = lstartptr; 222 tmpptr = buf; 223 n = llen; 224 while (n-- > 0) { 225 if (state) { 226 *tmpptr = *lptr; 227 ++tmpptr; 228 } 229 *lptr = *ptr; 230 ptr += twostep; 231 lptr += step; 232 state ^= 1; 233 } 234 /* Copy the appropriate samples into the highpass subband 235 signal. */ 236 /* Handle the nonoverwritten samples. */ 237 ptr = &startptr[((((llen + hlen) >> 1) << 1) - 1) * step]; 238 hptr = &hstartptr[(hlen - 1) * step]; 239 n = hlen - (tmpptr - buf); 240 while (n-- > 0) { 241 *hptr = *ptr; 242 ptr -= twostep; 243 hptr -= step; 244 } 245 /* Handle the overwritten samples. */ 246 n = tmpptr - buf; 247 while (n-- > 0) { 248 --tmpptr; 249 *hptr = *tmpptr; 250 hptr -= step; 251 } 252 } 253 254 #if !defined(HAVE_VLA) 255 /* If the split buffer was allocated on the heap, free this memory. */ 256 if (buf != splitbuf) { 257 jas_free(buf); 258 } 259 #endif 260 } 261 262 static void jpc_qmfb1d_join(jpc_fix_t *startptr, int startind, int endind, 263 register int step, jpc_fix_t *lstartptr, int lstartind, int lendind, 264 jpc_fix_t *hstartptr, int hstartind, int hendind) 265 { 266 int bufsize = JPC_CEILDIVPOW2(endind - startind, 2); 267 #if !defined(HAVE_VLA) 268 #define QMFB_JOINBUFSIZE 4096 269 jpc_fix_t joinbuf[QMFB_JOINBUFSIZE]; 270 #else 271 jpc_fix_t joinbuf[bufsize]; 272 #endif 273 jpc_fix_t *buf = joinbuf; 274 int llen; 275 int hlen; 276 int twostep; 277 jpc_fix_t *tmpptr; 278 register jpc_fix_t *ptr; 279 register jpc_fix_t *hptr; 280 register jpc_fix_t *lptr; 281 register int n; 282 int state; 283 284 #if !defined(HAVE_VLA) 285 /* Allocate memory for the join buffer from the heap. */ 286 if (bufsize > QMFB_JOINBUFSIZE) { 287 if (!(buf = jas_malloc(bufsize * sizeof(jpc_fix_t)))) { 288 /* We have no choice but to commit suicide. */ 289 abort(); 290 } 291 } 292 #endif 293 294 twostep = step << 1; 295 llen = lendind - lstartind; 296 hlen = hendind - hstartind; 297 298 if (hstartind < lstartind) { 299 /* The first sample in the highpass subband signal is to 300 appear first in the output signal. */ 301 /* Copy the appropriate samples into the first phase of the 302 output signal. */ 303 tmpptr = buf; 304 hptr = hstartptr; 305 ptr = startptr; 306 n = (llen + 1) >> 1; 307 while (n-- > 0) { 308 *tmpptr = *ptr; 309 *ptr = *hptr; 310 ++tmpptr; 311 ptr += twostep; 312 hptr += step; 313 } 314 n = hlen - ((llen + 1) >> 1); 315 while (n-- > 0) { 316 *ptr = *hptr; 317 ptr += twostep; 318 hptr += step; 319 } 320 /* Copy the appropriate samples into the second phase of 321 the output signal. */ 322 ptr -= (lendind > hendind) ? (step) : (step + twostep); 323 state = !((llen - 1) & 1); 324 lptr = &lstartptr[(llen - 1) * step]; 325 n = llen; 326 while (n-- > 0) { 327 if (state) { 328 --tmpptr; 329 *ptr = *tmpptr; 330 } else { 331 *ptr = *lptr; 332 } 333 lptr -= step; 334 ptr -= twostep; 335 state ^= 1; 336 } 337 } else { 338 /* The first sample in the lowpass subband signal is to 339 appear first in the output signal. */ 340 /* Copy the appropriate samples into the first phase of the 341 output signal (corresponding to even indexed samples). */ 342 lptr = &lstartptr[(llen - 1) * step]; 343 ptr = &startptr[((llen - 1) << 1) * step]; 344 n = llen >> 1; 345 tmpptr = buf; 346 while (n-- > 0) { 347 *tmpptr = *ptr; 348 *ptr = *lptr; 349 ++tmpptr; 350 ptr -= twostep; 351 lptr -= step; 352 } 353 n = llen - (llen >> 1); 354 while (n-- > 0) { 355 *ptr = *lptr; 356 ptr -= twostep; 357 lptr -= step; 358 } 359 /* Copy the appropriate samples into the second phase of 360 the output signal (corresponding to odd indexed 361 samples). */ 362 ptr = &startptr[step]; 363 hptr = hstartptr; 364 state = !(llen & 1); 365 n = hlen; 366 while (n-- > 0) { 367 if (state) { 368 --tmpptr; 369 *ptr = *tmpptr; 370 } else { 371 *ptr = *hptr; 372 } 373 hptr += step; 374 ptr += twostep; 375 state ^= 1; 376 } 377 } 378 379 #if !defined(HAVE_VLA) 380 /* If the join buffer was allocated on the heap, free this memory. */ 381 if (buf != joinbuf) { 382 jas_free(buf); 383 } 384 #endif 385 } 386 387 /******************************************************************************\ 388 * Code for 5/3 transform. 389 \******************************************************************************/ 390 391 static int jpc_ft_getnumchans(jpc_qmfb1d_t *qmfb) 392 { 393 /* Avoid compiler warnings about unused parameters. */ 394 qmfb = 0; 395 396 return 2; 397 } 398 399 static int jpc_ft_getanalfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters) 400 { 401 /* Avoid compiler warnings about unused parameters. */ 402 qmfb = 0; 403 len = 0; 404 filters = 0; 405 abort(); 406 return -1; 407 } 408 409 static int jpc_ft_getsynfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters) 410 { 411 jas_seq_t *lf; 412 jas_seq_t *hf; 413 414 /* Avoid compiler warnings about unused parameters. */ 415 qmfb = 0; 416 417 lf = 0; 418 hf = 0; 419 420 if (len > 1 || (!len)) { 421 if (!(lf = jas_seq_create(-1, 2))) { 422 goto error; 423 } 424 jas_seq_set(lf, -1, jpc_dbltofix(0.5)); 425 jas_seq_set(lf, 0, jpc_dbltofix(1.0)); 426 jas_seq_set(lf, 1, jpc_dbltofix(0.5)); 427 if (!(hf = jas_seq_create(-1, 4))) { 428 goto error; 429 } 430 jas_seq_set(hf, -1, jpc_dbltofix(-0.125)); 431 jas_seq_set(hf, 0, jpc_dbltofix(-0.25)); 432 jas_seq_set(hf, 1, jpc_dbltofix(0.75)); 433 jas_seq_set(hf, 2, jpc_dbltofix(-0.25)); 434 jas_seq_set(hf, 3, jpc_dbltofix(-0.125)); 435 } else if (len == 1) { 436 if (!(lf = jas_seq_create(0, 1))) { 437 goto error; 438 } 439 jas_seq_set(lf, 0, jpc_dbltofix(1.0)); 440 if (!(hf = jas_seq_create(0, 1))) { 441 goto error; 442 } 443 jas_seq_set(hf, 0, jpc_dbltofix(2.0)); 444 } else { 445 abort(); 446 } 447 448 filters[0] = lf; 449 filters[1] = hf; 450 451 return 0; 452 453 error: 454 if (lf) { 455 jas_seq_destroy(lf); 456 } 457 if (hf) { 458 jas_seq_destroy(hf); 459 } 460 return -1; 461 } 462 463 #define NFT_LIFT0(lstartptr, lstartind, lendind, hstartptr, hstartind, hendind, step, pluseq) \ 464 { \ 465 register jpc_fix_t *lptr = (lstartptr); \ 466 register jpc_fix_t *hptr = (hstartptr); \ 467 register int n = (hendind) - (hstartind); \ 468 if ((hstartind) < (lstartind)) { \ 469 pluseq(*hptr, *lptr); \ 470 hptr += (step); \ 471 --n; \ 472 } \ 473 if ((hendind) >= (lendind)) { \ 474 --n; \ 475 } \ 476 while (n-- > 0) { \ 477 pluseq(*hptr, jpc_fix_asr(jpc_fix_add(*lptr, lptr[(step)]), 1)); \ 478 hptr += (step); \ 479 lptr += (step); \ 480 } \ 481 if ((hendind) >= (lendind)) { \ 482 pluseq(*hptr, *lptr); \ 483 } \ 484 } 485 486 #define NFT_LIFT1(lstartptr, lstartind, lendind, hstartptr, hstartind, hendind, step, pluseq) \ 487 { \ 488 register jpc_fix_t *lptr = (lstartptr); \ 489 register jpc_fix_t *hptr = (hstartptr); \ 490 register int n = (lendind) - (lstartind); \ 491 if ((hstartind) >= (lstartind)) { \ 492 pluseq(*lptr, *hptr); \ 493 lptr += (step); \ 494 --n; \ 495 } \ 496 if ((lendind) > (hendind)) { \ 497 --n; \ 498 } \ 499 while (n-- > 0) { \ 500 pluseq(*lptr, jpc_fix_asr(jpc_fix_add(*hptr, hptr[(step)]), 2)); \ 501 lptr += (step); \ 502 hptr += (step); \ 503 } \ 504 if ((lendind) > (hendind)) { \ 505 pluseq(*lptr, *hptr); \ 506 } \ 507 } 508 509 #define RFT_LIFT0(lstartptr, lstartind, lendind, hstartptr, hstartind, hendind, step, pmeqop) \ 510 { \ 511 register jpc_fix_t *lptr = (lstartptr); \ 512 register jpc_fix_t *hptr = (hstartptr); \ 513 register int n = (hendind) - (hstartind); \ 514 if ((hstartind) < (lstartind)) { \ 515 *hptr pmeqop *lptr; \ 516 hptr += (step); \ 517 --n; \ 518 } \ 519 if ((hendind) >= (lendind)) { \ 520 --n; \ 521 } \ 522 while (n-- > 0) { \ 523 *hptr pmeqop (*lptr + lptr[(step)]) >> 1; \ 524 hptr += (step); \ 525 lptr += (step); \ 526 } \ 527 if ((hendind) >= (lendind)) { \ 528 *hptr pmeqop *lptr; \ 529 } \ 530 } 531 532 #define RFT_LIFT1(lstartptr, lstartind, lendind, hstartptr, hstartind, hendind, step, pmeqop) \ 533 { \ 534 register jpc_fix_t *lptr = (lstartptr); \ 535 register jpc_fix_t *hptr = (hstartptr); \ 536 register int n = (lendind) - (lstartind); \ 537 if ((hstartind) >= (lstartind)) { \ 538 *lptr pmeqop ((*hptr << 1) + 2) >> 2; \ 539 lptr += (step); \ 540 --n; \ 541 } \ 542 if ((lendind) > (hendind)) { \ 543 --n; \ 544 } \ 545 while (n-- > 0) { \ 546 *lptr pmeqop ((*hptr + hptr[(step)]) + 2) >> 2; \ 547 lptr += (step); \ 548 hptr += (step); \ 549 } \ 550 if ((lendind) > (hendind)) { \ 551 *lptr pmeqop ((*hptr << 1) + 2) >> 2; \ 552 } \ 553 } 554 555 static void jpc_ft_analyze(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x) 556 { 557 jpc_fix_t *startptr; 558 int startind; 559 int endind; 560 jpc_fix_t * lstartptr; 561 int lstartind; 562 int lendind; 563 jpc_fix_t * hstartptr; 564 int hstartind; 565 int hendind; 566 int interstep; 567 int intrastep; 568 int numseq; 569 570 /* Avoid compiler warnings about unused parameters. */ 571 qmfb = 0; 572 573 if (flags & JPC_QMFB1D_VERT) { 574 interstep = 1; 575 intrastep = jas_seq2d_rowstep(x); 576 numseq = jas_seq2d_width(x); 577 startind = jas_seq2d_ystart(x); 578 endind = jas_seq2d_yend(x); 579 } else { 580 interstep = jas_seq2d_rowstep(x); 581 intrastep = 1; 582 numseq = jas_seq2d_height(x); 583 startind = jas_seq2d_xstart(x); 584 endind = jas_seq2d_xend(x); 585 } 586 587 assert(startind < endind); 588 589 startptr = jas_seq2d_getref(x, jas_seq2d_xstart(x), jas_seq2d_ystart(x)); 590 if (flags & JPC_QMFB1D_RITIMODE) { 591 while (numseq-- > 0) { 592 jpc_qmfb1d_setup(startptr, startind, endind, intrastep, 593 &lstartptr, &lstartind, &lendind, &hstartptr, 594 &hstartind, &hendind); 595 if (endind - startind > 1) { 596 jpc_qmfb1d_split(startptr, startind, endind, 597 intrastep, lstartptr, lstartind, lendind, 598 hstartptr, hstartind, hendind); 599 RFT_LIFT0(lstartptr, lstartind, lendind, 600 hstartptr, hstartind, hendind, intrastep, -=); 601 RFT_LIFT1(lstartptr, lstartind, lendind, 602 hstartptr, hstartind, hendind, intrastep, +=); 603 } else { 604 if (lstartind == lendind) { 605 *startptr <<= 1; 606 } 607 } 608 startptr += interstep; 609 } 610 } else { 611 while (numseq-- > 0) { 612 jpc_qmfb1d_setup(startptr, startind, endind, intrastep, 613 &lstartptr, &lstartind, &lendind, &hstartptr, 614 &hstartind, &hendind); 615 if (endind - startind > 1) { 616 jpc_qmfb1d_split(startptr, startind, endind, 617 intrastep, lstartptr, lstartind, lendind, 618 hstartptr, hstartind, hendind); 619 NFT_LIFT0(lstartptr, lstartind, lendind, 620 hstartptr, hstartind, hendind, intrastep, 621 jpc_fix_minuseq); 622 NFT_LIFT1(lstartptr, lstartind, lendind, 623 hstartptr, hstartind, hendind, intrastep, 624 jpc_fix_pluseq); 625 } else { 626 if (lstartind == lendind) { 627 *startptr = jpc_fix_asl(*startptr, 1); 628 } 629 } 630 startptr += interstep; 631 } 632 } 633 } 634 635 static void jpc_ft_synthesize(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x) 636 { 637 jpc_fix_t *startptr; 638 int startind; 639 int endind; 640 jpc_fix_t *lstartptr; 641 int lstartind; 642 int lendind; 643 jpc_fix_t *hstartptr; 644 int hstartind; 645 int hendind; 646 int interstep; 647 int intrastep; 648 int numseq; 649 650 /* Avoid compiler warnings about unused parameters. */ 651 qmfb = 0; 652 653 if (flags & JPC_QMFB1D_VERT) { 654 interstep = 1; 655 intrastep = jas_seq2d_rowstep(x); 656 numseq = jas_seq2d_width(x); 657 startind = jas_seq2d_ystart(x); 658 endind = jas_seq2d_yend(x); 659 } else { 660 interstep = jas_seq2d_rowstep(x); 661 intrastep = 1; 662 numseq = jas_seq2d_height(x); 663 startind = jas_seq2d_xstart(x); 664 endind = jas_seq2d_xend(x); 665 } 666 667 assert(startind < endind); 668 669 startptr = jas_seq2d_getref(x, jas_seq2d_xstart(x), jas_seq2d_ystart(x)); 670 if (flags & JPC_QMFB1D_RITIMODE) { 671 while (numseq-- > 0) { 672 jpc_qmfb1d_setup(startptr, startind, endind, intrastep, 673 &lstartptr, &lstartind, &lendind, &hstartptr, 674 &hstartind, &hendind); 675 if (endind - startind > 1) { 676 RFT_LIFT1(lstartptr, lstartind, lendind, 677 hstartptr, hstartind, hendind, intrastep, -=); 678 RFT_LIFT0(lstartptr, lstartind, lendind, 679 hstartptr, hstartind, hendind, intrastep, +=); 680 jpc_qmfb1d_join(startptr, startind, endind, 681 intrastep, lstartptr, lstartind, lendind, 682 hstartptr, hstartind, hendind); 683 } else { 684 if (lstartind == lendind) { 685 *startptr >>= 1; 686 } 687 } 688 startptr += interstep; 689 } 690 } else { 691 while (numseq-- > 0) { 692 jpc_qmfb1d_setup(startptr, startind, endind, intrastep, 693 &lstartptr, &lstartind, &lendind, &hstartptr, 694 &hstartind, &hendind); 695 if (endind - startind > 1) { 696 NFT_LIFT1(lstartptr, lstartind, lendind, 697 hstartptr, hstartind, hendind, intrastep, 698 jpc_fix_minuseq); 699 NFT_LIFT0(lstartptr, lstartind, lendind, 700 hstartptr, hstartind, hendind, intrastep, 701 jpc_fix_pluseq); 702 jpc_qmfb1d_join(startptr, startind, endind, 703 intrastep, lstartptr, lstartind, lendind, 704 hstartptr, hstartind, hendind); 705 } else { 706 if (lstartind == lendind) { 707 *startptr = jpc_fix_asr(*startptr, 1); 708 } 709 } 710 startptr += interstep; 711 } 712 } 713 } 714 715 /******************************************************************************\ 716 * Code for 9/7 transform. 717 \******************************************************************************/ 718 719 static int jpc_ns_getnumchans(jpc_qmfb1d_t *qmfb) 720 { 721 /* Avoid compiler warnings about unused parameters. */ 722 qmfb = 0; 723 724 return 2; 725 } 726 727 static int jpc_ns_getanalfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters) 728 { 729 /* Avoid compiler warnings about unused parameters. */ 730 qmfb = 0; 731 len = 0; 732 filters = 0; 733 734 abort(); 735 return -1; 736 } 737 738 static int jpc_ns_getsynfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters) 739 { 740 jas_seq_t *lf; 741 jas_seq_t *hf; 742 743 /* Avoid compiler warnings about unused parameters. */ 744 qmfb = 0; 745 746 lf = 0; 747 hf = 0; 748 749 if (len > 1 || (!len)) { 750 if (!(lf = jas_seq_create(-3, 4))) { 751 goto error; 752 } 753 jas_seq_set(lf, -3, jpc_dbltofix(-0.09127176311424948)); 754 jas_seq_set(lf, -2, jpc_dbltofix(-0.05754352622849957)); 755 jas_seq_set(lf, -1, jpc_dbltofix(0.5912717631142470)); 756 jas_seq_set(lf, 0, jpc_dbltofix(1.115087052456994)); 757 jas_seq_set(lf, 1, jpc_dbltofix(0.5912717631142470)); 758 jas_seq_set(lf, 2, jpc_dbltofix(-0.05754352622849957)); 759 jas_seq_set(lf, 3, jpc_dbltofix(-0.09127176311424948)); 760 if (!(hf = jas_seq_create(-3, 6))) { 761 goto error; 762 } 763 jas_seq_set(hf, -3, jpc_dbltofix(-0.02674875741080976 * 2.0)); 764 jas_seq_set(hf, -2, jpc_dbltofix(-0.01686411844287495 * 2.0)); 765 jas_seq_set(hf, -1, jpc_dbltofix(0.07822326652898785 * 2.0)); 766 jas_seq_set(hf, 0, jpc_dbltofix(0.2668641184428723 * 2.0)); 767 jas_seq_set(hf, 1, jpc_dbltofix(-0.6029490182363579 * 2.0)); 768 jas_seq_set(hf, 2, jpc_dbltofix(0.2668641184428723 * 2.0)); 769 jas_seq_set(hf, 3, jpc_dbltofix(0.07822326652898785 * 2.0)); 770 jas_seq_set(hf, 4, jpc_dbltofix(-0.01686411844287495 * 2.0)); 771 jas_seq_set(hf, 5, jpc_dbltofix(-0.02674875741080976 * 2.0)); 772 } else if (len == 1) { 773 if (!(lf = jas_seq_create(0, 1))) { 774 goto error; 775 } 776 jas_seq_set(lf, 0, jpc_dbltofix(1.0)); 777 if (!(hf = jas_seq_create(0, 1))) { 778 goto error; 779 } 780 jas_seq_set(hf, 0, jpc_dbltofix(2.0)); 781 } else { 782 abort(); 783 } 784 785 filters[0] = lf; 786 filters[1] = hf; 787 788 return 0; 789 790 error: 791 if (lf) { 792 jas_seq_destroy(lf); 793 } 794 if (hf) { 795 jas_seq_destroy(hf); 796 } 797 return -1; 798 } 799 800 #define NNS_LIFT0(lstartptr, lstartind, lendind, hstartptr, hstartind, hendind, step, alpha) \ 801 { \ 802 register jpc_fix_t *lptr = (lstartptr); \ 803 register jpc_fix_t *hptr = (hstartptr); \ 804 register int n = (hendind) - (hstartind); \ 805 jpc_fix_t twoalpha = jpc_fix_mulbyint(alpha, 2); \ 806 if ((hstartind) < (lstartind)) { \ 807 jpc_fix_pluseq(*hptr, jpc_fix_mul(*lptr, (twoalpha))); \ 808 hptr += (step); \ 809 --n; \ 810 } \ 811 if ((hendind) >= (lendind)) { \ 812 --n; \ 813 } \ 814 while (n-- > 0) { \ 815 jpc_fix_pluseq(*hptr, jpc_fix_mul(jpc_fix_add(*lptr, lptr[(step)]), (alpha))); \ 816 hptr += (step); \ 817 lptr += (step); \ 818 } \ 819 if ((hendind) >= (lendind)) { \ 820 jpc_fix_pluseq(*hptr, jpc_fix_mul(*lptr, (twoalpha))); \ 821 } \ 822 } 823 824 #define NNS_LIFT1(lstartptr, lstartind, lendind, hstartptr, hstartind, hendind, step, alpha) \ 825 { \ 826 register jpc_fix_t *lptr = (lstartptr); \ 827 register jpc_fix_t *hptr = (hstartptr); \ 828 register int n = (lendind) - (lstartind); \ 829 int twoalpha = jpc_fix_mulbyint(alpha, 2); \ 830 if ((hstartind) >= (lstartind)) { \ 831 jpc_fix_pluseq(*lptr, jpc_fix_mul(*hptr, (twoalpha))); \ 832 lptr += (step); \ 833 --n; \ 834 } \ 835 if ((lendind) > (hendind)) { \ 836 --n; \ 837 } \ 838 while (n-- > 0) { \ 839 jpc_fix_pluseq(*lptr, jpc_fix_mul(jpc_fix_add(*hptr, hptr[(step)]), (alpha))); \ 840 lptr += (step); \ 841 hptr += (step); \ 842 } \ 843 if ((lendind) > (hendind)) { \ 844 jpc_fix_pluseq(*lptr, jpc_fix_mul(*hptr, (twoalpha))); \ 845 } \ 846 } 847 848 #define NNS_SCALE(startptr, startind, endind, step, alpha) \ 849 { \ 850 register jpc_fix_t *ptr = (startptr); \ 851 register int n = (endind) - (startind); \ 852 while (n-- > 0) { \ 853 jpc_fix_muleq(*ptr, alpha); \ 854 ptr += (step); \ 855 } \ 856 } 857 858 static void jpc_ns_analyze(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x) 859 { 860 jpc_fix_t *startptr; 861 int startind; 862 int endind; 863 jpc_fix_t *lstartptr; 864 int lstartind; 865 int lendind; 866 jpc_fix_t *hstartptr; 867 int hstartind; 868 int hendind; 869 int interstep; 870 int intrastep; 871 int numseq; 872 873 /* Avoid compiler warnings about unused parameters. */ 874 qmfb = 0; 875 876 if (flags & JPC_QMFB1D_VERT) { 877 interstep = 1; 878 intrastep = jas_seq2d_rowstep(x); 879 numseq = jas_seq2d_width(x); 880 startind = jas_seq2d_ystart(x); 881 endind = jas_seq2d_yend(x); 882 } else { 883 interstep = jas_seq2d_rowstep(x); 884 intrastep = 1; 885 numseq = jas_seq2d_height(x); 886 startind = jas_seq2d_xstart(x); 887 endind = jas_seq2d_xend(x); 888 } 889 890 assert(startind < endind); 891 892 startptr = jas_seq2d_getref(x, jas_seq2d_xstart(x), jas_seq2d_ystart(x)); 893 if (!(flags & JPC_QMFB1D_RITIMODE)) { 894 while (numseq-- > 0) { 895 jpc_qmfb1d_setup(startptr, startind, endind, intrastep, 896 &lstartptr, &lstartind, &lendind, &hstartptr, 897 &hstartind, &hendind); 898 if (endind - startind > 1) { 899 jpc_qmfb1d_split(startptr, startind, endind, 900 intrastep, lstartptr, lstartind, lendind, 901 hstartptr, hstartind, hendind); 902 NNS_LIFT0(lstartptr, lstartind, lendind, 903 hstartptr, hstartind, hendind, intrastep, 904 jpc_dbltofix(-1.586134342)); 905 NNS_LIFT1(lstartptr, lstartind, lendind, 906 hstartptr, hstartind, hendind, intrastep, 907 jpc_dbltofix(-0.052980118)); 908 NNS_LIFT0(lstartptr, lstartind, lendind, 909 hstartptr, hstartind, hendind, intrastep, 910 jpc_dbltofix(0.882911075)); 911 NNS_LIFT1(lstartptr, lstartind, lendind, 912 hstartptr, hstartind, hendind, intrastep, 913 jpc_dbltofix(0.443506852)); 914 NNS_SCALE(lstartptr, lstartind, lendind, 915 intrastep, jpc_dbltofix(1.0/1.23017410558578)); 916 NNS_SCALE(hstartptr, hstartind, hendind, 917 intrastep, jpc_dbltofix(1.0/1.62578613134411)); 918 } else { 919 #if 0 920 if (lstartind == lendind) { 921 *startptr = jpc_fix_asl(*startptr, 1); 922 } 923 #endif 924 } 925 startptr += interstep; 926 } 927 } else { 928 /* The reversible integer-to-integer mode is not supported 929 for this transform. */ 930 abort(); 931 } 932 } 933 934 static void jpc_ns_synthesize(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x) 935 { 936 jpc_fix_t *startptr; 937 int startind; 938 int endind; 939 jpc_fix_t *lstartptr; 940 int lstartind; 941 int lendind; 942 jpc_fix_t *hstartptr; 943 int hstartind; 944 int hendind; 945 int interstep; 946 int intrastep; 947 int numseq; 948 949 /* Avoid compiler warnings about unused parameters. */ 950 qmfb = 0; 951 952 if (flags & JPC_QMFB1D_VERT) { 953 interstep = 1; 954 intrastep = jas_seq2d_rowstep(x); 955 numseq = jas_seq2d_width(x); 956 startind = jas_seq2d_ystart(x); 957 endind = jas_seq2d_yend(x); 958 } else { 959 interstep = jas_seq2d_rowstep(x); 960 intrastep = 1; 961 numseq = jas_seq2d_height(x); 962 startind = jas_seq2d_xstart(x); 963 endind = jas_seq2d_xend(x); 964 } 965 966 assert(startind < endind); 967 968 startptr = jas_seq2d_getref(x, jas_seq2d_xstart(x), jas_seq2d_ystart(x)); 969 if (!(flags & JPC_QMFB1D_RITIMODE)) { 970 while (numseq-- > 0) { 971 jpc_qmfb1d_setup(startptr, startind, endind, intrastep, 972 &lstartptr, &lstartind, &lendind, &hstartptr, 973 &hstartind, &hendind); 974 if (endind - startind > 1) { 975 NNS_SCALE(lstartptr, lstartind, lendind, 976 intrastep, jpc_dbltofix(1.23017410558578)); 977 NNS_SCALE(hstartptr, hstartind, hendind, 978 intrastep, jpc_dbltofix(1.62578613134411)); 979 NNS_LIFT1(lstartptr, lstartind, lendind, 980 hstartptr, hstartind, hendind, intrastep, 981 jpc_dbltofix(-0.443506852)); 982 NNS_LIFT0(lstartptr, lstartind, lendind, 983 hstartptr, hstartind, hendind, intrastep, 984 jpc_dbltofix(-0.882911075)); 985 NNS_LIFT1(lstartptr, lstartind, lendind, 986 hstartptr, hstartind, hendind, intrastep, 987 jpc_dbltofix(0.052980118)); 988 NNS_LIFT0(lstartptr, lstartind, lendind, 989 hstartptr, hstartind, hendind, intrastep, 990 jpc_dbltofix(1.586134342)); 991 jpc_qmfb1d_join(startptr, startind, endind, 992 intrastep, lstartptr, lstartind, lendind, 993 hstartptr, hstartind, hendind); 994 } else { 995 #if 0 996 if (lstartind == lendind) { 997 *startptr = jpc_fix_asr(*startptr, 1); 998 } 999 #endif 1000 } 1001 startptr += interstep; 1002 } 1003 } else { 1004 /* The reversible integer-to-integer mode is not supported 1005 for this transform. */ 1006 abort(); 1007 } 1008 } 1009 1010 /******************************************************************************\ 1011 * 1012 \******************************************************************************/ 1013 1014 jpc_qmfb1d_t *jpc_qmfb1d_make(int qmfbid) 1015 { 1016 jpc_qmfb1d_t *qmfb; 1017 if (!(qmfb = jpc_qmfb1d_create())) { 1018 return 0; 1019 } 1020 switch (qmfbid) { 1021 case JPC_QMFB1D_FT: 1022 qmfb->ops = &jpc_ft_ops; 1023 break; 1024 case JPC_QMFB1D_NS: 1025 qmfb->ops = &jpc_ns_ops; 1026 break; 1027 default: 1028 jpc_qmfb1d_destroy(qmfb); 1029 return 0; 1030 break; 1031 } 1032 return qmfb; 1033 } 1034 1035 static jpc_qmfb1d_t *jpc_qmfb1d_create() 1036 { 1037 jpc_qmfb1d_t *qmfb; 1038 if (!(qmfb = jas_malloc(sizeof(jpc_qmfb1d_t)))) { 1039 return 0; 1040 } 1041 qmfb->ops = 0; 1042 return qmfb; 1043 } 1044 1045 jpc_qmfb1d_t *jpc_qmfb1d_copy(jpc_qmfb1d_t *qmfb) 1046 { 1047 jpc_qmfb1d_t *newqmfb; 1048 1049 if (!(newqmfb = jpc_qmfb1d_create())) { 1050 return 0; 1051 } 1052 newqmfb->ops = qmfb->ops; 1053 return newqmfb; 1054 } 1055 1056 void jpc_qmfb1d_destroy(jpc_qmfb1d_t *qmfb) 1057 { 1058 jas_free(qmfb); 1059 } 1060 1061 /******************************************************************************\ 1062 * 1063 \******************************************************************************/ 1064 1065 void jpc_qmfb1d_getbands(jpc_qmfb1d_t *qmfb, int flags, uint_fast32_t xstart, 1066 uint_fast32_t ystart, uint_fast32_t xend, uint_fast32_t yend, int maxbands, 1067 int *numbandsptr, jpc_qmfb1dband_t *bands) 1068 { 1069 int start; 1070 int end; 1071 1072 assert(maxbands >= 2); 1073 1074 if (flags & JPC_QMFB1D_VERT) { 1075 start = ystart; 1076 end = yend; 1077 } else { 1078 start = xstart; 1079 end = xend; 1080 } 1081 assert(jpc_qmfb1d_getnumchans(qmfb) == 2); 1082 assert(start <= end); 1083 bands[0].start = JPC_CEILDIVPOW2(start, 1); 1084 bands[0].end = JPC_CEILDIVPOW2(end, 1); 1085 bands[0].locstart = start; 1086 bands[0].locend = start + bands[0].end - bands[0].start; 1087 bands[1].start = JPC_FLOORDIVPOW2(start, 1); 1088 bands[1].end = JPC_FLOORDIVPOW2(end, 1); 1089 bands[1].locstart = bands[0].locend; 1090 bands[1].locend = bands[1].locstart + bands[1].end - bands[1].start; 1091 assert(bands[1].locend == end); 1092 *numbandsptr = 2; 1093 } 1094 1095 /******************************************************************************\ 1096 * 1097 \******************************************************************************/ 1098 1099 int jpc_qmfb1d_getnumchans(jpc_qmfb1d_t *qmfb) 1100 { 1101 return (*qmfb->ops->getnumchans)(qmfb); 1102 } 1103 1104 int jpc_qmfb1d_getanalfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters) 1105 { 1106 return (*qmfb->ops->getanalfilters)(qmfb, len, filters); 1107 } 1108 1109 int jpc_qmfb1d_getsynfilters(jpc_qmfb1d_t *qmfb, int len, jas_seq2d_t **filters) 1110 { 1111 return (*qmfb->ops->getsynfilters)(qmfb, len, filters); 1112 } 1113 1114 void jpc_qmfb1d_analyze(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x) 1115 { 1116 (*qmfb->ops->analyze)(qmfb, flags, x); 1117 } 1118 1119 void jpc_qmfb1d_synthesize(jpc_qmfb1d_t *qmfb, int flags, jas_seq2d_t *x) 1120 { 1121 (*qmfb->ops->synthesize)(qmfb, flags, x); 1122 }