tcd.c

关于JPEG2000压缩标准的模板程序的内容,是我们导师给的,据说是国外的版本.

                                continue;
                            }
                            if (dd/dr>thresh) {
                                n=passno+1;
                            }
                        }
                        layer->numpasses=n-cblk->numpassesinlayers;
                        if (!layer->numpasses) {
                            continue;
                        }
                        if (cblk->numpassesinlayers==0) {
                            layer->len=cblk->passes[n-1].rate;
                            layer->data=cblk->data;
                        } else {
                            layer->len=cblk->passes[n-1].rate-cblk->passes[cblk->numpassesinlayers-1].rate;
                            layer->data=cblk->data+cblk->passes[cblk->numpassesinlayers-1].rate;
                        }
                        if (final) {
                            cblk->numpassesinlayers=n;
                        }
                    }
                }
            }
        }
    }
}

void tcd_rateallocate(unsigned char *dest, int len) {
    int compno, resno, bandno, precno, cblkno, passno, layno;
    double min, max;
    min=DBL_MAX;
    max=0;
    for (compno=0; compno<tcd_tile->numcomps; compno++) {
        tcd_tilecomp_t *tilec=&tcd_tile->comps[compno];
        for (resno=0; resno<tilec->numresolutions; resno++) {
            tcd_resolution_t *res=&tilec->resolutions[resno];
            for (bandno=0; bandno<res->numbands; bandno++) {
                tcd_band_t *band=&res->bands[bandno];
                for (precno=0; precno<res->pw*res->ph; precno++) {
                    tcd_precinct_t *prc=&band->precincts[precno];
                    for (cblkno=0; cblkno<prc->cw*prc->ch; cblkno++) {
                        tcd_cblk_t *cblk=&prc->cblks[cblkno];
                        for (passno=0; passno<cblk->totalpasses; passno++) {
                            tcd_pass_t *pass=&cblk->passes[passno];
                            int dr;
                            double dd, rdslope;
                            if (passno==0) {
                                dr=pass->rate;
                                dd=pass->distortiondec;
                            } else {
                                dr=pass->rate-cblk->passes[passno-1].rate;
                                dd=pass->distortiondec-cblk->passes[passno-1].distortiondec;
                            }
                            if (dr==0) {
                                continue;
                            }
                            rdslope=dd/dr;
                            if (rdslope<min) {
                                min=rdslope;
                            }
                            if (rdslope>max) {
                                max=rdslope;
                            }
                        }
                    }
                }
            }
        }
    }
    for (layno=0; layno<tcd_tcp->numlayers; layno++) {
        volatile double lo=min;
        volatile double hi=max;
        volatile int success=0;
        volatile int maxlen=int_min(tcd_tcp->rates[layno], len);
        volatile double goodthresh;
        volatile int goodlen;
        volatile jmp_buf oldenv;
        volatile int i;
        memcpy((void*)oldenv, j2k_error, sizeof(jmp_buf));
        for (i=0; i<32; i++) {
            volatile double thresh=(lo+hi)/2;
            int l;
            tcd_makelayer(layno, thresh, 0);
            if (setjmp(j2k_error)) {
                lo=thresh;
                continue;
            }
            l=t2_encode_packets(tcd_img, tcd_cp, tcd_tileno, tcd_tile, layno+1, dest, maxlen);
            fprintf(stderr, "rate alloc: len=%d, max=%d\n", l, maxlen);
            hi=thresh;
            success=1;
            goodthresh=thresh;
            goodlen=l;
        }
        memcpy(j2k_error, (void*)oldenv, sizeof(jmp_buf));
        if (!success) {
            longjmp(j2k_error, 1);
        }
        tcd_makelayer(layno, goodthresh, 1);
    }
}

int tcd_encode_tile(int tileno, unsigned char *dest, int len) {
    int compno;
    int l;
    clock_t time1, time2, time3, time4, time5, time6, time7;
    tcd_tile_t *tile;
    tcd_tileno=tileno;
    tcd_tile=&tcd_image.tiles[tileno];
    tcd_tcp=&tcd_cp->tcps[tileno];
    tile=tcd_tile;

    time7=clock();

    time1=clock();
    for (compno=0; compno<tile->numcomps; compno++) {
        int i, j;
        int tw, w;
        tcd_tilecomp_t *tilec=&tile->comps[compno];
        int adjust=tcd_img->comps[compno].sgnd?0:1<<(tcd_img->comps[compno].prec-1);
        tw=tilec->x1-tilec->x0;
        w=int_ceildiv(tcd_img->x1-tcd_img->x0, tcd_img->comps[compno].dx);
        for (j=tilec->y0; j<tilec->y1; j++) {
            for (i=tilec->x0; i<tilec->x1; i++) {
                if (tcd_tcp->tccps[compno].qmfbid==1) {
                    tilec->data[i-tilec->x0+(j-tilec->y0)*tw]=tcd_img->comps[compno].data[i+j*w]-adjust;
                } else if (tcd_tcp->tccps[compno].qmfbid==0) {
                    tilec->data[i-tilec->x0+(j-tilec->y0)*tw]=(tcd_img->comps[compno].data[i+j*w]-adjust)<<13;
                }
            }
        }
    }
    time1=clock()-time1;

    time2=clock();
    if (tcd_tcp->mct) {
        if (tcd_tcp->tccps[0].qmfbid==0) {
            mct_encode_real(tile->comps[0].data, tile->comps[1].data, tile->comps[2].data, (tile->comps[0].x1-tile->comps[0].x0)*(tile->comps[0].y1-tile->comps[0].y0));
        } else {
            mct_encode(tile->comps[0].data, tile->comps[1].data, tile->comps[2].data, (tile->comps[0].x1-tile->comps[0].x0)*(tile->comps[0].y1-tile->comps[0].y0));
        }
    }
    time2=clock()-time2;

    time3=clock();
    for (compno=0; compno<tile->numcomps; compno++) {
        tcd_tilecomp_t *tilec=&tile->comps[compno];
        if (tcd_tcp->tccps[compno].qmfbid==1) {
            dwt_encode(tilec->data, tilec->x1-tilec->x0, tilec->y1-tilec->y0, tilec->numresolutions-1);
        } else if (tcd_tcp->tccps[compno].qmfbid==0) {
            dwt_encode_real(tilec->data, tilec->x1-tilec->x0, tilec->y1-tilec->y0, tilec->numresolutions-1);
        }
    }
    time3=clock()-time3;

    time4=clock();
    t1_init_luts();

    t1_encode_cblks(tile, tcd_tcp);
    time4=clock()-time4;

    time5=clock();
    tcd_rateallocate(dest, len);
    time5=clock()-time5;

    time6=clock();
    l=t2_encode_packets(tcd_img, tcd_cp, tileno, tile, tcd_tcp->numlayers, dest, len);
    time6=clock()-time6;

    time7=clock()-time7;

    printf("tile encoding times:\n");
    printf("img->tile: %d.%.3d s\n", time1/CLOCKS_PER_SEC, (time1%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("mct:       %d.%.3d s\n", time2/CLOCKS_PER_SEC, (time2%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("dwt:       %d.%.3d s\n", time3/CLOCKS_PER_SEC, (time3%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("tier 1:    %d.%.3d s\n", time4/CLOCKS_PER_SEC, (time4%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("ratealloc: %d.%.3d s\n", time5/CLOCKS_PER_SEC, (time5%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("tier 2:    %d.%.3d s\n", time6/CLOCKS_PER_SEC, (time6%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("total:     %d.%.3d s\n", time7/CLOCKS_PER_SEC, (time7%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);

    return l;
}

int tcd_decode_tile(unsigned char *src, int len, int tileno) {
    int l;
    int compno;
    int eof=0;
    jmp_buf oldenv;
    clock_t time1, time2, time3, time4, time5, time6;

    tcd_tile_t *tile;
    tcd_tileno=tileno;
    tcd_tile=&tcd_image.tiles[tileno];
    tcd_tcp=&tcd_cp->tcps[tileno];
    tile=tcd_tile;

    time6=clock();

    time1=clock();
    memcpy(oldenv, j2k_error, sizeof(jmp_buf));
    if (setjmp(j2k_error)) {
        eof=1;
        fprintf(stderr, "tcd_decode: incomplete bistream\n");
    } else {
        l=t2_decode_packets(src, len, tcd_img, tcd_cp, tileno, tile);
    }
    memcpy(j2k_error, oldenv, sizeof(jmp_buf));
    time1=clock()-time1;

    time2=clock();
    t1_init_luts();
    t1_decode_cblks(tile, tcd_tcp);
    time2=clock()-time2;

    time3=clock();
    for (compno=0; compno<tile->numcomps; compno++) {
        tcd_tilecomp_t *tilec=&tile->comps[compno];
        if (tcd_tcp->tccps[compno].qmfbid==1) {
            dwt_decode(tilec->data, tilec->x1-tilec->x0, tilec->y1-tilec->y0, tilec->numresolutions-1);
        } else if (tcd_tcp->tccps[compno].qmfbid==0) {
            dwt_decode_real(tilec->data, tilec->x1-tilec->x0, tilec->y1-tilec->y0, tilec->numresolutions-1);
        }
    }
    time3=clock()-time3;

    time4=clock();
    if (tcd_tcp->mct) {
        if (tcd_tcp->tccps[0].qmfbid==0) {
            mct_decode_real(tile->comps[0].data, tile->comps[1].data, tile->comps[2].data, (tile->comps[0].x1-tile->comps[0].x0)*(tile->comps[0].y1-tile->comps[0].y0));
        } else {
            mct_decode(tile->comps[0].data, tile->comps[1].data, tile->comps[2].data, (tile->comps[0].x1-tile->comps[0].x0)*(tile->comps[0].y1-tile->comps[0].y0));
        }
    }
    time4=clock()-time4;

    time5=clock();
    for (compno=0; compno<tile->numcomps; compno++) {
        tcd_tilecomp_t *tilec=&tile->comps[compno];
        int adjust=tcd_img->comps[compno].sgnd?0:1<<(tcd_img->comps[compno].prec-1);
        int min=tcd_img->comps[compno].sgnd?-(1<<(tcd_img->comps[compno].prec-1)):0;
        int max=tcd_img->comps[compno].sgnd?(1<<(tcd_img->comps[compno].prec-1))-1:(1<<tcd_img->comps[compno].prec)-1;
        int tw=tilec->x1-tilec->x0;
        int w=int_ceildiv(tcd_img->x1-tcd_img->x0, tcd_img->comps[compno].dx);
        int i, j;
        for (j=tilec->y0; j<tilec->y1; j++) {
            for (i=tilec->x0; i<tilec->x1; i++) {
                int v;
                if (tcd_tcp->tccps[compno].qmfbid==1) {
                    v=tilec->data[i-tilec->x0+(j-tilec->y0)*tw];
                } else if (tcd_tcp->tccps[compno].qmfbid==0) {
                    v=tilec->data[i-tilec->x0+(j-tilec->y0)*tw]>>13;
                }
                v+=adjust;
                tcd_img->comps[compno].data[i+j*w]=int_clamp(v, min, max);
            }
        }
    }
    time5=clock()-time5;

    time6=clock()-time6;

    printf("tile decoding times:\n");
    printf("tier 2:    %d.%.3d s\n", time1/CLOCKS_PER_SEC, (time1%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("tier 1:    %d.%.3d s\n", time2/CLOCKS_PER_SEC, (time2%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("dwt:       %d.%.3d s\n", time3/CLOCKS_PER_SEC, (time3%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("mct:       %d.%.3d s\n", time4/CLOCKS_PER_SEC, (time4%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("tile->img: %d.%.3d s\n", time5/CLOCKS_PER_SEC, (time5%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);
    printf("total:     %d.%.3d s\n", time6/CLOCKS_PER_SEC, (time6%CLOCKS_PER_SEC)*1000/CLOCKS_PER_SEC);

    if (eof) {
        longjmp(j2k_error, 1);
    }
    return l;
}