/***************************************************************************
 *   Copyright (C) 2006 by Nicolas PASCAL   *
 *   pascal@icare-pc12   *
 *                                                                         *
 *   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.             *
 ***************************************************************************/

#include "parasolfiledata.h"

int main(int argc, char *argv[]) {

    // ***** print sat constants *****
    cout<<"--- Constants ---"<<endl;
    cout << "satellite_altitude " << PARASOLFileData::satellite_altitude << " m" << endl;
    cout << "pix_sz_nadir " << PARASOLFileData::pix_sz_nadir << " m" << endl;
    cout << "pix_obs_nadir_angle " << PARASOLFileData::pix_obs_nadir_angle << " rad " << degrees(PARASOLFileData::pix_obs_nadir_angle) << " deg" << endl;
    cout << "max_sat_nadir_angle " << PARASOLFileData::max_sat_nadir_angle << " rad " << degrees(PARASOLFileData::max_sat_nadir_angle) << " deg" << endl;
    cout << "pix_sz_max " << PARASOLFileData::pix_sz_max << " m" << endl;
    cout<<"-------------------------------------------------------"<<endl;

    string filename("/DATA/LIENS/PARASOL/L1_B/2009/2009_01_01/P3L1TBG1094004KL");
    PARASOLFileData pfd(filename);

//     pfd.print_spatio_temp_char();
//     pfd.print_techno_param();

//     int nrec=pfd.get_nb_data();

    // ***** test some attributes building *****
//     cout<<"-------------------------------------------------------"<<endl;
//     cout<<"Instrument "<<pfd.get_instrument()<<endl;
//     cout<<"Level "<<pfd.get_level()<<endl;
//     cout<<"Processing line "<<pfd.get_processing_line_string()<<endl;
//     cout<<"Product "<<pfd.get_product()<<endl;
//     cout<<"Orbit cycle number "<<pfd.get_orbit_cycle_nb()<<endl;
//     cout<<"Orbit number "<<pfd.get_orbit_nb()<<endl;
//     cout<<"Reprocessing number "<<pfd.get_reprocessing_nb()<<endl;
//     cout<<"Radix "<<pfd.get_radix()<<endl;
//     cout<<"Software Version Number "<<pfd.get_software_version_number()<<endl;
    cout << "satellite_altitude " << PARASOLFileData::satellite_altitude << " m" << endl;
    cout << "pix_sz_nadir " << PARASOLFileData::pix_sz_nadir << " m" << endl;
    cout << "pix_obs_nadir_angle " << PARASOLFileData::pix_obs_nadir_angle << " rad " << degrees(PARASOLFileData::pix_obs_nadir_angle) << " deg" << endl;
    cout << "max_sat_nadir_angle " << PARASOLFileData::max_sat_nadir_angle << " rad " << degrees(PARASOLFileData::max_sat_nadir_angle) << " deg" << endl;
    cout << "pix_sz_max " << PARASOLFileData::pix_sz_max << " m" << endl;
    cout<<"-------------------------------------------------------"<<endl;

//     /* test neighbour pixels methods */
//     cout << "--- test neighbour pixels indices methods ---" << endl;
//     vector <int> v_irow_near, v_icol_near;
//     int irow, icol, dicol_min = -5, dicol_max = 5;
// 
//     cout << "-> columns" << endl;
//     irow = 1, icol = 3240, dicol_min = 0, dicol_max = 0;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     irow = 1, icol = 3240, dicol_min = -1, dicol_max = 1;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     irow = 1, icol = 3240, dicol_min = -5, dicol_max = 5;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     irow = 1, icol = 3240, dicol_min = -2, dicol_max = 0;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     irow = 1, icol = 3241, dicol_min = 0, dicol_max = 2;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     irow = 2, icol = 3240, dicol_min = -2, dicol_max = 2;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     irow = 2, icol = 3240, dicol_min = -6, dicol_max = 1;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     irow = 2, icol = 3240, dicol_min = -1, dicol_max = 6;
//     pfd.get_v_icol_neighbours (irow, icol, dicol_min, dicol_max, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " << MyTools::vec2str(v_icol_near) << endl;
// 
//     cout << "-> rows + columns" << endl;
//     int dirow_min = -2, dirow_max = 2;
// 
//     irow = 1, icol = 3240, dirow_min = -2, dirow_max = 2, dicol_min = -2, dicol_max = 2;
//     pfd.get_neighbours_coord (irow, icol, dirow_min, dirow_max, dicol_min, dicol_max, v_irow_near, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " <<  " dirow [" << dirow_min << ", " << dirow_max << "]" << endl;
//     cout << "v_irow " << MyTools::vec2str(v_irow_near) << endl;
//     cout << "v_icol " << MyTools::vec2str(v_icol_near) << endl;
//     cout << endl;
// 
//     irow = 1, icol = 3240, dirow_min = 0, dirow_max = 3, dicol_min = -10, dicol_max = 10;
//     pfd.get_neighbours_coord (irow, icol, dirow_min, dirow_max, dicol_min, dicol_max, v_irow_near, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " <<  " dirow [" << dirow_min << ", " << dirow_max << "]" << endl;
//     cout << "v_irow " << MyTools::vec2str(v_irow_near) << endl;
//     cout << "v_icol " << MyTools::vec2str(v_icol_near) << endl;
//     cout << endl;
// 
//     irow = 3087, icol = 2787, dirow_min = -9, dirow_max = 9, dicol_min = -13, dicol_max = 15;
//     pfd.get_neighbours_coord (irow, icol, dirow_min, dirow_max, dicol_min, dicol_max, v_irow_near, v_icol_near);
//     cout << "irow " << irow << " icol " << icol << " dicol [" << dicol_min << ", " << dicol_max << "] " <<  " dirow [" << dirow_min << ", " << dirow_max << "]" << endl;
//     cout << "v_irow " << MyTools::vec2str(v_irow_near) << endl;
//     cout << "v_icol " << MyTools::vec2str(v_icol_near) << endl;
//     cout << endl;

    // ***** test geolocation data reading *****
    /*
    Init the geolocation data : it loads the lat,lon indexes in memory. To mandatory, but fasten methods that map lat,lon indexes to a record index like get_index(lat,lon,rec_idx).
    */
//     pfd.load_geolocation_data();
//     pfd.print_pix2data_map();
//     unsigned short* lin=const_cast<unsigned short*>(pfd.get_line_data());
//     unsigned short* col=const_cast<unsigned short*>(pfd.get_col_data());
//     if (true) {
//         cout<<"line\tcolumn"<<endl;
//         for (int i = 0 ; i<nrec ; ++i)
//             cout<<lin[i]<<"\t"<<col[i]<<endl;
//     }
//     cout<<"-------------------------------------------------------"<<endl; // OK

    // ***** test grid/record/(lat, lon) mapping methods *****
    cout << "--- check grid/record/(lat, lon) mapping methods ---" << endl;

    // igrid = [3218, 3195] pos = [-88.750, -115.352] irec = 0 (2 in the data because starts at 2)
    float pos[] = {0.,0.};
    int   ipix[2];
    int  irec;

    // pos -> igrid
//     pos[0] = -88.750, pos[1] = -115.352;
    pos[0] = -56.446396, pos[1] = 127.66864;

    pfd.geolocation_to_grid(pos, ipix);
    cout << "point (" << pos[0] << ", " << pos[1] << ") -> ipix [" << ipix[0] << ", " << ipix[1] << "]" << endl;
    // ipix -> pos
    pfd.grid_to_geolocation(ipix, pos);
    cout << "ipix [" << ipix[0] << ", " << ipix[1] << "] -> point (" << pos[0] << ", " << pos[1] << ")" << endl;
    // pos -> irec
    pfd.get_index(pos[0], pos[1], irec);
    cout << "point (" << pos[0] << "," << pos[1] << ") -> record " << irec << endl;
    // irec -> pos
    pfd.get_record_coord(irec, pos[0], pos[1]);
    cout << "irec " << irec << " -> point ( " << pos[0] << "," << pos[1] << ")" << endl;
    // irec -> ipix
    pfd.idata2igrid(irec, ipix);
    cout << "irec " << irec << " -> ipix (" << ipix[0] << ", " << ipix[1] << ")" << endl;
    // ipix -> irec
    pfd.igrid2idata(ipix, irec);
    cout << "ipix [" << ipix[0] << ", " << ipix[1] << "] -> irec " << irec << "" << endl;

//     /* --- read some data --- */
//     cout << "--- check data read ---" << endl;
//     cout << "-> load whole dataset" << endl;
//     int sz_dir = 16;
//     unsigned short* v_irow = const_cast<unsigned short*>(pfd.get_line_data());
//     unsigned short* v_icol = const_cast<unsigned short*>(pfd.get_col_data());
//     // number of available directions. 1 per record
//     char * v_ni = (char *)(pfd.read_count_data(NULL, "ni"));
//     // solar azimuth angles. 1 per record
//     unsigned char * v_saa = (unsigned char *)(pfd.read_count_data(NULL, "phis"));
//     double * v_saa_scaled = (double *)(pfd.read_data(NULL, "phis"));
//     // dqx. N_directions per record
//     unsigned short * v_dqx = (unsigned short *)(pfd.read_count_data(NULL, "dqx"));
//     cout << "irow\ticol\tirec\tni\tphis_raw\tphis\tdqx_1\tdqx_2..." << endl;
// //     for (int irec = 0 ; irec < nrec ; ++irec) {
//     for (int irec = 97552 ; irec < 97552 + 1 ; ++irec) {
//         // non directionnal params
//         cout << v_irow[irec] << "\t" << v_icol[irec] << "\t" << irec << "\t" << (int)(v_ni[irec]) << "\t" << (int)(v_saa[irec]) << "\t" << v_saa_scaled[irec];
//         for (int idir = 0 ; idir < sz_dir ; ++idir) {
//             if (! MyTools::get_endianness())
//                 MyTools::byte_swap(v_dqx [irec * sz_dir + idir]);
//             cout << "\t" << v_dqx [irec * sz_dir + idir];
//         }
//         cout << endl;
//     }
//     delete [] v_saa, v_saa = NULL;
//     delete [] v_saa_scaled, v_saa_scaled = NULL;
//     delete [] v_dqx, v_dqx = NULL;
// 
//     // viewing zenith angles. N_directions per record
//     unsigned short * v_vza = (unsigned short *)(pfd.read_count_data(NULL, "thetav"));
//     double * v_vza_scaled  = (double *)(pfd.read_data(NULL, "thetav"));
//     cout << "irow\ticol\tirec\tni\t\tthetav_1_raw\tthetav_1\tthetav_2_raw..." << endl;
// //     for (int irec = 0 ; irec < nrec ; ++irec) {
//     for (int irec = 97552 ; irec < 97552 + 1 ; ++irec) {
//         // non directionnal params
//         cout << v_irow[irec] << "\t" << v_icol[irec] << "\t" << irec << "\t" << (int)(v_ni[irec]) << "\t";
//         for (int idir = 0 ; idir < sz_dir ; ++idir) {
//            cout << "\t" << v_vza [irec * sz_dir + idir] << "\t" << v_vza_scaled [irec * sz_dir + idir];
//         }
//         cout << endl;
//     }
//     delete [] v_vza, v_vza = NULL;
//     delete [] v_vza_scaled, v_vza_scaled = NULL;
// 
//     /* --- read a unique data value --- */
//     cout << "-> read a unique value" << endl;
// 
//     int start_2D[] = {97552, 0}; // [irecord, idirection]
//     int edges_2D[] = {1, 1};
//     unsigned short ushort_val;
// 
//     cout << "* thetav *" << endl;
//     for (int idir = 0 ; idir < 16 ; ++idir) {
//         start_2D[1] = idir;
//         pfd.read_count_data(&ushort_val, "thetav", start_2D, edges_2D, 2);
//         cout << ushort_val << " ";
//     }
//     cout << endl;
// 
//     cout << "* dqx *" << endl;
//     for (int idir = 0 ; idir < 16 ; ++idir) {
//         start_2D[1] = idir;
//         pfd.read_count_data(&ushort_val, "dqx", start_2D, edges_2D, 2);
//         if (! MyTools::get_endianness())
//             MyTools::byte_swap(ushort_val);
//         cout << ushort_val << " ";
//     }
//     cout << endl;

    /* --- check pixels viewing directions --- */
    cout << "--- check pixels viewing directions ---" << endl;
    vector <Observation> v_obs;
    pfd.load_viewing_directions_data();

    // get all directions
//     cout << "-> all directions" << endl;
//     vector <int> v_ipix2(2);
//     v_ipix2 [0] = 711, v_ipix2 [1] = 4358;
//     pfd.get_viewing_directions(v_ipix2, v_obs);
//     for (vector <Observation>::iterator it_obs = v_obs.begin() ; it_obs != v_obs.end() ; ++it_obs) {
//         cout << *it_obs << endl;
//         cout << "Geodetic " << it_obs->viewing.get_p1().to_geodetic() << " -> " << it_obs->viewing.get_p2().to_geodetic() << endl;
//     }

    // get one direction
    cout << "-> one by one" << endl;
    vector <int> v_ipix3(3);
    v_ipix3 [0] = 711, v_ipix3 [1] = 4358;
    for (int idir = 10 ; idir < 11 ; ++idir) {
        v_ipix3 [2] = idir;
        cout << "idir = " << idir << endl;
        pfd.get_viewing_directions(v_ipix3, v_obs);
        for (vector <Observation>::iterator it_obs = v_obs.begin() ; it_obs != v_obs.end() ; ++it_obs) {
            cout << *it_obs << endl;
            cout << "Geodetic " << it_obs->viewing.get_p1().to_geodetic() << " -> " << it_obs->viewing.get_p2().to_geodetic() << endl;
        }
    }

    pfd.free_viewing_directions_data();

    return 0;
}