Load an image to a GSL matrix

Question

Does anyone know about some function that loads a grayscale image into a GSL matrix?

Something like:

gsl_matrix *M;
load_image("./image.jpg", M); // any image extension would also be fine

Answer 1

Something like this works on my PC: allows you to load a grayscale JPG image using libjpeg.

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <jpeglib.h>
#include <gsl/gsl_matrix.h>

gsl_matrix *load_jpg_image(const char *pFileName)
{
    FILE *pFile;
    long jpegSize;
    unsigned char *pJpegBytes, *pPixels;
    struct jpeg_decompress_struct cinfo;
    struct jpeg_error_mgr jerr;
    int status, w, h, numComponents, stride, x, y;
    gsl_matrix *pMatrix;

    pFile = fopen(pFileName, "rb");
    if (!pFile)
    {
        fprintf(stderr, "Can't open file\n");
        return NULL;
    }

    // Get the size of the file
    fseek(pFile, 0, SEEK_END);
    jpegSize = ftell(pFile);
    rewind(pFile);

    if (jpegSize == 0)
    {
        fclose(pFile);
        fprintf(stderr, "Empty file\n");
        return NULL;
    }

    // Read it into memory and close file
    pJpegBytes = (unsigned char *)malloc(jpegSize);
    fread(pJpegBytes, 1, jpegSize, pFile);
    fclose(pFile);

    // Jpeg decompression starts here
    memset(&cinfo, 0, sizeof(struct jpeg_decompress_struct));
    cinfo.err = jpeg_std_error(&jerr);
    jpeg_create_decompress(&cinfo);

    jpeg_mem_src(&cinfo, pJpegBytes, jpegSize);
    status = jpeg_read_header(&cinfo, TRUE);

    if (status != 1)
    {
        free(pJpegBytes);
        fprintf(stderr, "Invalid JPEG header\n");
        return NULL;
    }

    jpeg_start_decompress(&cinfo);

    w = cinfo.output_width;
    h = cinfo.output_height;
    numComponents = cinfo.output_components;

    if (numComponents != 1)
    {
        free(pJpegBytes);
        fprintf(stderr, "Can only handle 1 color component\n");
        return NULL;
    }

    pPixels = (unsigned char *)malloc(w*h);
    stride = w*numComponents;

    // perhaps this can de done much faster by processing
    // multiple lines at once
    while (cinfo.output_scanline < cinfo.output_height) 
    {
        unsigned char *buffer_array[1];
        buffer_array[0] = pPixels + cinfo.output_scanline * stride;

        jpeg_read_scanlines(&cinfo, buffer_array, 1);
    }

    jpeg_finish_decompress(&cinfo);
    jpeg_destroy_decompress(&cinfo);

    free(pJpegBytes);

    // Now, create and fill in the matrix
    pMatrix = gsl_matrix_alloc(h, w);    
    for (y = 0 ; y < h ; y++)
        for (x = 0 ; x < w ; x++)
            gsl_matrix_set(pMatrix, y, x, pPixels[x+y*stride]);

    return pMatrix;
}

int main(void)
{
    gsl_matrix *pMatrix;
    int rows, cols;
    int i, j;

    pMatrix = load_jpg_image("test.jpg");
    if (pMatrix == NULL)
    {
        fprintf(stderr, "Can't load matrix\n");
        return -1;
    }

    //
    // Use the matrix
    // 

    gsl_matrix_free(pMatrix);
    return 0;
}