jpc_qmfb.c - vx32 - Local 9vx git repository for patches.
(HTM) git clone git://r-36.net/vx32
(DIR) Log
(DIR) Files
(DIR) Refs
---
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 }