rhkspm32topng.C

Go to the documentation of this file.

/* Gxsm - Gnome X Scanning Microscopy
 * universal STM/AFM/SARLS/SPALEED/... controlling and
 * data analysis software
 * 
 * Copyright (C) 1999,2000,2001 Percy Zahl
 *
 * Authors: Percy Zahl <zahl@users.sf.net>
 * additional features: Andreas Klust <klust@users.sf.net>
 * WWW Home: http://gxsm.sf.net
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 *
 * This program 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 General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307, USA.
 */

/* -*- Mode: C++; indent-tabs-mode: nil; c-basic-offset: 8 c-style: "K&R" -*- */

#include <cstdio>
#include <cstdlib>
#include <fstream>
#include <cmath>

#include <glib.h>
#include <netcdf.hh>
#include <png.h>
#include <popt.h>
#include "config.h"
#include "thumbpal.h"
#include "writepng.h"

using namespace std;

// ==============================================================
//
// based on the nctopng code
//
// Read SPM32 image/data generated by STM-1000 systems from 
// RHK Technology Inc.
// and convert them to png image/thumbnail
// - read only required subgrid data (scaled down)
// - use quick, "more intelligent" autoscale
// - use palette
//
// ==============================================================

#ifndef WORDS_BIGENDIAN
# define WORDS_BIGENDIAN 0
#endif

typedef enum RHK_SPM32_STATUS { SPM32_READ_OK, SPM32_OPEN_FAILED, SPM32_FILE_NOT_VALID };

#define THUMB_X 96 // Thumbnail max X size
#define THUMB_Y 91 // Thumbnail max Y size

class raw_image{
public:
        raw_image(){ rgb=NULL; };
        virtual ~raw_image(){ 
                if(rgb){
                        for (int i=0; i<ny; delete[] rgb[i++]);
                        delete[] rgb; 
                }
        };

        int width() { return nx; };
        int height() { return ny; };
        void calc_line_regress(int y, double &slope, double &offset){
                //
                // OutputFeld[i] = InputFeld[i] - (slope * i + offset)
                //
                double sumi, sumy, sumiy, sumyq, sumiq;
                int i, istep;
                int Xmax, Xmin;
                double a, n, imean, ymean;
                        
                sumi = sumiq = sumy = sumyq = sumiy = 0.0;
                Xmax = nx-nx/10; // etwas Abstand (10%) vom Rand !
                Xmin = nx/10;
                istep = 1;
                for (i=Xmin; i < Xmax; i+=istep) { /* Rev 2 */
                        sumi   += (double)i;
                        sumiq  += (double)i*(double)i;
                        a =  rowdata[y][i];
                        sumy   += a;
                        sumyq  += a * a;
                        sumiy  += a * (double)i;
                }
                n = (double)((Xmax-Xmin)/istep);
                if ((Xmax-Xmin) % istep != 0) n += 1.0;
                imean = sumi / n;
                ymean = sumy / n;
                slope  = (sumiy- n * imean * ymean ) / (sumiq - n * imean * imean);
                offset = ymean - slope * imean;
        };
        double soft(int i, int j, double lim=0){
/*
                int a,b,c,d;
                a = max (0, i-1);
                b = min (i+1, nx-1);
                c = max (0, j-1);
                d = min (j+1, ny-1);
*/
                return rowdata[j][i];
        }
        void find_soft_min_max (double &min, double &max){
                min = max = soft(0,0);
                for (int i=0; i<ny; ++i)
                        for (int j=0; j<nx; ++j){
                                double v = soft (j,i);
                                if (min > v) min = v;
                                if (max < v) max = v;
                        }
        };
        void quick_rgb(int linreg=TRUE) {
                        if (rgb){
                                for (int i=0; i<ny; delete[] rgb[i++]);
                                delete[] rgb;
                        }
                        rgb = new unsigned char* [ny];
                        for (int i=0; i<ny; rgb[i++] = new unsigned char[3*nx]);

                        if (linreg){
                                        // Calc and Apply "Quick"
                                        for (int i=0; i<ny; ++i){
                                                        double a,b;
                                                        calc_line_regress (i, a, b);
                                                        for (int j=0; j<nx; ++j)
                                                                        rowdata[i][j] -= j*a+b;
                                        }
                        }

                        // Find Scale
                        double min, max, range;
                        find_soft_min_max (min, max);
                        range = max-min;
                        for (int i=0; i<ny; ++i)
                                for (int j=0, k=0; j<nx; ++j){
                                        int idx = 3*(int)((rowdata[i][j]-min)/range*PALETTE_ENTRIES+0.5);
                                        rgb[i][k++] = thumb_palette[idx++];
                                        rgb[i][k++] = thumb_palette[idx++];
                                        rgb[i][k++] = thumb_palette[idx];
                                }
                        
        };
        unsigned char **row_rgb_pointers() { return rgb; };

protected:
        int nx, ny; // image size
        int x0, y0, w0; // offset, data width
        int onx, ony; // original data size of NC-file
        double **rowdata;
        unsigned char **rgb;
};



class spm32_image : public raw_image {
public:
        spm32_image (const gchar *file_name, int thumb, int new_x, int x_off, int y_off, int width) {

                char header[512];
                long dummy, type, data_offset, size;

                x0 = x_off; y0 = y_off;

                // read file header
                std::ifstream f(file_name, ios::in);
                if (!f.good()){
                        state = SPM32_OPEN_FAILED;
                        return;
                }

                // read header
                f.read(header, 512);

                // check file type
                // check for Magic No./ID at file header start:
                if (strncmp (header, "STiMage 3.1", 11) != 0){
                        f.close ();
                        state = SPM32_FILE_NOT_VALID;
                        return;
                }

                // need: onx, ony, data_offset
                sscanf (header + 0x20, "%ld %ld %ld %d %d %ld %ld",
                        &type, &dummy, &dummy, &onx, &ony, &size, &dummy);
                sscanf (header + 0x100, "id %ld %ld",
                        &dummy, &data_offset);

                // check if header information valid
                if ( type != 0 || size != onx*ony*2 ){
                        f.close ();
                        state = SPM32_FILE_NOT_VALID;
                        return;
                }
                
                //
                // compute scale and offsets
                //
                if (width>0){
                        w0=width;
                        nx=new_x;
                        ny=new_x;
                }else{
                        w0=onx;
                        if (thumb){
                                if(onx < ony*THUMB_X/THUMB_Y){
                                        ny=THUMB_Y;
                                        nx=onx*ny/ony;
                                }else{
                                        nx=THUMB_X;
                                        ny=ony*nx/onx;
                                }
                        }else{
                                if (new_x){ // scale to new X w aspect
                                        nx=new_x;
                                        ny=ony*nx/onx;
                                }else{
                                        nx=onx; 
                                        ny=ony;
                                }
                        }
                }
                if (x0+nx >= onx) x0=0; // fallback
                if (y0+ny >= ony) y0=0; // fallback


                //prepare memory
                rowdata = new double* [ny];
                for (int i=0; i<ny; rowdata[i++] = new double [nx]);


                //
                // read data
                //
        
                double scale = (double)nx/(double)w0;
                short *row = new short[onx];
        
                union order {short s; char c[2];};
                order swap;
                char low;
        
                for (int y=0; y<ny; y++) {
                        f.seekg(data_offset + (y0+(int)((double)y/scale+.5))*onx*sizeof(short), ios::beg);
                        f.read((char*)row, onx*sizeof(short));
        
                        for (int i=0; i<nx; ++i) {
                                swap.s = row[x0+(int)((double)i/scale+.5)];
                                
                                if (WORDS_BIGENDIAN) {
                                        low = swap.c[0];
                                        swap.c[0] = swap.c[1];
                                        swap.c[1] =low;
                                }

                                rowdata[y][i] = (double)swap.s;
                        }
                }

                f.close();

                state = SPM32_READ_OK;
        };

        const RHK_SPM32_STATUS status(){
                return state;
        };
        

        virtual ~spm32_image() {
                for (int i=0; i<ny; delete [] rowdata[i++]);
                delete [] rowdata;
        };

private:
        RHK_SPM32_STATUS state;
};


int write_png(gchar *file_name, raw_image *img){
// see here for pnglib docu:
// http://www.libpng.org/pub/png/libpng-manual.html

        mainprog_info m;

        m.gamma   = 1.;
        m.width   = img->width();
        m.height  = img->height();
        m.modtime = time(NULL);
        m.infile  = NULL;
        if (! (m.outfile = fopen(file_name, "wb")))
                return -1;
        m.row_pointers = img->row_rgb_pointers();
        m.title  = "rhkspm32topng";
        m.author = "A. Klust";
        m.desc   = "RHK SPM32 data to png";
        m.copyright = "GPL";
        m.email   = "klust@users.sourceforge.net";
        m.url     = "http://gxsm.sf.net";
        m.have_bg   = FALSE;
        m.have_time = FALSE;
        m.have_text = FALSE;
        m.pnmtype   = 6; // TYPE_RGB
        m.sample_depth = 8;
        m.interlaced= FALSE;

        writepng_init (&m);
        writepng_encode_image (&m);
        writepng_encode_finish (&m);
        writepng_cleanup (&m);

        return 0;
}


/* ****************************
              main
   ***************************/
int main(int argc, const char *argv[]) {
        int thumb = 1;
        int new_x = 0;
        int x_off = 0, y_off = 0, width = 0;    
        int verbose = 0;
        int noquick = 0;
        int help = 0;
        gchar *filename;
        gchar *destinationfilename;
        spm32_image *img = NULL;
        poptContext optCon; 

        struct poptOption optionsTable[] = {
                {"x-offset",'x',POPT_ARG_INT,&x_off,0,"x-offset",NULL},
                {"y-offset",'y',POPT_ARG_INT,&y_off,0,"y-offset",NULL},
                {"width",'w',POPT_ARG_INT,&width,0,"width",NULL},
                {"size",'s',POPT_ARG_INT,&new_x,0,"Size in x-direction.",NULL},
                {"verbose",'v',POPT_ARG_NONE,&verbose,1,"Display information.",NULL},
                {"noquick",'n',POPT_ARG_NONE,&noquick,1,"No-quick",NULL},
                {"help",'h',POPT_ARG_NONE,&help,1,"Print help.",NULL},
                { NULL, 0, 0, NULL, 0 }
        };

        optCon = poptGetContext(NULL, argc, argv, optionsTable, 0);
        poptSetOtherOptionHelp(optCon, "spm32file.sm2? [outfile.png]");
        while (poptGetNextOpt(optCon) >= 0) {;} //Now parse until no more options.

        if ((argc < 2 )||(help)) { 
                poptPrintHelp(optCon, stderr, 0);
                exit(1);
        }
                
        filename = g_strdup(poptGetArg(optCon));
        destinationfilename = g_strdup(poptGetArg(optCon));

        if (destinationfilename == NULL){
                destinationfilename = g_strjoin(NULL, filename, ".png", NULL);
                // using simple join. if you need more sophisticated
                // have a look at 'mmv' for suffix handling.
        }

        if(verbose){
                if (new_x == 0) 
                        std::cout << "Thumbnail-size" << std::endl;
                else
                        std::cout << "Rescaling to new x-size = " << new_x << std::endl;
                std::cout << "Sourcefile: " << filename << std::endl;
                std::cout << "Destinationfile: " << destinationfilename << std::endl;
        }

        if (new_x > 0)
                thumb = 0;
        
        if ( x_off + y_off + width > 0){
                if ( (x_off==0) || (y_off==0) || (width==0) ) {
                        std::cout << "Please use -x and -y and -w together."<< std::endl;
                        exit(-1);
                }
                if(verbose){
                        std::cout << "Offset: " << x_off << " x " << y_off << std::endl;
                        std::cout << "Width set to: " << width << std::endl;
                }
        }
                
        // load & rescale image
        img = new spm32_image(filename, thumb, new_x, x_off, y_off, width);

        switch (img->status()){
                case SPM32_READ_OK: break;
                case SPM32_OPEN_FAILED: 
                        std::cerr << "Error opening SPM32 file >" << filename << "<" << std::endl; 
                        exit(-1);
                        break;
                case SPM32_FILE_NOT_VALID:
                        std::cerr << "invalid SPM32 file >" << filename << "<" << std::endl; 
                        exit(-1);
                        break;
        }
                
        if  (!img){
                std::cerr << "Error while creating image from SPM32 file." << std::endl;
                exit(-1);
        }
                
        if (verbose)
                std::cout << "Converting ..." << std::endl; 
                
        if (noquick)
                img->quick_rgb(FALSE);
        else
                img->quick_rgb();

        if (verbose)
                std::cout << "Writing >" << destinationfilename << "<" << std::endl; 
                
        write_png(destinationfilename, img);

        if(verbose)
                std::cout << "Writing complete." << std::endl;
                
        g_free(filename);
        g_free(destinationfilename);
        poptFreeContext(optCon);
        exit(0);
}

---