source: fact/tools/pyscripts/sandbox/dneise/fact_compress/c++/readfits.cpp@ 14259

// reading a complete binary file
#include <iostream>
#include <fstream>
#include <vector>

// I will use the fits class by TB for reading the interesting data from the files
#define HAVE_ZLIB
#include "fits"

using namespace std;

ifstream::pos_type size;
char * memblock;

int main (int argc, char * argv[]) {
    
    if (argc < 3)
    {
        cout << "Usage: " << argv[0] << " data-file-name calib-file-name [ouput-file-name]" << endl;
        cout << "" << endl;
    }
    
    char * data_file_name = argv[1];
    char * calib_file_name = argv[2];
    char * out_file_name = 0;
    if (argc == 4)
    {
        out_file_name = argv[3];
    }
    else
    {
        out_file_name = new char[strlen(data_file_name)+4];
        strcpy(out_file_name, data_file_name);
        strcpy(out_file_name+strlen(data_file_name),".fc");
        out_file_name[strlen(data_file_name)+3] = 0;
    }
    
    
    
    
    //=========================================================================
    // CALIBRATION CONSTANTS
    //=========================================================================
    fits * calib =new fits(calib_file_name);
    calib->PrintKeys();
    calib->PrintColums();
    offset_mV = new float[]
    
    // I need the offset calibration constants from the calibration files
    // they are stored as floats at a known position inside the calibration file...
    long position_of_calibration_constants = 0x1000;
    // The calibration constants are stored in units of pseudo-mV. 
    // but I need them in ADC units
    // I will first read offset_mV from the file, then multiply (or divide) 
    // with the conversion factor 2000/4096.;
    // then I will convert it to shorts.
    // From that point on I will only use the shorts in *offset*, so I can
    // free the memory for *offset_mV* already.
    int size_of_offset = 1440*1024;
    int size_of_offset_memblock = 1440*1024*sizeof(float);
    float * offset_mV = new float[size_of_offset];
    short * offset = new short[size_of_offset];
    char * memblock = new char[size_of_offset_memblock];
    
    ifstream calib (calib_file_name, ios::in|ios::binary);
    if ( !calib.is_open() )
    {
        cerr << "Could not open Calibration File:" << calib_file_name << ".. ABORT." << endl;
        return 1;
    }
    else // file was opened, I can go on...
    {
        calib.seekg(position_of_calibration_constants, ios::beg);
        calib.read(memblock, size_of_offset_memblock);
        offset_mV = (float*)memblock;
        for (int i=0; i<size_of_offset; i++)
        {
            // -0.5 is for rounding correctly negative integers.
            // in all cases where it worked, the offset should be negative for FACT.
            offset[i] = short(offset_mV / 2000. * 4096 - 0.5);
        }
    }
    delete[] offset_mV;
    delete[] memblock;
    calib.close();
    //=========================================================================
    // END OF                CALIBRATION CONSTANTS
    //=========================================================================
    
    
    ifstream data (data_file_name, ios::in|ios::binary);
    ofstream out (out_file_name, ios::out|ios::binary|ios::trunc);
    if (data.is_open() && out.is_open())
    {
        // create our own header
        //out.write( "FACT_COMPRESS                                                                   ", 80);
        //out.write( "END.                                                                            ", 80);
        
        // copy first 0x2d00 bytes to new file
        const int ascii_header_size = 0x2d00;
        char * memblock = new char [ascii_header_size];
        data.read(memblock, ascii_header_size);
        out.write(memblock, ascii_header_size);
        delete[] memblock;
        
        for ( int event_id = 0 ; event_id < 1000; event_id++)
        //while ( data.good() )
        {
            // copy binary header to new file
            const int bin_header_size = 3390;
            char * memblock = new char [bin_header_size];
            data.read(memblock, bin_header_size);
            out.write(memblock, bin_header_size);
            delete[] memblock;
            
            // read 300 shorts out of the file
            for (int chid = 0 ; chid < 1440; chid++)
            {
                const int data_size = 300;
                const int diff_size = data_size;
                char * memblock = new char [data_size*sizeof(short)];
                short * data = (short *)memblock;
                short * diffs = new short [diff_size];
                unsigned char * sizes = new unsigned char [diff_size];
                
                data.read(memblock, data_size*sizeof(short) );
                
                for ( int i = 0; i<diff_size; i++)
                {
                    diffs[i] = data[i]-data[i];
                    //~ if (diffs[i] >= -32 && diffs[i] <= 31)
                        //~ sizes[i] = 6;
                    //~ else if (diffs[i] >= -128 && diffs[i] <= 127)
                        //~ sizes[i] = 8;
                    //~ else if (diffs[i] >= -512 && diffs[i] <= 511)
                        //~ sizes[i] = 10;
                    //~ else
                        //~ sizes[i] = 16;
                    if (diffs[i] >= -128 && diffs[i] <= 127)
                        sizes[i] = 8;
                    else
                        sizes[i] = 16;
                }
                
                // calculate group sizes
                int counter = 0;
                unsigned char last_size = sizes[0];
                vector<int> group_sizes;
                vector<unsigned char> groups;
                
                for (int i=0 ; i < diff_size; i++)
                {
                    
                    if (sizes[i] != last_size)
                    {
                        group_sizes.push_back( counter );
                        groups.push_back( last_size );
                        last_size = sizes[i];
                    }
                    counter++;
                }
                groups.push_back( last_size );
                group_sizes.push_back( counter );
                
                
                // write the first short
                out.write(memblock, 1*sizeof(short) );
                // for all groups write header and group
                short header = 0;
                int diff_index = 0;
                for ( int i = 0 ; i < (int)groups.size() ; i++)
                {
                    // write header
                    if (groups[i] == 8)
                    {
                        header = short(group_sizes[i]);
                    }
                    else
                    {
                        header = 0x8000 | short(group_sizes[i] );
                    }
                    out.write( (char*)&header, sizeof(short));
                    
                    // write group
                    if (groups[i] == 8)
                    {
                        for (int j = 0; j<group_sizes[i]; j++)
                        {
                            out.write( (char*)&(diffs[diff_index++]), 1);
                        }
                    }
                    else
                    {
                        for (int j = 0; j<group_sizes[i]; j++)
                        {
                            out.write( (short*)&(diffs[diff_index++]), 2);
                        }
                    }
                }
                
                //out.write(memblock, data_size*sizeof(short) );
                
                delete[] memblock;
                delete[] diffs;
                delete[] sizes;
            }
        }
        cout << "finished with 1000 events." << endl;
        long after_address = data.tellg();
        data.seekg (0, ios::end);
        long end = data.tellg();
        data.seekg (after_address, ios::beg);
        
        cout << "between last event and end:" << end - after_address << endl;
        // read the last piece...
        const int rest_size = end-after_address;
        char * memblock2 = new char [rest_size];
        data.read(memblock2, rest_size);
        cout << "first char in memblock: " << int(memblock2[0]) << endl;
        char lastchar = memblock2[0];
        for (int i =0 ; i<rest_size; i++)
        {
            if (memblock2[i] != lastchar)
            {
                cout << "new char at: " << i << " with value:" << int(memblock[0]) << endl;
                lastchar = memblock2[i];
            }
        }
        
        out.write(memblock2, rest_size);
        delete[] memblock2;
        
        
        data.close();
        out.close();
  } //end of if file open
} //end of main