diff mbox

[FFmpeg-devel] Add FITS Decoder

Message ID 1497361980-30310-1-git-send-email-paraschadha18@gmail.com
State Superseded
Headers show

Commit Message

Paras June 13, 2017, 1:53 p.m. UTC
Above mentioned changes have been done.

Signed-off-by: Paras Chadha <paraschadha18@gmail.com>
---
 Changelog               |   1 +
 doc/general.texi        |   2 +
 libavcodec/Makefile     |   1 +
 libavcodec/allcodecs.c  |   1 +
 libavcodec/avcodec.h    |   1 +
 libavcodec/codec_desc.c |   8 +
 libavcodec/fitsdec.c    | 515 ++++++++++++++++++++++++++++++++++++++++++++++++
 libavcodec/version.h    |   2 +-
 libavformat/img2.c      |   1 +
 9 files changed, 531 insertions(+), 1 deletion(-)
 create mode 100644 libavcodec/fitsdec.c

--
2.4.11

Comments

Paul B Mahol June 13, 2017, 2:15 p.m. UTC | #1
On 6/13/17, Paras Chadha <paraschadha18@gmail.com> wrote:
> Above mentioned changes have been done.
>
> Signed-off-by: Paras Chadha <paraschadha18@gmail.com>
> ---
>  Changelog               |   1 +
>  doc/general.texi        |   2 +
>  libavcodec/Makefile     |   1 +
>  libavcodec/allcodecs.c  |   1 +
>  libavcodec/avcodec.h    |   1 +
>  libavcodec/codec_desc.c |   8 +
>  libavcodec/fitsdec.c    | 515
> ++++++++++++++++++++++++++++++++++++++++++++++++
>  libavcodec/version.h    |   2 +-
>  libavformat/img2.c      |   1 +
>  9 files changed, 531 insertions(+), 1 deletion(-)
>  create mode 100644 libavcodec/fitsdec.c
>
> diff --git a/Changelog b/Changelog
> index cd91f63..db495d9 100644
> --- a/Changelog
> +++ b/Changelog
> @@ -20,6 +20,7 @@ version <next>:
>  - sofalizer filter switched to libmysofa
>  - Gremlin Digital Video demuxer and decoder
>  - headphone audio filter
> +- FITS decoder
>
>  version 3.3:
>  - CrystalHD decoder moved to new decode API
> diff --git a/doc/general.texi b/doc/general.texi
> index 8f582d5..213e50e 100644
> --- a/doc/general.texi
> +++ b/doc/general.texi
> @@ -591,6 +591,8 @@ following image formats are supported:
>      @tab Digital Picture Exchange
>  @item EXR          @tab   @tab X
>      @tab OpenEXR
> +@item FITS         @tab   @tab X
> +    @tab FITS image format
>  @item JPEG         @tab X @tab X
>      @tab Progressive JPEG is not supported.
>  @item JPEG 2000    @tab X @tab X
> diff --git a/libavcodec/Makefile b/libavcodec/Makefile
> index 2e7f19d..2cfd1fe 100644
> --- a/libavcodec/Makefile
> +++ b/libavcodec/Makefile
> @@ -291,6 +291,7 @@ OBJS-$(CONFIG_FFV1_DECODER)            += ffv1dec.o
> ffv1.o
>  OBJS-$(CONFIG_FFV1_ENCODER)            += ffv1enc.o ffv1.o
>  OBJS-$(CONFIG_FFWAVESYNTH_DECODER)     += ffwavesynth.o
>  OBJS-$(CONFIG_FIC_DECODER)             += fic.o
> +OBJS-$(CONFIG_FITS_DECODER)            += fitsdec.o
>  OBJS-$(CONFIG_FLAC_DECODER)            += flacdec.o flacdata.o flac.o
>  OBJS-$(CONFIG_FLAC_ENCODER)            += flacenc.o flacdata.o flac.o
> vorbis_data.o
>  OBJS-$(CONFIG_FLASHSV_DECODER)         += flashsv.o
> diff --git a/libavcodec/allcodecs.c b/libavcodec/allcodecs.c
> index 27110e1..8a89264 100644
> --- a/libavcodec/allcodecs.c
> +++ b/libavcodec/allcodecs.c
> @@ -186,6 +186,7 @@ static void register_all(void)
>      REGISTER_ENCDEC (FFV1,              ffv1);
>      REGISTER_ENCDEC (FFVHUFF,           ffvhuff);
>      REGISTER_DECODER(FIC,               fic);
> +    REGISTER_DECODER(FITS,              fits);
>      REGISTER_ENCDEC (FLASHSV,           flashsv);
>      REGISTER_ENCDEC (FLASHSV2,          flashsv2);
>      REGISTER_DECODER(FLIC,              flic);
> diff --git a/libavcodec/avcodec.h b/libavcodec/avcodec.h
> index dcdcfe0..e0c1089 100644
> --- a/libavcodec/avcodec.h
> +++ b/libavcodec/avcodec.h
> @@ -447,6 +447,7 @@ enum AVCodecID {
>      AV_CODEC_ID_SRGC,
>      AV_CODEC_ID_SVG,
>      AV_CODEC_ID_GDV,
> +    AV_CODEC_ID_FITS,
>
>      /* various PCM "codecs" */
>      AV_CODEC_ID_FIRST_AUDIO = 0x10000,     ///< A dummy id pointing at the
> start of audio codecs
> diff --git a/libavcodec/codec_desc.c b/libavcodec/codec_desc.c
> index cf1246e..0112517 100644
> --- a/libavcodec/codec_desc.c
> +++ b/libavcodec/codec_desc.c
> @@ -1464,6 +1464,14 @@ static const AVCodecDescriptor codec_descriptors[] =
> {
>                       AV_CODEC_PROP_LOSSLESS,
>      },
>      {
> +        .id        = AV_CODEC_ID_FITS,
> +        .type      = AVMEDIA_TYPE_VIDEO,
> +        .name      = "fits",
> +        .long_name = NULL_IF_CONFIG_SMALL("Flexible Image Transport
> System"),
> +        .props     = AV_CODEC_PROP_INTRA_ONLY | AV_CODEC_PROP_LOSSY |
> +                     AV_CODEC_PROP_LOSSLESS,
> +    },
> +    {
>          .id        = AV_CODEC_ID_GIF,
>          .type      = AVMEDIA_TYPE_VIDEO,
>          .name      = "gif",
> diff --git a/libavcodec/fitsdec.c b/libavcodec/fitsdec.c
> new file mode 100644
> index 0000000..dfc4995
> --- /dev/null
> +++ b/libavcodec/fitsdec.c
> @@ -0,0 +1,515 @@
> +/*
> + * FITS image decoder
> + * Copyright (c) 2017 Paras Chadha
> + *
> + * This file is part of FFmpeg.
> + *
> + * FFmpeg is free software; you can redistribute it and/or
> + * modify it under the terms of the GNU Lesser General Public
> + * License as published by the Free Software Foundation; either
> + * version 2.1 of the License, or (at your option) any later version.
> + *
> + * FFmpeg is distributed in the hope that it will be useful,
> + * but WITHOUT ANY WARRANTY; without even the implied warranty of
> + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
> + * Lesser General Public License for more details.
> + *
> + * You should have received a copy of the GNU Lesser General Public
> + * License along with FFmpeg; if not, write to the Free Software
> + * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301
> USA
> + */
> +
> +/**
> + * @file
> + * FITS image decoder
> + * It supports all 2-d images alongwith, bzero, bscale and blank keywords.
> + * RGBA images are supported as NAXIS3 = 3 or 4 i.e. Planes in RGBA order.
> Also CTYPE = 'RGB ' should be present.
> + * It currently does not support XTENSION keyword.
> + * Also to interpret data, values are linearly scaled using min-max scaling
> but not RGB images.
> + */
> +
> +#include "avcodec.h"
> +#include "internal.h"
> +#include <float.h>
> +#include "libavutil/intreadwrite.h"
> +
> +/**
> + * Structure to store the header keywords in FITS file
> + */
> +typedef struct FITSContext {
> +    char simple;
> +    int bitpix;
> +    int blank;
> +    int naxis;
> +    int naxisn[3];
> +    int rgb; /**< 1 if file contains RGB image, 0 otherwise */
> +    double bscale;
> +    double bzero;
> +    double data_min;
> +    double data_max;
> +} FITSDecContext;
> +
> +/**
> + * function calculates the data_min and data_max values from the data.
> + * This is called if the values are not present in the header.
> + * @param ptr8 - pointer to the data
> + * @param header - pointer to the header
> + * @return 1, if calculated successfully, otherwise AVERROR_INVALIDDATA
> + */
> +static int fill_data_min_max(const uint8_t * ptr8, FITSDecContext * header,
> const uint8_t * end)
> +{
> +    int16_t t16;
> +    int32_t t32;
> +    int64_t t64;
> +    float tflt;
> +    double tdbl;
> +    int i, j;
> +
> +    header->data_min = DBL_MAX;
> +    header->data_max = DBL_MIN;
> +    switch (header->bitpix) {
> +        case -64:
> +            t64 = AV_RB64(ptr8);
> +            memcpy(&tdbl, &t64, 8);
> +            for (i = 0; i < header->naxisn[1]; i++) {
> +                for (j = 0; j < header->naxisn[0]; j++) {
> +                    t64 = AV_RB64(ptr8);
> +                    memcpy(&tdbl, &t64, 8);

Why are you using memcpy calls for this? Why not assign directly?

> +                    if (tdbl > header->data_max)
> +                        header->data_max = tdbl;
> +                    if (tdbl < header->data_min)
> +                        header->data_min = tdbl;
> +                    ptr8 += 8;
> +                }
> +            }
> +            break;
> +        case -32:
> +            t32 = AV_RB32(ptr8);
> +            memcpy(&tflt, &t32, 4);
> +            for (i = 0; i < header->naxisn[1]; i++) {
> +                for (j = 0; j < header->naxisn[0]; j++) {
> +                    t32 = AV_RB32(ptr8);
> +                    memcpy(&tflt, &t32, 4);
> +                    if (tflt > header->data_max)
> +                        header->data_max = tflt;
> +                    if (tflt < header->data_min)
> +                        header->data_min = tflt;
> +                    ptr8 += 4;
> +                }
> +            }
> +            break;
> +        case 8:
> +            for (i = 0; i < header->naxisn[1]; i++) {
> +                for (j = 0; j < header->naxisn[0]; j++) {
> +                    if (ptr8[0] != header->blank) {
> +                        if (ptr8[0] > header->data_max)
> +                            header->data_max = ptr8[0];
> +                        if (ptr8[0] < header->data_min)
> +                            header->data_min = ptr8[0];
> +                    }
> +                    ptr8++;
> +                }
> +            }
> +            break;
> +        case 16:
> +            t16 = AV_RB16(ptr8);
> +            for (i = 0; i < header->naxisn[1]; i++) {
> +                for (j = 0; j < header->naxisn[0]; j++) {
> +                    t16 = AV_RB16(ptr8);
> +                    if (t16 != header->blank) {
> +                        if (t16 > header->data_max)
> +                            header->data_max = t16;
> +                        if (t16 < header->data_min)
> +                            header->data_min = t16;
> +                    }
> +                    ptr8 += 2;
> +                }
> +            }
> +            break;
> +        case 32:
> +            t32 = AV_RB32(ptr8);
> +            for (i = 0; i < header->naxisn[1]; i++) {
> +                for (j = 0; j < header->naxisn[0]; j++) {
> +                    t32 = AV_RB32(ptr8);
> +                    if (t32 != header->blank) {
> +                        if (t32 > header->data_max)
> +                            header->data_max = t32;
> +                        if (t32 < header->data_min)
> +                            header->data_min = t32;
> +                    }
> +                    ptr8 += 4;
> +                }
> +            }
> +            break;
> +        case 64:
> +            t64 = AV_RB64(ptr8);
> +            for (i = 0; i < header->naxisn[1]; i++) {
> +                for (j = 0; j < header->naxisn[0]; j++) {
> +                    t64 = AV_RB64(ptr8);
> +                    if (t64 != header->blank) {
> +                        if (t64 > header->data_max)
> +                            header->data_max = t64;
> +                        if (t64 < header->data_min)
> +                            header->data_min = t64;
> +                    }
> +                    ptr8 += 8;
> +                }
> +            }
> +            break;
> +        default:
> +            return AVERROR_INVALIDDATA;
> +    }
> +    return 1;
> +}
> +
> +/**
> + * function reads the fits header and stores the values in FITSDecContext
> pointed by header
> + * @param avctx - AVCodec context
> + * @param ptr - pointer to pointer to the data
> + * @param header - pointer to the FITSDecContext
> + * @param end - pointer to end of packet
> + * @return 1, if calculated successfully, otherwise AVERROR_INVALIDDATA
> + */
> +static int fits_read_header(AVCodecContext *avctx, const uint8_t **ptr,
> FITSDecContext * header, const uint8_t * end)
> +{
> +    const uint8_t *ptr8 = *ptr;
> +    int lines_read = 0, i, dim_no, t, data_min_found = 0, data_max_found =
> 0, ret;
> +    uint64_t size=1;
> +    double d;
> +
> +    header->blank = 0;
> +    header->bscale = 1.0;
> +    header->bzero = 0;
> +    header->rgb = 0;
> +
> +    if (end - ptr8 < 80)
> +        return AVERROR_INVALIDDATA;
> +
> +    if (sscanf(ptr8, "SIMPLE = %c", &header->simple) != 1) {
> +        av_log(avctx, AV_LOG_ERROR, "missing SIMPLE keyword\n");
> +        return AVERROR_INVALIDDATA;
> +    }
> +
> +    if (header->simple == 'F') {
> +        av_log(avctx, AV_LOG_WARNING, "not a standard FITS file\n");
> +    } else if (header->simple != 'T') {
> +        av_log(avctx, AV_LOG_ERROR, "invalid SIMPLE value, SIMPLE = %c\n",
> header->simple);
> +        return AVERROR_INVALIDDATA;
> +    }
> +
> +    ptr8 += 80;
> +    lines_read++;
> +
> +    if (end - ptr8 < 80)
> +        return AVERROR_INVALIDDATA;
> +
> +    if (sscanf(ptr8, "BITPIX = %d", &header->bitpix) != 1) {
> +        av_log(avctx, AV_LOG_ERROR, "missing BITPIX keyword\n");
> +        return AVERROR_INVALIDDATA;
> +    }
> +
> +    size = abs(header->bitpix) >> 3;
> +    ptr8 += 80;
> +    lines_read++;
> +
> +    if (end - ptr8 < 80)
> +        return AVERROR_INVALIDDATA;
> +
> +    if (sscanf(ptr8, "NAXIS = %d", &header->naxis) != 1) {
> +        av_log(avctx, AV_LOG_ERROR, "missing NAXIS keyword\n");
> +        return AVERROR_INVALIDDATA;
> +    }
> +
> +    if (!header->naxis) {
> +        av_log(avctx, AV_LOG_ERROR, "No image data found, NAXIS = %d\n",
> header->naxis);
> +        return AVERROR_INVALIDDATA;
> +    }
> +
> +    if (header->naxis != 2 && header->naxis != 3) {
> +        av_log(avctx, AV_LOG_ERROR, "unsupported number of dimensions,
> NAXIS = %d\n", header->naxis);
> +        return AVERROR_INVALIDDATA;
> +    }
> +
> +    ptr8 += 80;
> +    lines_read++;
> +
> +    for (i = 0; i < header->naxis; i++) {
> +        if (end - ptr8 < 80)
> +            return AVERROR_INVALIDDATA;
> +
> +        if (sscanf(ptr8, "NAXIS%d = %d", &dim_no, &header->naxisn[i]) != 2
> || dim_no != i+1) {
> +            av_log(avctx, AV_LOG_ERROR, "missing NAXIS%d keyword\n", i+1);
> +            return AVERROR_INVALIDDATA;
> +        }
> +
> +        size *= header->naxisn[i];
> +        ptr8 += 80;
> +        lines_read++;
> +    }
> +
> +    if (end - ptr8 < 80)
> +            return AVERROR_INVALIDDATA;
> +
> +    while (strncmp(ptr8, "END", 3)) {
> +        if (sscanf(ptr8, "BLANK = %d", &t) == 1) {
> +            header->blank = t;
> +        } else if (sscanf(ptr8, "BSCALE = %lf", &d) == 1) {
> +            header->bscale = d;
> +        } else if (sscanf(ptr8, "BZERO = %lf", &d) == 1) {
> +            header->bzero = d;
> +        } else if (sscanf(ptr8, "DATAMAX = %lf", &d) == 1) {
> +            data_max_found = 1;
> +            header->data_max = d;
> +        } else if (sscanf(ptr8, "DATAMIN = %lf", &d) == 1) {
> +            data_min_found = 1;
> +            header->data_min = d;
> +        } else if (!strncmp(ptr8, "CTYPE3  = 'RGB", 14)) {
> +            header->rgb = 1;
> +            if (header->naxis != 3 || (header->naxisn[2] != 3 &&
> header->naxisn[2] != 4)) {
> +                av_log(avctx, AV_LOG_ERROR, "File contains RGB image but
> NAXIS = %d and NAXIS3 = %d\n", header->naxis, header->naxisn[2]);
> +                return AVERROR_INVALIDDATA;
> +            }
> +        }
> +        ptr8 += 80;
> +        lines_read++;
> +
> +        if (end - ptr8 < 80)
> +            return AVERROR_INVALIDDATA;
> +    }
> +
> +    if (!header->rgb && header->naxis != 2){
> +        av_log(avctx, AV_LOG_ERROR, "unsupported number of dimensions,
> NAXIS = %d\n", header->naxis);
> +        return AVERROR_INVALIDDATA;
> +    }
> +
> +    ptr8 += 80;
> +    lines_read++;
> +    lines_read %= 36;
> +
> +    t = ((36 - lines_read) % 36) * 80;
> +    if (end - ptr8 < t)
> +        return AVERROR_INVALIDDATA;
> +    ptr8 += t;
> +    *ptr = ptr8;
> +
> +    if (end - ptr8 < size)
> +        return AVERROR_INVALIDDATA;
> +
> +    if (!header->rgb && (!data_min_found || !data_max_found)) {
> +        if ((ret = fill_data_min_max(ptr8, header, end)) < 0) {
> +            av_log(avctx, AV_LOG_ERROR, "invalid BITPIX, %d\n",
> header->bitpix);
> +            return AVERROR_INVALIDDATA;
> +        }
> +    } else {
> +        /*
> +         * instead of applying bscale and bzero to every element, we can do
> inverse transformation on data_min and
> +         * data_max
> +         */
> +        header->data_min = (header->data_min - header->bzero) /
> header->bscale;
> +        header->data_max = (header->data_max - header->bzero) /
> header->bscale;
> +    }
> +    return 1;
> +}
> +
> +static int fits_decode_frame(AVCodecContext *avctx, void *data, int
> *got_frame, AVPacket *avpkt)
> +{
> +    AVFrame *p=data;
> +    const uint8_t *ptr8 = avpkt->data, *end;
> +    int16_t t16;
> +    int32_t t32;
> +    int64_t t64;
> +    float   tflt;
> +    double  tdbl;
> +    int ret, i, j;
> +    uint8_t *dst8;
> +    uint16_t *dst16;
> +    uint32_t *dst32;
> +    uint64_t *dst64, size, r, g, b, a, t;
> +    FITSDecContext * header = avctx->priv_data;
> +
> +    end = ptr8 + avpkt->size;
> +    if (ret = fits_read_header(avctx, &ptr8, header, end) < 0)
> +        return ret;
> +
> +    size = (header->naxisn[0]) * (header->naxisn[1]);
> +
> +    if (header->rgb) {
> +        if (header->bitpix == 8) {
> +            avctx->pix_fmt = AV_PIX_FMT_RGB32;
> +        } else if (header->bitpix == 16) {
> +            avctx->pix_fmt = AV_PIX_FMT_RGBA64;
> +        } else {
> +            av_log(avctx, AV_LOG_ERROR, "unsupported BITPIX = %d\n",
> header->bitpix);
> +            return AVERROR_INVALIDDATA;
> +        }
> +    } else {
> +        if (header->bitpix == 8) {
> +            avctx->pix_fmt = AV_PIX_FMT_GRAY8;
> +        } else {
> +            avctx->pix_fmt = AV_PIX_FMT_GRAY16;
> +        }
> +    }
> +
> +    if ((ret = ff_set_dimensions(avctx, header->naxisn[0],
> header->naxisn[1])) < 0)
> +        return ret;
> +
> +    if ((ret = ff_get_buffer(avctx, p, 0)) < 0)
> +        return ret;
> +
> +    if (header->rgb) {
> +        if (header->bitpix == 8) {
> +            for (i = 0; i < avctx->height; i++) {
> +                /*
> +                 * FITS stores images with bottom row first. Therefore we
> have
> +                 * to fill the image from bottom to top.
> +                 */
> +                dst32 = (uint32_t *)(p->data[0] + (avctx->height-i-1)*
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +                    if (header->naxisn[2] == 4) {
> +                        if (ptr8[size * 3] != header->blank)
> +                            t = ptr8[size * 3] * header->bscale +
> header->bzero;
> +                        a = t << 24;
> +                    } else {
> +                        a = (255 << 24);
> +                    }
> +
> +                    if (ptr8[0] != header->blank)
> +                        t = ptr8[0] * header->bscale + header->bzero;
> +                    r = t << 16;
> +
> +                    if (ptr8[size] != header->blank)
> +                        t = ptr8[size] * header->bscale + header->bzero;
> +                    g = t << 8;
> +
> +                    if (ptr8[size * 2] != header->blank)
> +                        t = ptr8[size * 2] * header->bscale +
> header->bzero;
> +                    b = t;
> +
> +                    *dst32++ = ((uint32_t)a) | ((uint32_t)r) |
> ((uint32_t)g) | ((uint32_t)b);
> +                    ptr8++;
> +                }
> +            }
> +        } else if (header->bitpix == 16) {
> +            // not tested ....
> +            for (i = 0; i < avctx->height; i++) {
> +                dst64 = (uint64_t *)(p->data[0] + (avctx->height-i-1) *
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +
> +                    if (header->naxisn[2] == 4) {
> +                        t = ((ptr8[size * 3] << 8) | ptr8[size * 3 + 1]);
> +                        if (t != header->blank)
> +                            t = t*header->bscale + header->bzero;
> +                        a = t << 48;
> +                    } else {
> +                        a = 65535ULL << 48;
> +                    }
> +
> +                    t = ptr8[0] << 8 | ptr8[1];
> +                    if (t != header->blank)
> +                        t = t*header->bscale + header->bzero;
> +                    r = t << 32;
> +
> +                    t = ptr8[size] << 8 | ptr8[size + 1];
> +                    if (t != header->blank)
> +                        t = t*header->bscale + header->bzero;
> +                    g = t << 16;
> +
> +                    t = ptr8[size * 2] << 8 | ptr8[size * 2 + 1];
> +                    if (t != header->blank)
> +                        t = t*header->bscale + header->bzero;
> +                    b = t;
> +
> +                    *dst64++ = a | r | g | b;
> +                    ptr8 += 2;
> +                }
> +            }
> +        }
> +    } else {
> +        if (header->bitpix == 8) {
> +            for (i = 0; i < avctx->height; i++) {
> +                dst8 = (uint8_t *) (p->data[0] + (avctx->height-i-1)*
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +                    if (ptr8[0] != header->blank) {
> +                        *dst8++ = ((ptr8[0] - header->data_min) * 255) /
> (header->data_max - header->data_min);
> +                    } else {
> +                        *dst8++ = ptr8[0];
> +                    }
> +                    ptr8++;
> +                }
> +            }
> +        } else if (header->bitpix == 16) {
> +            for (i = 0; i < avctx->height; i++) {
> +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +                    t16 = AV_RB16(ptr8);
> +                    if (t16 != header->blank)
> +                        t16 = ((t16 - header->data_min) * 65535) /
> (header->data_max - header->data_min);
> +                    *dst16++ = t16;
> +                    ptr8 += 2;
> +                }
> +            }
> +        } else if (header->bitpix == 32) {
> +            for (i = 0; i < avctx->height; i++) {
> +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +                    t32 = AV_RB32(ptr8);
> +                    if (t32 != header->blank)
> +                        t16 = ((t32 - header->data_min) * 65535) /
> (header->data_max - header->data_min);
> +                    *dst16++ = t16;
> +                    ptr8 += 4;
> +                }
> +            }
> +        } else if (header->bitpix == 64) {
> +            for (i = 0; i < avctx->height; i++) {
> +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +                    t64 = AV_RB64(ptr8);
> +                    if (t64 != header->blank)
> +                        t16 = ((t64 - header->data_min) * 65535) /
> (header->data_max - header->data_min);
> +                    *dst16++ = t16;
> +                    ptr8 += 8;
> +                }
> +            }
> +        } else if (header->bitpix == -32) {
> +            for (i = 0; i < avctx->height; i++) {
> +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +                    t32 = AV_RB32(ptr8);
> +                    memcpy(&tflt, &t32, 4);
> +                    *dst16++ = ((tflt - header->data_min) * 65535) /
> (header->data_max - header->data_min);
> +                    ptr8 += 4;
> +                }
> +            }
> +        } else if (header->bitpix == -64) {
> +            for (i = 0; i < avctx->height; i++) {
> +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> p->linesize[0]);
> +                for (j = 0; j < avctx->width; j++) {
> +                    t64 = AV_RB64(ptr8);
> +                    memcpy(&tdbl, &t64, 8);
> +                    *dst16++ = ((tdbl - header->data_min) * 65535) /
> (header->data_max - header->data_min);
> +                    ptr8 += 8;
> +                }
> +            }
> +        } else {
> +            av_log(avctx, AV_LOG_ERROR, "invalid BITPIX, %d\n",
> header->bitpix);
> +            return AVERROR_INVALIDDATA;
> +        }
> +    }
> +
> +    p->key_frame = 1;
> +    p->pict_type = AV_PICTURE_TYPE_I;
> +
> +    *got_frame = 1;
> +
> +    return avpkt->size;
> +}
> +
> +AVCodec ff_fits_decoder = {
> +    .name           = "fits",
> +    .type           = AVMEDIA_TYPE_VIDEO,
> +    .id             = AV_CODEC_ID_FITS,
> +    .priv_data_size = sizeof(FITSDecContext),
> +    .decode         = fits_decode_frame,
> +    .capabilities   = AV_CODEC_CAP_DR1,
> +    .long_name      = NULL_IF_CONFIG_SMALL("Flexible Image Transport
> System")
> +};
> diff --git a/libavcodec/version.h b/libavcodec/version.h
> index c934872..a44a888 100644
> --- a/libavcodec/version.h
> +++ b/libavcodec/version.h
> @@ -28,7 +28,7 @@
>  #include "libavutil/version.h"
>
>  #define LIBAVCODEC_VERSION_MAJOR  57
> -#define LIBAVCODEC_VERSION_MINOR  98
> +#define LIBAVCODEC_VERSION_MINOR  99
>  #define LIBAVCODEC_VERSION_MICRO 100
>
>  #define LIBAVCODEC_VERSION_INT  AV_VERSION_INT(LIBAVCODEC_VERSION_MAJOR, \
> diff --git a/libavformat/img2.c b/libavformat/img2.c
> index 8432cc0..e405df8 100644
> --- a/libavformat/img2.c
> +++ b/libavformat/img2.c
> @@ -80,6 +80,7 @@ const IdStrMap ff_img_tags[] = {
>      { AV_CODEC_ID_XPM,        "xpm"      },
>      { AV_CODEC_ID_XFACE,      "xface"    },
>      { AV_CODEC_ID_XWD,        "xwd"      },
> +    { AV_CODEC_ID_FITS,       "fits"     },
>      { AV_CODEC_ID_NONE,       NULL       }
>  };
>
> --
> 2.4.11
>
> _______________________________________________
> ffmpeg-devel mailing list
> ffmpeg-devel@ffmpeg.org
> http://ffmpeg.org/mailman/listinfo/ffmpeg-devel
>
Paras June 14, 2017, 7:59 a.m. UTC | #2
On Tue, Jun 13, 2017 at 7:45 PM, Paul B Mahol <onemda@gmail.com> wrote:

> On 6/13/17, Paras Chadha <paraschadha18@gmail.com> wrote:
> > Above mentioned changes have been done.
> >
> > Signed-off-by: Paras Chadha <paraschadha18@gmail.com>
> > ---
> >  Changelog               |   1 +
> >  doc/general.texi        |   2 +
> >  libavcodec/Makefile     |   1 +
> >  libavcodec/allcodecs.c  |   1 +
> >  libavcodec/avcodec.h    |   1 +
> >  libavcodec/codec_desc.c |   8 +
> >  libavcodec/fitsdec.c    | 515
> > ++++++++++++++++++++++++++++++++++++++++++++++++
> >  libavcodec/version.h    |   2 +-
> >  libavformat/img2.c      |   1 +
> >  9 files changed, 531 insertions(+), 1 deletion(-)
> >  create mode 100644 libavcodec/fitsdec.c
> >
> > diff --git a/Changelog b/Changelog
> > index cd91f63..db495d9 100644
> > --- a/Changelog
> > +++ b/Changelog
> > @@ -20,6 +20,7 @@ version <next>:
> >  - sofalizer filter switched to libmysofa
> >  - Gremlin Digital Video demuxer and decoder
> >  - headphone audio filter
> > +- FITS decoder
> >
> >  version 3.3:
> >  - CrystalHD decoder moved to new decode API
> > diff --git a/doc/general.texi b/doc/general.texi
> > index 8f582d5..213e50e 100644
> > --- a/doc/general.texi
> > +++ b/doc/general.texi
> > @@ -591,6 +591,8 @@ following image formats are supported:
> >      @tab Digital Picture Exchange
> >  @item EXR          @tab   @tab X
> >      @tab OpenEXR
> > +@item FITS         @tab   @tab X
> > +    @tab FITS image format
> >  @item JPEG         @tab X @tab X
> >      @tab Progressive JPEG is not supported.
> >  @item JPEG 2000    @tab X @tab X
> > diff --git a/libavcodec/Makefile b/libavcodec/Makefile
> > index 2e7f19d..2cfd1fe 100644
> > --- a/libavcodec/Makefile
> > +++ b/libavcodec/Makefile
> > @@ -291,6 +291,7 @@ OBJS-$(CONFIG_FFV1_DECODER)            += ffv1dec.o
> > ffv1.o
> >  OBJS-$(CONFIG_FFV1_ENCODER)            += ffv1enc.o ffv1.o
> >  OBJS-$(CONFIG_FFWAVESYNTH_DECODER)     += ffwavesynth.o
> >  OBJS-$(CONFIG_FIC_DECODER)             += fic.o
> > +OBJS-$(CONFIG_FITS_DECODER)            += fitsdec.o
> >  OBJS-$(CONFIG_FLAC_DECODER)            += flacdec.o flacdata.o flac.o
> >  OBJS-$(CONFIG_FLAC_ENCODER)            += flacenc.o flacdata.o flac.o
> > vorbis_data.o
> >  OBJS-$(CONFIG_FLASHSV_DECODER)         += flashsv.o
> > diff --git a/libavcodec/allcodecs.c b/libavcodec/allcodecs.c
> > index 27110e1..8a89264 100644
> > --- a/libavcodec/allcodecs.c
> > +++ b/libavcodec/allcodecs.c
> > @@ -186,6 +186,7 @@ static void register_all(void)
> >      REGISTER_ENCDEC (FFV1,              ffv1);
> >      REGISTER_ENCDEC (FFVHUFF,           ffvhuff);
> >      REGISTER_DECODER(FIC,               fic);
> > +    REGISTER_DECODER(FITS,              fits);
> >      REGISTER_ENCDEC (FLASHSV,           flashsv);
> >      REGISTER_ENCDEC (FLASHSV2,          flashsv2);
> >      REGISTER_DECODER(FLIC,              flic);
> > diff --git a/libavcodec/avcodec.h b/libavcodec/avcodec.h
> > index dcdcfe0..e0c1089 100644
> > --- a/libavcodec/avcodec.h
> > +++ b/libavcodec/avcodec.h
> > @@ -447,6 +447,7 @@ enum AVCodecID {
> >      AV_CODEC_ID_SRGC,
> >      AV_CODEC_ID_SVG,
> >      AV_CODEC_ID_GDV,
> > +    AV_CODEC_ID_FITS,
> >
> >      /* various PCM "codecs" */
> >      AV_CODEC_ID_FIRST_AUDIO = 0x10000,     ///< A dummy id pointing at
> the
> > start of audio codecs
> > diff --git a/libavcodec/codec_desc.c b/libavcodec/codec_desc.c
> > index cf1246e..0112517 100644
> > --- a/libavcodec/codec_desc.c
> > +++ b/libavcodec/codec_desc.c
> > @@ -1464,6 +1464,14 @@ static const AVCodecDescriptor
> codec_descriptors[] =
> > {
> >                       AV_CODEC_PROP_LOSSLESS,
> >      },
> >      {
> > +        .id        = AV_CODEC_ID_FITS,
> > +        .type      = AVMEDIA_TYPE_VIDEO,
> > +        .name      = "fits",
> > +        .long_name = NULL_IF_CONFIG_SMALL("Flexible Image Transport
> > System"),
> > +        .props     = AV_CODEC_PROP_INTRA_ONLY | AV_CODEC_PROP_LOSSY |
> > +                     AV_CODEC_PROP_LOSSLESS,
> > +    },
> > +    {
> >          .id        = AV_CODEC_ID_GIF,
> >          .type      = AVMEDIA_TYPE_VIDEO,
> >          .name      = "gif",
> > diff --git a/libavcodec/fitsdec.c b/libavcodec/fitsdec.c
> > new file mode 100644
> > index 0000000..dfc4995
> > --- /dev/null
> > +++ b/libavcodec/fitsdec.c
> > @@ -0,0 +1,515 @@
> > +/*
> > + * FITS image decoder
> > + * Copyright (c) 2017 Paras Chadha
> > + *
> > + * This file is part of FFmpeg.
> > + *
> > + * FFmpeg is free software; you can redistribute it and/or
> > + * modify it under the terms of the GNU Lesser General Public
> > + * License as published by the Free Software Foundation; either
> > + * version 2.1 of the License, or (at your option) any later version.
> > + *
> > + * FFmpeg is distributed in the hope that it will be useful,
> > + * but WITHOUT ANY WARRANTY; without even the implied warranty of
> > + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
> > + * Lesser General Public License for more details.
> > + *
> > + * You should have received a copy of the GNU Lesser General Public
> > + * License along with FFmpeg; if not, write to the Free Software
> > + * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
> 02110-1301
> > USA
> > + */
> > +
> > +/**
> > + * @file
> > + * FITS image decoder
> > + * It supports all 2-d images alongwith, bzero, bscale and blank
> keywords.
> > + * RGBA images are supported as NAXIS3 = 3 or 4 i.e. Planes in RGBA
> order.
> > Also CTYPE = 'RGB ' should be present.
> > + * It currently does not support XTENSION keyword.
> > + * Also to interpret data, values are linearly scaled using min-max
> scaling
> > but not RGB images.
> > + */
> > +
> > +#include "avcodec.h"
> > +#include "internal.h"
> > +#include <float.h>
> > +#include "libavutil/intreadwrite.h"
> > +
> > +/**
> > + * Structure to store the header keywords in FITS file
> > + */
> > +typedef struct FITSContext {
> > +    char simple;
> > +    int bitpix;
> > +    int blank;
> > +    int naxis;
> > +    int naxisn[3];
> > +    int rgb; /**< 1 if file contains RGB image, 0 otherwise */
> > +    double bscale;
> > +    double bzero;
> > +    double data_min;
> > +    double data_max;
> > +} FITSDecContext;
> > +
> > +/**
> > + * function calculates the data_min and data_max values from the data.
> > + * This is called if the values are not present in the header.
> > + * @param ptr8 - pointer to the data
> > + * @param header - pointer to the header
> > + * @return 1, if calculated successfully, otherwise AVERROR_INVALIDDATA
> > + */
> > +static int fill_data_min_max(const uint8_t * ptr8, FITSDecContext *
> header,
> > const uint8_t * end)
> > +{
> > +    int16_t t16;
> > +    int32_t t32;
> > +    int64_t t64;
> > +    float tflt;
> > +    double tdbl;
> > +    int i, j;
> > +
> > +    header->data_min = DBL_MAX;
> > +    header->data_max = DBL_MIN;
> > +    switch (header->bitpix) {
> > +        case -64:
> > +            t64 = AV_RB64(ptr8);
> > +            memcpy(&tdbl, &t64, 8);
> > +            for (i = 0; i < header->naxisn[1]; i++) {
> > +                for (j = 0; j < header->naxisn[0]; j++) {
> > +                    t64 = AV_RB64(ptr8);
> > +                    memcpy(&tdbl, &t64, 8);
>
> Why are you using memcpy calls for this? Why not assign directly?
>

Actually the data portion of FITS file contains raw bytes. For double
values, it contains bits in IEEE double precision floating point format.
So, i am reading the number in a 64-bit integer and then interpreting the
bits as a double value. So, memcpy is used. If I assign it directly, it
will do typeconversion and convert the integer value to double value which
is wrong. So, to interpret the underlying bits as double, memcpy is used.


>
> > +                    if (tdbl > header->data_max)
> > +                        header->data_max = tdbl;
> > +                    if (tdbl < header->data_min)
> > +                        header->data_min = tdbl;
> > +                    ptr8 += 8;
> > +                }
> > +            }
> > +            break;
> > +        case -32:
> > +            t32 = AV_RB32(ptr8);
> > +            memcpy(&tflt, &t32, 4);
> > +            for (i = 0; i < header->naxisn[1]; i++) {
> > +                for (j = 0; j < header->naxisn[0]; j++) {
> > +                    t32 = AV_RB32(ptr8);
> > +                    memcpy(&tflt, &t32, 4);
> > +                    if (tflt > header->data_max)
> > +                        header->data_max = tflt;
> > +                    if (tflt < header->data_min)
> > +                        header->data_min = tflt;
> > +                    ptr8 += 4;
> > +                }
> > +            }
> > +            break;
> > +        case 8:
> > +            for (i = 0; i < header->naxisn[1]; i++) {
> > +                for (j = 0; j < header->naxisn[0]; j++) {
> > +                    if (ptr8[0] != header->blank) {
> > +                        if (ptr8[0] > header->data_max)
> > +                            header->data_max = ptr8[0];
> > +                        if (ptr8[0] < header->data_min)
> > +                            header->data_min = ptr8[0];
> > +                    }
> > +                    ptr8++;
> > +                }
> > +            }
> > +            break;
> > +        case 16:
> > +            t16 = AV_RB16(ptr8);
> > +            for (i = 0; i < header->naxisn[1]; i++) {
> > +                for (j = 0; j < header->naxisn[0]; j++) {
> > +                    t16 = AV_RB16(ptr8);
> > +                    if (t16 != header->blank) {
> > +                        if (t16 > header->data_max)
> > +                            header->data_max = t16;
> > +                        if (t16 < header->data_min)
> > +                            header->data_min = t16;
> > +                    }
> > +                    ptr8 += 2;
> > +                }
> > +            }
> > +            break;
> > +        case 32:
> > +            t32 = AV_RB32(ptr8);
> > +            for (i = 0; i < header->naxisn[1]; i++) {
> > +                for (j = 0; j < header->naxisn[0]; j++) {
> > +                    t32 = AV_RB32(ptr8);
> > +                    if (t32 != header->blank) {
> > +                        if (t32 > header->data_max)
> > +                            header->data_max = t32;
> > +                        if (t32 < header->data_min)
> > +                            header->data_min = t32;
> > +                    }
> > +                    ptr8 += 4;
> > +                }
> > +            }
> > +            break;
> > +        case 64:
> > +            t64 = AV_RB64(ptr8);
> > +            for (i = 0; i < header->naxisn[1]; i++) {
> > +                for (j = 0; j < header->naxisn[0]; j++) {
> > +                    t64 = AV_RB64(ptr8);
> > +                    if (t64 != header->blank) {
> > +                        if (t64 > header->data_max)
> > +                            header->data_max = t64;
> > +                        if (t64 < header->data_min)
> > +                            header->data_min = t64;
> > +                    }
> > +                    ptr8 += 8;
> > +                }
> > +            }
> > +            break;
> > +        default:
> > +            return AVERROR_INVALIDDATA;
> > +    }
> > +    return 1;
> > +}
> > +
> > +/**
> > + * function reads the fits header and stores the values in
> FITSDecContext
> > pointed by header
> > + * @param avctx - AVCodec context
> > + * @param ptr - pointer to pointer to the data
> > + * @param header - pointer to the FITSDecContext
> > + * @param end - pointer to end of packet
> > + * @return 1, if calculated successfully, otherwise AVERROR_INVALIDDATA
> > + */
> > +static int fits_read_header(AVCodecContext *avctx, const uint8_t **ptr,
> > FITSDecContext * header, const uint8_t * end)
> > +{
> > +    const uint8_t *ptr8 = *ptr;
> > +    int lines_read = 0, i, dim_no, t, data_min_found = 0,
> data_max_found =
> > 0, ret;
> > +    uint64_t size=1;
> > +    double d;
> > +
> > +    header->blank = 0;
> > +    header->bscale = 1.0;
> > +    header->bzero = 0;
> > +    header->rgb = 0;
> > +
> > +    if (end - ptr8 < 80)
> > +        return AVERROR_INVALIDDATA;
> > +
> > +    if (sscanf(ptr8, "SIMPLE = %c", &header->simple) != 1) {
> > +        av_log(avctx, AV_LOG_ERROR, "missing SIMPLE keyword\n");
> > +        return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    if (header->simple == 'F') {
> > +        av_log(avctx, AV_LOG_WARNING, "not a standard FITS file\n");
> > +    } else if (header->simple != 'T') {
> > +        av_log(avctx, AV_LOG_ERROR, "invalid SIMPLE value, SIMPLE =
> %c\n",
> > header->simple);
> > +        return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    ptr8 += 80;
> > +    lines_read++;
> > +
> > +    if (end - ptr8 < 80)
> > +        return AVERROR_INVALIDDATA;
> > +
> > +    if (sscanf(ptr8, "BITPIX = %d", &header->bitpix) != 1) {
> > +        av_log(avctx, AV_LOG_ERROR, "missing BITPIX keyword\n");
> > +        return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    size = abs(header->bitpix) >> 3;
> > +    ptr8 += 80;
> > +    lines_read++;
> > +
> > +    if (end - ptr8 < 80)
> > +        return AVERROR_INVALIDDATA;
> > +
> > +    if (sscanf(ptr8, "NAXIS = %d", &header->naxis) != 1) {
> > +        av_log(avctx, AV_LOG_ERROR, "missing NAXIS keyword\n");
> > +        return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    if (!header->naxis) {
> > +        av_log(avctx, AV_LOG_ERROR, "No image data found, NAXIS = %d\n",
> > header->naxis);
> > +        return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    if (header->naxis != 2 && header->naxis != 3) {
> > +        av_log(avctx, AV_LOG_ERROR, "unsupported number of dimensions,
> > NAXIS = %d\n", header->naxis);
> > +        return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    ptr8 += 80;
> > +    lines_read++;
> > +
> > +    for (i = 0; i < header->naxis; i++) {
> > +        if (end - ptr8 < 80)
> > +            return AVERROR_INVALIDDATA;
> > +
> > +        if (sscanf(ptr8, "NAXIS%d = %d", &dim_no, &header->naxisn[i])
> != 2
> > || dim_no != i+1) {
> > +            av_log(avctx, AV_LOG_ERROR, "missing NAXIS%d keyword\n",
> i+1);
> > +            return AVERROR_INVALIDDATA;
> > +        }
> > +
> > +        size *= header->naxisn[i];
> > +        ptr8 += 80;
> > +        lines_read++;
> > +    }
> > +
> > +    if (end - ptr8 < 80)
> > +            return AVERROR_INVALIDDATA;
> > +
> > +    while (strncmp(ptr8, "END", 3)) {
> > +        if (sscanf(ptr8, "BLANK = %d", &t) == 1) {
> > +            header->blank = t;
> > +        } else if (sscanf(ptr8, "BSCALE = %lf", &d) == 1) {
> > +            header->bscale = d;
> > +        } else if (sscanf(ptr8, "BZERO = %lf", &d) == 1) {
> > +            header->bzero = d;
> > +        } else if (sscanf(ptr8, "DATAMAX = %lf", &d) == 1) {
> > +            data_max_found = 1;
> > +            header->data_max = d;
> > +        } else if (sscanf(ptr8, "DATAMIN = %lf", &d) == 1) {
> > +            data_min_found = 1;
> > +            header->data_min = d;
> > +        } else if (!strncmp(ptr8, "CTYPE3  = 'RGB", 14)) {
> > +            header->rgb = 1;
> > +            if (header->naxis != 3 || (header->naxisn[2] != 3 &&
> > header->naxisn[2] != 4)) {
> > +                av_log(avctx, AV_LOG_ERROR, "File contains RGB image but
> > NAXIS = %d and NAXIS3 = %d\n", header->naxis, header->naxisn[2]);
> > +                return AVERROR_INVALIDDATA;
> > +            }
> > +        }
> > +        ptr8 += 80;
> > +        lines_read++;
> > +
> > +        if (end - ptr8 < 80)
> > +            return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    if (!header->rgb && header->naxis != 2){
> > +        av_log(avctx, AV_LOG_ERROR, "unsupported number of dimensions,
> > NAXIS = %d\n", header->naxis);
> > +        return AVERROR_INVALIDDATA;
> > +    }
> > +
> > +    ptr8 += 80;
> > +    lines_read++;
> > +    lines_read %= 36;
> > +
> > +    t = ((36 - lines_read) % 36) * 80;
> > +    if (end - ptr8 < t)
> > +        return AVERROR_INVALIDDATA;
> > +    ptr8 += t;
> > +    *ptr = ptr8;
> > +
> > +    if (end - ptr8 < size)
> > +        return AVERROR_INVALIDDATA;
> > +
> > +    if (!header->rgb && (!data_min_found || !data_max_found)) {
> > +        if ((ret = fill_data_min_max(ptr8, header, end)) < 0) {
> > +            av_log(avctx, AV_LOG_ERROR, "invalid BITPIX, %d\n",
> > header->bitpix);
> > +            return AVERROR_INVALIDDATA;
> > +        }
> > +    } else {
> > +        /*
> > +         * instead of applying bscale and bzero to every element, we
> can do
> > inverse transformation on data_min and
> > +         * data_max
> > +         */
> > +        header->data_min = (header->data_min - header->bzero) /
> > header->bscale;
> > +        header->data_max = (header->data_max - header->bzero) /
> > header->bscale;
> > +    }
> > +    return 1;
> > +}
> > +
> > +static int fits_decode_frame(AVCodecContext *avctx, void *data, int
> > *got_frame, AVPacket *avpkt)
> > +{
> > +    AVFrame *p=data;
> > +    const uint8_t *ptr8 = avpkt->data, *end;
> > +    int16_t t16;
> > +    int32_t t32;
> > +    int64_t t64;
> > +    float   tflt;
> > +    double  tdbl;
> > +    int ret, i, j;
> > +    uint8_t *dst8;
> > +    uint16_t *dst16;
> > +    uint32_t *dst32;
> > +    uint64_t *dst64, size, r, g, b, a, t;
> > +    FITSDecContext * header = avctx->priv_data;
> > +
> > +    end = ptr8 + avpkt->size;
> > +    if (ret = fits_read_header(avctx, &ptr8, header, end) < 0)
> > +        return ret;
> > +
> > +    size = (header->naxisn[0]) * (header->naxisn[1]);
> > +
> > +    if (header->rgb) {
> > +        if (header->bitpix == 8) {
> > +            avctx->pix_fmt = AV_PIX_FMT_RGB32;
> > +        } else if (header->bitpix == 16) {
> > +            avctx->pix_fmt = AV_PIX_FMT_RGBA64;
> > +        } else {
> > +            av_log(avctx, AV_LOG_ERROR, "unsupported BITPIX = %d\n",
> > header->bitpix);
> > +            return AVERROR_INVALIDDATA;
> > +        }
> > +    } else {
> > +        if (header->bitpix == 8) {
> > +            avctx->pix_fmt = AV_PIX_FMT_GRAY8;
> > +        } else {
> > +            avctx->pix_fmt = AV_PIX_FMT_GRAY16;
> > +        }
> > +    }
> > +
> > +    if ((ret = ff_set_dimensions(avctx, header->naxisn[0],
> > header->naxisn[1])) < 0)
> > +        return ret;
> > +
> > +    if ((ret = ff_get_buffer(avctx, p, 0)) < 0)
> > +        return ret;
> > +
> > +    if (header->rgb) {
> > +        if (header->bitpix == 8) {
> > +            for (i = 0; i < avctx->height; i++) {
> > +                /*
> > +                 * FITS stores images with bottom row first. Therefore
> we
> > have
> > +                 * to fill the image from bottom to top.
> > +                 */
> > +                dst32 = (uint32_t *)(p->data[0] + (avctx->height-i-1)*
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +                    if (header->naxisn[2] == 4) {
> > +                        if (ptr8[size * 3] != header->blank)
> > +                            t = ptr8[size * 3] * header->bscale +
> > header->bzero;
> > +                        a = t << 24;
> > +                    } else {
> > +                        a = (255 << 24);
> > +                    }
> > +
> > +                    if (ptr8[0] != header->blank)
> > +                        t = ptr8[0] * header->bscale + header->bzero;
> > +                    r = t << 16;
> > +
> > +                    if (ptr8[size] != header->blank)
> > +                        t = ptr8[size] * header->bscale + header->bzero;
> > +                    g = t << 8;
> > +
> > +                    if (ptr8[size * 2] != header->blank)
> > +                        t = ptr8[size * 2] * header->bscale +
> > header->bzero;
> > +                    b = t;
> > +
> > +                    *dst32++ = ((uint32_t)a) | ((uint32_t)r) |
> > ((uint32_t)g) | ((uint32_t)b);
> > +                    ptr8++;
> > +                }
> > +            }
> > +        } else if (header->bitpix == 16) {
> > +            // not tested ....
> > +            for (i = 0; i < avctx->height; i++) {
> > +                dst64 = (uint64_t *)(p->data[0] + (avctx->height-i-1) *
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +
> > +                    if (header->naxisn[2] == 4) {
> > +                        t = ((ptr8[size * 3] << 8) | ptr8[size * 3 +
> 1]);
> > +                        if (t != header->blank)
> > +                            t = t*header->bscale + header->bzero;
> > +                        a = t << 48;
> > +                    } else {
> > +                        a = 65535ULL << 48;
> > +                    }
> > +
> > +                    t = ptr8[0] << 8 | ptr8[1];
> > +                    if (t != header->blank)
> > +                        t = t*header->bscale + header->bzero;
> > +                    r = t << 32;
> > +
> > +                    t = ptr8[size] << 8 | ptr8[size + 1];
> > +                    if (t != header->blank)
> > +                        t = t*header->bscale + header->bzero;
> > +                    g = t << 16;
> > +
> > +                    t = ptr8[size * 2] << 8 | ptr8[size * 2 + 1];
> > +                    if (t != header->blank)
> > +                        t = t*header->bscale + header->bzero;
> > +                    b = t;
> > +
> > +                    *dst64++ = a | r | g | b;
> > +                    ptr8 += 2;
> > +                }
> > +            }
> > +        }
> > +    } else {
> > +        if (header->bitpix == 8) {
> > +            for (i = 0; i < avctx->height; i++) {
> > +                dst8 = (uint8_t *) (p->data[0] + (avctx->height-i-1)*
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +                    if (ptr8[0] != header->blank) {
> > +                        *dst8++ = ((ptr8[0] - header->data_min) * 255) /
> > (header->data_max - header->data_min);
> > +                    } else {
> > +                        *dst8++ = ptr8[0];
> > +                    }
> > +                    ptr8++;
> > +                }
> > +            }
> > +        } else if (header->bitpix == 16) {
> > +            for (i = 0; i < avctx->height; i++) {
> > +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +                    t16 = AV_RB16(ptr8);
> > +                    if (t16 != header->blank)
> > +                        t16 = ((t16 - header->data_min) * 65535) /
> > (header->data_max - header->data_min);
> > +                    *dst16++ = t16;
> > +                    ptr8 += 2;
> > +                }
> > +            }
> > +        } else if (header->bitpix == 32) {
> > +            for (i = 0; i < avctx->height; i++) {
> > +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +                    t32 = AV_RB32(ptr8);
> > +                    if (t32 != header->blank)
> > +                        t16 = ((t32 - header->data_min) * 65535) /
> > (header->data_max - header->data_min);
> > +                    *dst16++ = t16;
> > +                    ptr8 += 4;
> > +                }
> > +            }
> > +        } else if (header->bitpix == 64) {
> > +            for (i = 0; i < avctx->height; i++) {
> > +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +                    t64 = AV_RB64(ptr8);
> > +                    if (t64 != header->blank)
> > +                        t16 = ((t64 - header->data_min) * 65535) /
> > (header->data_max - header->data_min);
> > +                    *dst16++ = t16;
> > +                    ptr8 += 8;
> > +                }
> > +            }
> > +        } else if (header->bitpix == -32) {
> > +            for (i = 0; i < avctx->height; i++) {
> > +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +                    t32 = AV_RB32(ptr8);
> > +                    memcpy(&tflt, &t32, 4);
> > +                    *dst16++ = ((tflt - header->data_min) * 65535) /
> > (header->data_max - header->data_min);
> > +                    ptr8 += 4;
> > +                }
> > +            }
> > +        } else if (header->bitpix == -64) {
> > +            for (i = 0; i < avctx->height; i++) {
> > +                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) *
> > p->linesize[0]);
> > +                for (j = 0; j < avctx->width; j++) {
> > +                    t64 = AV_RB64(ptr8);
> > +                    memcpy(&tdbl, &t64, 8);
> > +                    *dst16++ = ((tdbl - header->data_min) * 65535) /
> > (header->data_max - header->data_min);
> > +                    ptr8 += 8;
> > +                }
> > +            }
> > +        } else {
> > +            av_log(avctx, AV_LOG_ERROR, "invalid BITPIX, %d\n",
> > header->bitpix);
> > +            return AVERROR_INVALIDDATA;
> > +        }
> > +    }
> > +
> > +    p->key_frame = 1;
> > +    p->pict_type = AV_PICTURE_TYPE_I;
> > +
> > +    *got_frame = 1;
> > +
> > +    return avpkt->size;
> > +}
> > +
> > +AVCodec ff_fits_decoder = {
> > +    .name           = "fits",
> > +    .type           = AVMEDIA_TYPE_VIDEO,
> > +    .id             = AV_CODEC_ID_FITS,
> > +    .priv_data_size = sizeof(FITSDecContext),
> > +    .decode         = fits_decode_frame,
> > +    .capabilities   = AV_CODEC_CAP_DR1,
> > +    .long_name      = NULL_IF_CONFIG_SMALL("Flexible Image Transport
> > System")
> > +};
> > diff --git a/libavcodec/version.h b/libavcodec/version.h
> > index c934872..a44a888 100644
> > --- a/libavcodec/version.h
> > +++ b/libavcodec/version.h
> > @@ -28,7 +28,7 @@
> >  #include "libavutil/version.h"
> >
> >  #define LIBAVCODEC_VERSION_MAJOR  57
> > -#define LIBAVCODEC_VERSION_MINOR  98
> > +#define LIBAVCODEC_VERSION_MINOR  99
> >  #define LIBAVCODEC_VERSION_MICRO 100
> >
> >  #define LIBAVCODEC_VERSION_INT  AV_VERSION_INT(LIBAVCODEC_VERSION_MAJOR,
> \
> > diff --git a/libavformat/img2.c b/libavformat/img2.c
> > index 8432cc0..e405df8 100644
> > --- a/libavformat/img2.c
> > +++ b/libavformat/img2.c
> > @@ -80,6 +80,7 @@ const IdStrMap ff_img_tags[] = {
> >      { AV_CODEC_ID_XPM,        "xpm"      },
> >      { AV_CODEC_ID_XFACE,      "xface"    },
> >      { AV_CODEC_ID_XWD,        "xwd"      },
> > +    { AV_CODEC_ID_FITS,       "fits"     },
> >      { AV_CODEC_ID_NONE,       NULL       }
> >  };
> >
> > --
> > 2.4.11
> >
> > _______________________________________________
> > ffmpeg-devel mailing list
> > ffmpeg-devel@ffmpeg.org
> > http://ffmpeg.org/mailman/listinfo/ffmpeg-devel
> >
>
diff mbox

Patch

diff --git a/Changelog b/Changelog
index cd91f63..db495d9 100644
--- a/Changelog
+++ b/Changelog
@@ -20,6 +20,7 @@  version <next>:
 - sofalizer filter switched to libmysofa
 - Gremlin Digital Video demuxer and decoder
 - headphone audio filter
+- FITS decoder

 version 3.3:
 - CrystalHD decoder moved to new decode API
diff --git a/doc/general.texi b/doc/general.texi
index 8f582d5..213e50e 100644
--- a/doc/general.texi
+++ b/doc/general.texi
@@ -591,6 +591,8 @@  following image formats are supported:
     @tab Digital Picture Exchange
 @item EXR          @tab   @tab X
     @tab OpenEXR
+@item FITS         @tab   @tab X
+    @tab FITS image format
 @item JPEG         @tab X @tab X
     @tab Progressive JPEG is not supported.
 @item JPEG 2000    @tab X @tab X
diff --git a/libavcodec/Makefile b/libavcodec/Makefile
index 2e7f19d..2cfd1fe 100644
--- a/libavcodec/Makefile
+++ b/libavcodec/Makefile
@@ -291,6 +291,7 @@  OBJS-$(CONFIG_FFV1_DECODER)            += ffv1dec.o ffv1.o
 OBJS-$(CONFIG_FFV1_ENCODER)            += ffv1enc.o ffv1.o
 OBJS-$(CONFIG_FFWAVESYNTH_DECODER)     += ffwavesynth.o
 OBJS-$(CONFIG_FIC_DECODER)             += fic.o
+OBJS-$(CONFIG_FITS_DECODER)            += fitsdec.o
 OBJS-$(CONFIG_FLAC_DECODER)            += flacdec.o flacdata.o flac.o
 OBJS-$(CONFIG_FLAC_ENCODER)            += flacenc.o flacdata.o flac.o vorbis_data.o
 OBJS-$(CONFIG_FLASHSV_DECODER)         += flashsv.o
diff --git a/libavcodec/allcodecs.c b/libavcodec/allcodecs.c
index 27110e1..8a89264 100644
--- a/libavcodec/allcodecs.c
+++ b/libavcodec/allcodecs.c
@@ -186,6 +186,7 @@  static void register_all(void)
     REGISTER_ENCDEC (FFV1,              ffv1);
     REGISTER_ENCDEC (FFVHUFF,           ffvhuff);
     REGISTER_DECODER(FIC,               fic);
+    REGISTER_DECODER(FITS,              fits);
     REGISTER_ENCDEC (FLASHSV,           flashsv);
     REGISTER_ENCDEC (FLASHSV2,          flashsv2);
     REGISTER_DECODER(FLIC,              flic);
diff --git a/libavcodec/avcodec.h b/libavcodec/avcodec.h
index dcdcfe0..e0c1089 100644
--- a/libavcodec/avcodec.h
+++ b/libavcodec/avcodec.h
@@ -447,6 +447,7 @@  enum AVCodecID {
     AV_CODEC_ID_SRGC,
     AV_CODEC_ID_SVG,
     AV_CODEC_ID_GDV,
+    AV_CODEC_ID_FITS,

     /* various PCM "codecs" */
     AV_CODEC_ID_FIRST_AUDIO = 0x10000,     ///< A dummy id pointing at the start of audio codecs
diff --git a/libavcodec/codec_desc.c b/libavcodec/codec_desc.c
index cf1246e..0112517 100644
--- a/libavcodec/codec_desc.c
+++ b/libavcodec/codec_desc.c
@@ -1464,6 +1464,14 @@  static const AVCodecDescriptor codec_descriptors[] = {
                      AV_CODEC_PROP_LOSSLESS,
     },
     {
+        .id        = AV_CODEC_ID_FITS,
+        .type      = AVMEDIA_TYPE_VIDEO,
+        .name      = "fits",
+        .long_name = NULL_IF_CONFIG_SMALL("Flexible Image Transport System"),
+        .props     = AV_CODEC_PROP_INTRA_ONLY | AV_CODEC_PROP_LOSSY |
+                     AV_CODEC_PROP_LOSSLESS,
+    },
+    {
         .id        = AV_CODEC_ID_GIF,
         .type      = AVMEDIA_TYPE_VIDEO,
         .name      = "gif",
diff --git a/libavcodec/fitsdec.c b/libavcodec/fitsdec.c
new file mode 100644
index 0000000..dfc4995
--- /dev/null
+++ b/libavcodec/fitsdec.c
@@ -0,0 +1,515 @@ 
+/*
+ * FITS image decoder
+ * Copyright (c) 2017 Paras Chadha
+ *
+ * This file is part of FFmpeg.
+ *
+ * FFmpeg is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public
+ * License as published by the Free Software Foundation; either
+ * version 2.1 of the License, or (at your option) any later version.
+ *
+ * FFmpeg is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with FFmpeg; if not, write to the Free Software
+ * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
+ */
+
+/**
+ * @file
+ * FITS image decoder
+ * It supports all 2-d images alongwith, bzero, bscale and blank keywords.
+ * RGBA images are supported as NAXIS3 = 3 or 4 i.e. Planes in RGBA order. Also CTYPE = 'RGB ' should be present.
+ * It currently does not support XTENSION keyword.
+ * Also to interpret data, values are linearly scaled using min-max scaling but not RGB images.
+ */
+
+#include "avcodec.h"
+#include "internal.h"
+#include <float.h>
+#include "libavutil/intreadwrite.h"
+
+/**
+ * Structure to store the header keywords in FITS file
+ */
+typedef struct FITSContext {
+    char simple;
+    int bitpix;
+    int blank;
+    int naxis;
+    int naxisn[3];
+    int rgb; /**< 1 if file contains RGB image, 0 otherwise */
+    double bscale;
+    double bzero;
+    double data_min;
+    double data_max;
+} FITSDecContext;
+
+/**
+ * function calculates the data_min and data_max values from the data.
+ * This is called if the values are not present in the header.
+ * @param ptr8 - pointer to the data
+ * @param header - pointer to the header
+ * @return 1, if calculated successfully, otherwise AVERROR_INVALIDDATA
+ */
+static int fill_data_min_max(const uint8_t * ptr8, FITSDecContext * header, const uint8_t * end)
+{
+    int16_t t16;
+    int32_t t32;
+    int64_t t64;
+    float tflt;
+    double tdbl;
+    int i, j;
+
+    header->data_min = DBL_MAX;
+    header->data_max = DBL_MIN;
+    switch (header->bitpix) {
+        case -64:
+            t64 = AV_RB64(ptr8);
+            memcpy(&tdbl, &t64, 8);
+            for (i = 0; i < header->naxisn[1]; i++) {
+                for (j = 0; j < header->naxisn[0]; j++) {
+                    t64 = AV_RB64(ptr8);
+                    memcpy(&tdbl, &t64, 8);
+                    if (tdbl > header->data_max)
+                        header->data_max = tdbl;
+                    if (tdbl < header->data_min)
+                        header->data_min = tdbl;
+                    ptr8 += 8;
+                }
+            }
+            break;
+        case -32:
+            t32 = AV_RB32(ptr8);
+            memcpy(&tflt, &t32, 4);
+            for (i = 0; i < header->naxisn[1]; i++) {
+                for (j = 0; j < header->naxisn[0]; j++) {
+                    t32 = AV_RB32(ptr8);
+                    memcpy(&tflt, &t32, 4);
+                    if (tflt > header->data_max)
+                        header->data_max = tflt;
+                    if (tflt < header->data_min)
+                        header->data_min = tflt;
+                    ptr8 += 4;
+                }
+            }
+            break;
+        case 8:
+            for (i = 0; i < header->naxisn[1]; i++) {
+                for (j = 0; j < header->naxisn[0]; j++) {
+                    if (ptr8[0] != header->blank) {
+                        if (ptr8[0] > header->data_max)
+                            header->data_max = ptr8[0];
+                        if (ptr8[0] < header->data_min)
+                            header->data_min = ptr8[0];
+                    }
+                    ptr8++;
+                }
+            }
+            break;
+        case 16:
+            t16 = AV_RB16(ptr8);
+            for (i = 0; i < header->naxisn[1]; i++) {
+                for (j = 0; j < header->naxisn[0]; j++) {
+                    t16 = AV_RB16(ptr8);
+                    if (t16 != header->blank) {
+                        if (t16 > header->data_max)
+                            header->data_max = t16;
+                        if (t16 < header->data_min)
+                            header->data_min = t16;
+                    }
+                    ptr8 += 2;
+                }
+            }
+            break;
+        case 32:
+            t32 = AV_RB32(ptr8);
+            for (i = 0; i < header->naxisn[1]; i++) {
+                for (j = 0; j < header->naxisn[0]; j++) {
+                    t32 = AV_RB32(ptr8);
+                    if (t32 != header->blank) {
+                        if (t32 > header->data_max)
+                            header->data_max = t32;
+                        if (t32 < header->data_min)
+                            header->data_min = t32;
+                    }
+                    ptr8 += 4;
+                }
+            }
+            break;
+        case 64:
+            t64 = AV_RB64(ptr8);
+            for (i = 0; i < header->naxisn[1]; i++) {
+                for (j = 0; j < header->naxisn[0]; j++) {
+                    t64 = AV_RB64(ptr8);
+                    if (t64 != header->blank) {
+                        if (t64 > header->data_max)
+                            header->data_max = t64;
+                        if (t64 < header->data_min)
+                            header->data_min = t64;
+                    }
+                    ptr8 += 8;
+                }
+            }
+            break;
+        default:
+            return AVERROR_INVALIDDATA;
+    }
+    return 1;
+}
+
+/**
+ * function reads the fits header and stores the values in FITSDecContext pointed by header
+ * @param avctx - AVCodec context
+ * @param ptr - pointer to pointer to the data
+ * @param header - pointer to the FITSDecContext
+ * @param end - pointer to end of packet
+ * @return 1, if calculated successfully, otherwise AVERROR_INVALIDDATA
+ */
+static int fits_read_header(AVCodecContext *avctx, const uint8_t **ptr, FITSDecContext * header, const uint8_t * end)
+{
+    const uint8_t *ptr8 = *ptr;
+    int lines_read = 0, i, dim_no, t, data_min_found = 0, data_max_found = 0, ret;
+    uint64_t size=1;
+    double d;
+
+    header->blank = 0;
+    header->bscale = 1.0;
+    header->bzero = 0;
+    header->rgb = 0;
+
+    if (end - ptr8 < 80)
+        return AVERROR_INVALIDDATA;
+
+    if (sscanf(ptr8, "SIMPLE = %c", &header->simple) != 1) {
+        av_log(avctx, AV_LOG_ERROR, "missing SIMPLE keyword\n");
+        return AVERROR_INVALIDDATA;
+    }
+
+    if (header->simple == 'F') {
+        av_log(avctx, AV_LOG_WARNING, "not a standard FITS file\n");
+    } else if (header->simple != 'T') {
+        av_log(avctx, AV_LOG_ERROR, "invalid SIMPLE value, SIMPLE = %c\n", header->simple);
+        return AVERROR_INVALIDDATA;
+    }
+
+    ptr8 += 80;
+    lines_read++;
+
+    if (end - ptr8 < 80)
+        return AVERROR_INVALIDDATA;
+
+    if (sscanf(ptr8, "BITPIX = %d", &header->bitpix) != 1) {
+        av_log(avctx, AV_LOG_ERROR, "missing BITPIX keyword\n");
+        return AVERROR_INVALIDDATA;
+    }
+
+    size = abs(header->bitpix) >> 3;
+    ptr8 += 80;
+    lines_read++;
+
+    if (end - ptr8 < 80)
+        return AVERROR_INVALIDDATA;
+
+    if (sscanf(ptr8, "NAXIS = %d", &header->naxis) != 1) {
+        av_log(avctx, AV_LOG_ERROR, "missing NAXIS keyword\n");
+        return AVERROR_INVALIDDATA;
+    }
+
+    if (!header->naxis) {
+        av_log(avctx, AV_LOG_ERROR, "No image data found, NAXIS = %d\n", header->naxis);
+        return AVERROR_INVALIDDATA;
+    }
+
+    if (header->naxis != 2 && header->naxis != 3) {
+        av_log(avctx, AV_LOG_ERROR, "unsupported number of dimensions, NAXIS = %d\n", header->naxis);
+        return AVERROR_INVALIDDATA;
+    }
+
+    ptr8 += 80;
+    lines_read++;
+
+    for (i = 0; i < header->naxis; i++) {
+        if (end - ptr8 < 80)
+            return AVERROR_INVALIDDATA;
+
+        if (sscanf(ptr8, "NAXIS%d = %d", &dim_no, &header->naxisn[i]) != 2 || dim_no != i+1) {
+            av_log(avctx, AV_LOG_ERROR, "missing NAXIS%d keyword\n", i+1);
+            return AVERROR_INVALIDDATA;
+        }
+
+        size *= header->naxisn[i];
+        ptr8 += 80;
+        lines_read++;
+    }
+
+    if (end - ptr8 < 80)
+            return AVERROR_INVALIDDATA;
+
+    while (strncmp(ptr8, "END", 3)) {
+        if (sscanf(ptr8, "BLANK = %d", &t) == 1) {
+            header->blank = t;
+        } else if (sscanf(ptr8, "BSCALE = %lf", &d) == 1) {
+            header->bscale = d;
+        } else if (sscanf(ptr8, "BZERO = %lf", &d) == 1) {
+            header->bzero = d;
+        } else if (sscanf(ptr8, "DATAMAX = %lf", &d) == 1) {
+            data_max_found = 1;
+            header->data_max = d;
+        } else if (sscanf(ptr8, "DATAMIN = %lf", &d) == 1) {
+            data_min_found = 1;
+            header->data_min = d;
+        } else if (!strncmp(ptr8, "CTYPE3  = 'RGB", 14)) {
+            header->rgb = 1;
+            if (header->naxis != 3 || (header->naxisn[2] != 3 && header->naxisn[2] != 4)) {
+                av_log(avctx, AV_LOG_ERROR, "File contains RGB image but NAXIS = %d and NAXIS3 = %d\n", header->naxis, header->naxisn[2]);
+                return AVERROR_INVALIDDATA;
+            }
+        }
+        ptr8 += 80;
+        lines_read++;
+
+        if (end - ptr8 < 80)
+            return AVERROR_INVALIDDATA;
+    }
+
+    if (!header->rgb && header->naxis != 2){
+        av_log(avctx, AV_LOG_ERROR, "unsupported number of dimensions, NAXIS = %d\n", header->naxis);
+        return AVERROR_INVALIDDATA;
+    }
+
+    ptr8 += 80;
+    lines_read++;
+    lines_read %= 36;
+
+    t = ((36 - lines_read) % 36) * 80;
+    if (end - ptr8 < t)
+        return AVERROR_INVALIDDATA;
+    ptr8 += t;
+    *ptr = ptr8;
+
+    if (end - ptr8 < size)
+        return AVERROR_INVALIDDATA;
+
+    if (!header->rgb && (!data_min_found || !data_max_found)) {
+        if ((ret = fill_data_min_max(ptr8, header, end)) < 0) {
+            av_log(avctx, AV_LOG_ERROR, "invalid BITPIX, %d\n", header->bitpix);
+            return AVERROR_INVALIDDATA;
+        }
+    } else {
+        /*
+         * instead of applying bscale and bzero to every element, we can do inverse transformation on data_min and
+         * data_max
+         */
+        header->data_min = (header->data_min - header->bzero) / header->bscale;
+        header->data_max = (header->data_max - header->bzero) / header->bscale;
+    }
+    return 1;
+}
+
+static int fits_decode_frame(AVCodecContext *avctx, void *data, int *got_frame, AVPacket *avpkt)
+{
+    AVFrame *p=data;
+    const uint8_t *ptr8 = avpkt->data, *end;
+    int16_t t16;
+    int32_t t32;
+    int64_t t64;
+    float   tflt;
+    double  tdbl;
+    int ret, i, j;
+    uint8_t *dst8;
+    uint16_t *dst16;
+    uint32_t *dst32;
+    uint64_t *dst64, size, r, g, b, a, t;
+    FITSDecContext * header = avctx->priv_data;
+
+    end = ptr8 + avpkt->size;
+    if (ret = fits_read_header(avctx, &ptr8, header, end) < 0)
+        return ret;
+
+    size = (header->naxisn[0]) * (header->naxisn[1]);
+
+    if (header->rgb) {
+        if (header->bitpix == 8) {
+            avctx->pix_fmt = AV_PIX_FMT_RGB32;
+        } else if (header->bitpix == 16) {
+            avctx->pix_fmt = AV_PIX_FMT_RGBA64;
+        } else {
+            av_log(avctx, AV_LOG_ERROR, "unsupported BITPIX = %d\n", header->bitpix);
+            return AVERROR_INVALIDDATA;
+        }
+    } else {
+        if (header->bitpix == 8) {
+            avctx->pix_fmt = AV_PIX_FMT_GRAY8;
+        } else {
+            avctx->pix_fmt = AV_PIX_FMT_GRAY16;
+        }
+    }
+
+    if ((ret = ff_set_dimensions(avctx, header->naxisn[0], header->naxisn[1])) < 0)
+        return ret;
+
+    if ((ret = ff_get_buffer(avctx, p, 0)) < 0)
+        return ret;
+
+    if (header->rgb) {
+        if (header->bitpix == 8) {
+            for (i = 0; i < avctx->height; i++) {
+                /*
+                 * FITS stores images with bottom row first. Therefore we have
+                 * to fill the image from bottom to top.
+                 */
+                dst32 = (uint32_t *)(p->data[0] + (avctx->height-i-1)* p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+                    if (header->naxisn[2] == 4) {
+                        if (ptr8[size * 3] != header->blank)
+                            t = ptr8[size * 3] * header->bscale + header->bzero;
+                        a = t << 24;
+                    } else {
+                        a = (255 << 24);
+                    }
+
+                    if (ptr8[0] != header->blank)
+                        t = ptr8[0] * header->bscale + header->bzero;
+                    r = t << 16;
+
+                    if (ptr8[size] != header->blank)
+                        t = ptr8[size] * header->bscale + header->bzero;
+                    g = t << 8;
+
+                    if (ptr8[size * 2] != header->blank)
+                        t = ptr8[size * 2] * header->bscale + header->bzero;
+                    b = t;
+
+                    *dst32++ = ((uint32_t)a) | ((uint32_t)r) | ((uint32_t)g) | ((uint32_t)b);
+                    ptr8++;
+                }
+            }
+        } else if (header->bitpix == 16) {
+            // not tested ....
+            for (i = 0; i < avctx->height; i++) {
+                dst64 = (uint64_t *)(p->data[0] + (avctx->height-i-1) * p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+
+                    if (header->naxisn[2] == 4) {
+                        t = ((ptr8[size * 3] << 8) | ptr8[size * 3 + 1]);
+                        if (t != header->blank)
+                            t = t*header->bscale + header->bzero;
+                        a = t << 48;
+                    } else {
+                        a = 65535ULL << 48;
+                    }
+
+                    t = ptr8[0] << 8 | ptr8[1];
+                    if (t != header->blank)
+                        t = t*header->bscale + header->bzero;
+                    r = t << 32;
+
+                    t = ptr8[size] << 8 | ptr8[size + 1];
+                    if (t != header->blank)
+                        t = t*header->bscale + header->bzero;
+                    g = t << 16;
+
+                    t = ptr8[size * 2] << 8 | ptr8[size * 2 + 1];
+                    if (t != header->blank)
+                        t = t*header->bscale + header->bzero;
+                    b = t;
+
+                    *dst64++ = a | r | g | b;
+                    ptr8 += 2;
+                }
+            }
+        }
+    } else {
+        if (header->bitpix == 8) {
+            for (i = 0; i < avctx->height; i++) {
+                dst8 = (uint8_t *) (p->data[0] + (avctx->height-i-1)* p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+                    if (ptr8[0] != header->blank) {
+                        *dst8++ = ((ptr8[0] - header->data_min) * 255) / (header->data_max - header->data_min);
+                    } else {
+                        *dst8++ = ptr8[0];
+                    }
+                    ptr8++;
+                }
+            }
+        } else if (header->bitpix == 16) {
+            for (i = 0; i < avctx->height; i++) {
+                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) * p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+                    t16 = AV_RB16(ptr8);
+                    if (t16 != header->blank)
+                        t16 = ((t16 - header->data_min) * 65535) / (header->data_max - header->data_min);
+                    *dst16++ = t16;
+                    ptr8 += 2;
+                }
+            }
+        } else if (header->bitpix == 32) {
+            for (i = 0; i < avctx->height; i++) {
+                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) * p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+                    t32 = AV_RB32(ptr8);
+                    if (t32 != header->blank)
+                        t16 = ((t32 - header->data_min) * 65535) / (header->data_max - header->data_min);
+                    *dst16++ = t16;
+                    ptr8 += 4;
+                }
+            }
+        } else if (header->bitpix == 64) {
+            for (i = 0; i < avctx->height; i++) {
+                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) * p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+                    t64 = AV_RB64(ptr8);
+                    if (t64 != header->blank)
+                        t16 = ((t64 - header->data_min) * 65535) / (header->data_max - header->data_min);
+                    *dst16++ = t16;
+                    ptr8 += 8;
+                }
+            }
+        } else if (header->bitpix == -32) {
+            for (i = 0; i < avctx->height; i++) {
+                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) * p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+                    t32 = AV_RB32(ptr8);
+                    memcpy(&tflt, &t32, 4);
+                    *dst16++ = ((tflt - header->data_min) * 65535) / (header->data_max - header->data_min);
+                    ptr8 += 4;
+                }
+            }
+        } else if (header->bitpix == -64) {
+            for (i = 0; i < avctx->height; i++) {
+                dst16 = (uint16_t *)(p->data[0] + (avctx->height-i-1) * p->linesize[0]);
+                for (j = 0; j < avctx->width; j++) {
+                    t64 = AV_RB64(ptr8);
+                    memcpy(&tdbl, &t64, 8);
+                    *dst16++ = ((tdbl - header->data_min) * 65535) / (header->data_max - header->data_min);
+                    ptr8 += 8;
+                }
+            }
+        } else {
+            av_log(avctx, AV_LOG_ERROR, "invalid BITPIX, %d\n", header->bitpix);
+            return AVERROR_INVALIDDATA;
+        }
+    }
+
+    p->key_frame = 1;
+    p->pict_type = AV_PICTURE_TYPE_I;
+
+    *got_frame = 1;
+
+    return avpkt->size;
+}
+
+AVCodec ff_fits_decoder = {
+    .name           = "fits",
+    .type           = AVMEDIA_TYPE_VIDEO,
+    .id             = AV_CODEC_ID_FITS,
+    .priv_data_size = sizeof(FITSDecContext),
+    .decode         = fits_decode_frame,
+    .capabilities   = AV_CODEC_CAP_DR1,
+    .long_name      = NULL_IF_CONFIG_SMALL("Flexible Image Transport System")
+};
diff --git a/libavcodec/version.h b/libavcodec/version.h
index c934872..a44a888 100644
--- a/libavcodec/version.h
+++ b/libavcodec/version.h
@@ -28,7 +28,7 @@ 
 #include "libavutil/version.h"

 #define LIBAVCODEC_VERSION_MAJOR  57
-#define LIBAVCODEC_VERSION_MINOR  98
+#define LIBAVCODEC_VERSION_MINOR  99
 #define LIBAVCODEC_VERSION_MICRO 100

 #define LIBAVCODEC_VERSION_INT  AV_VERSION_INT(LIBAVCODEC_VERSION_MAJOR, \
diff --git a/libavformat/img2.c b/libavformat/img2.c
index 8432cc0..e405df8 100644
--- a/libavformat/img2.c
+++ b/libavformat/img2.c
@@ -80,6 +80,7 @@  const IdStrMap ff_img_tags[] = {
     { AV_CODEC_ID_XPM,        "xpm"      },
     { AV_CODEC_ID_XFACE,      "xface"    },
     { AV_CODEC_ID_XWD,        "xwd"      },
+    { AV_CODEC_ID_FITS,       "fits"     },
     { AV_CODEC_ID_NONE,       NULL       }
 };