/home/pascal/depot/filedata/src/parasolfiledata.h

/***************************************************************************
 *   Copyright (C) 2005 by Nicolas PASCAL   *
 *   nicolas.pascal@icare.univ-lille1.fr   *
 *                                                                         *
 *   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.             *
 ***************************************************************************/
#ifndef PARASOLFILEDATA_H
#define PARASOLFILEDATA_H

#include "filedata.h"
#include "parasolleader.h"
#include "parasoldata.h"
#include "satellitefiledata.h"

#include <cmath>
#include <cfloat>
#include "hdfi.h"

class PARASOLFileData : public SatelliteFileData
{
public:
    static const int nb_parasol_direction;
    // altitude of the satellite in m. It is a virtual mean altitude. The real movement is an ellipse, not a perfect circle
    static const double satellite_altitude;
    // size in m of a NADIR pixel border (6km)
    static const double pix_sz_nadir;
    // maximum size in m of an observed pixel border, ie a pixel observed with a zenithal angle of 60 degrees
    static const double pix_sz_max;
    // maximum angle between viewibng direction and satellite nadir : 60 degrees
    static const double max_sat_nadir_angle;
    // angle in radians at satellite position on an observed nadir pixel at full resolution (= 1/16 degree earth angle)
    static const double pix_obs_nadir_angle;
    // maximum number of acquisition sequences
    static const int nb_acq_seq_max;
    // number of images during an acquisition sequence
    static const int nb_acq_img_max;

protected:
    PARASOLLeader* leader;
    PARASOLData* data;
private:
    template <typename NumType>
    static bool greater_comparator (NumType i,NumType j) { return (i > j); };
    enum ProcessingLine{
        UNDEFINED,
        BASIC,
        LAND_SURFACE,
        OCEAN_COLOR,
        RADIATION_BUDGET
    };

    static const double lat2time_coef;

    int instrument;
    int level;
    ProcessingLine processing_line;
    char product;
    int orbit_cycle_nb;
    int orbit_nb;
    char reprocessing_nb;
    int grid_factor;
    int sz_grid [2];
    double sz_cell;
    string software_version_number;
    unsigned short *line_data;
    unsigned short *col_data;
    float lat_min,lat_max,lon_min,lon_max;
    // inner data buffers used to speed up viewing directions computations
    short * v_alt_view;
    double * v_saa;
    double * v_vza;
    double * v_raa;
    unsigned char * v_ni;
    unsigned char * v_sn;
    double * v_x_sat;
    double * v_y_sat;
    double * v_z_sat;

    void init_viewing_direction();

    void init();
    void parse_filename(const string& short_filename);
    bool check_filename(const string& short_filename) const;
    void check_grid_coord (const int irow, const int icol) const;
    void set_grid_factor();
    bool init_read_data_range(const int entry_idx, int * &start, int * &edges, int &rank);
    void* read_count_data(void* data, const int var_idx = -1, int * start = NULL, int * edges = NULL, int rank = -1);
    void get_scaling(const int var_idx, double &slope, double &offset, unsigned short& nb_bytes);

    const EntryFormat *get_entry_format(const int &var_idx);
protected:
    vector<int> get_data_dimension(const int var_idx);
public:
    PARASOLFileData(const string& name, const string& mode="r");
    ~PARASOLFileData();
    const string get_radix();
    const int get_instrument() const {
        return instrument;
    }
    const int get_level() const {
        return level;
    }
    const ProcessingLine get_processing_line() const {
        return processing_line;
    }
    const string get_processing_line_string() const;
    const char get_processing_line_char() const;
    const char get_product() const {
        return product;
    }
    const int get_orbit_cycle_nb() const {
        return orbit_cycle_nb;
    }
    const int get_orbit_nb() const {
        return orbit_nb;
    }
    const char get_reprocessing_nb() const {
        return reprocessing_nb;
    }
    void geolocation_to_grid(const float * lat_lon, int * grid_idx);
    void geolocation_to_grid(const float lat, const float lon, int & irow, int & icol);
    void grid_to_geolocation (const int igrid[2], float lat_minmax[2], float lon_minmax[2]);
    void grid_to_geolocation (const int igrid[2], float pos[2]);
    void grid_to_geolocation (const int igrid[2], float & lat, float & lon);
    void grid_to_geolocation (const int irow, const int icol, float & lat, float & lon);
    const int get_nrow() const;
    const int get_ncol(const int irow = -1) const;
    const int get_ncol(const double lat) const;
    bool is_earth_grid (const int irow, const int icol) const;
    const double get_grid_lat (const int irow) const;
    void get_record_coord(const int irec, float &lat, float &lon);
    void idata2igrid (const int irec, int igrid[2]);
    vector<int> get_data_dimension(const string &var_name);
    void* read_count_data(void * data, const char* var_name = "alt", int * start = NULL, int * edges = NULL, int rank = -1);
    double* read_scaled_data(double * data, const char * var_name = "alt", int * start = NULL, int * edges = NULL, int rank = -1);
    void * read_data(void * data, const char * var_name, int * start = NULL, int * stride = NULL, int * edges = NULL, int rank = -1);
    const int get_nb_data();
    virtual void load_geolocation_data();
    virtual void free_geolocation_data();
    void load_v_pixel() ;
    virtual const bool is_geolocation_data_loaded() const;
    const bool contain_data(const float &lat,const float &lon, const double &time);
    const bool contain_data(const float &lat,const float &lon, const double &time, const double &tolerance);
    const bool contain_location(const float &lat,const float &lon);
    const bool contain_location(const float &lat,const float &lon, const double &tolerance);
    const bool contain_pixel(const int* pix_idx);
    void get_scaling(const string &var_name, double &slope, double &offset, unsigned short& nb_bytes);
    const EntryFormat *get_entry_format(const string &var_name);
    void open_data_file(){
        assert(data!=NULL);
        data->open_file();
    };
    void close_data_file(){
        assert(data!=NULL);
        data->close_file();
    };
    const bool is_data_file_opened(){
        return (data!=NULL && data->is_file_loaded());
    };
    unsigned short* get_line_data() const {
        assert(line_data!=NULL); // TODO do load geolocation_data before if NULL
        return line_data;
    }
    unsigned short* get_col_data() const {
        assert(col_data!=NULL); // TODO do load geolocation_data before if NULL
        return col_data;
    }
    const bool get_index(const float& lat, const float& lon, int &irec);
    const bool get_index(const float &lat, const float& lon, int &irec, const float colocation_tolerance);
    const bool get_index(const int ipix[2], int & irec);
    const bool get_index(const vector <int> & ipix, int & irec);
    void igrid2idata (const int igrid[2], int & irec);
    void igrid2idata (const vector<int> & igrid, int & irec);
     virtual void get_vindex(vector < vector < int > > &v_index, const float &lat, const float& lon,
                                            const float colocation_tolerance) {
        UnimplementedMethod e(__FILE__,__LINE__,__PRETTY_FUNCTION__);
        throw e;
     };

    const float get_nearest_point_distance(const float & lat, const float & lon, const float coloc_tolerance);

    const bool is_file_loaded();
    void close_file();
    static const int get_max_direction();
    void get_pixel(const int &irec, int ipix[2]);
    void print_scaling_factors(){
        leader->print_scaling_factors();
    };
    void print_spatio_temp_char(){
        leader->print_spatio_temp_char();
    };
    void print_techno_param(){
        leader->print_techno_param();
    };
    inline const double get_record_time(const int &irec) {
      assert(leader->is_record_loaded(SPATIO_TEMPORAL_CHARACTERISTICS));
      int pixel[]={-1,-1};
      // record to pixel mapping
      get_pixel(irec,pixel);
      // pixel to geolocation mapping
      float latminmax[]={-1.,-1.};
      float lonminmax[]={-1.,-1.};
      grid_to_geolocation(pixel,latminmax,lonminmax);
      float lat=(latminmax[1]-latminmax[0])/2.;
      // Ascending node time
      string s_asc_node_time(leader->spatio_temp_char->desc_node_date);
      s_asc_node_time=s_asc_node_time.substr(0,s_asc_node_time.size()-2);
      Date d; d.set_date_str(s_asc_node_time,"%Y%m%d%H%M%S");
      return d.get_TAI93_time()+lat*lat2time_coef;
    };

    static inline const double get_filter_view_zenith_angle(const double &view_zenith_angle_filter8,
                                                             const double &relative_azimuth_angle_filter8,
                                                             const double &delta_thetav_cosphi,
                                                             const double &delta_thetav_sinphi,
                                                             const int    xj) {
        double a1=view_zenith_angle_filter8*cos(relative_azimuth_angle_filter8) +
                      (double)(xj)*delta_thetav_cosphi;
        double a2=view_zenith_angle_filter8*sin(relative_azimuth_angle_filter8) +
                      (double)(xj)*delta_thetav_sinphi;
        return sqrt(pow (a1,2) + pow(a2,2));
    }
    static inline const double get_filter_relative_azimuth_angle(const double &view_zenith_angle_filter8,
                                                                  const double &relative_azimuth_angle_filter8,
                                                                  const double &delta_thetav_cosphi,
                                                                  const double &delta_thetav_sinphi,
                                                                  const int    xj) {
        double a1=view_zenith_angle_filter8*sin(relative_azimuth_angle_filter8) +
                      (double)(xj)*delta_thetav_sinphi;
        double a2=view_zenith_angle_filter8*cos(relative_azimuth_angle_filter8) +
                      (double)(xj)*delta_thetav_cosphi;
        double offset=0.;
        if (a2<0)
            offset+=180.;
        return offset + atan (a1/a2) ;
    }
    virtual void get_pixel_coord (const vector < int > & ipix, float &lat, float &lon, double &time) {
        UnimplementedMethod e(__FILE__,__LINE__,__PRETTY_FUNCTION__);
        throw e;
    };
    const string get_software_version_number() const {
        return software_version_number;
    }
//     virtual void load_viewing_directions ();
    virtual void get_viewing_directions (const vector <int> & ipix, vector <Observation> & v_obs);
    virtual bool is_viewing_directions_data_loaded();
    virtual void load_viewing_directions_data ();
    virtual void free_viewing_directions_data ();
    virtual void load_pix2data_map ();
//     void get_neighbours_coord (const int irow, const int icol, const int d_irow, const int d_icol, vector<int> & v_irow, vector<int> & v_icol);
    void get_neighbours_coord (const int irow, const int icol, int dirow_min, int dirow_max, int dicol_min, int dicol_max, vector<int> & v_irow, vector<int> & v_icol);
    void get_v_icol_neighbours (const int irow, const int icol, int dicol_min, int dicol_max, vector <int> & v_icol);

    unsigned char * get_ndir_data() const;
private:
    void _get_viewing_directions_l1(const vector< int > & ipix, const int irec, vector< Observation > & v_obs);
    void _get_viewing_directions_l2(const vector< int > & ipix, const int irec, vector< Observation > & v_obs);

//     virtual void _load_viewing_directions_l1 ();
//     virtual void _load_viewing_directions_l2 ();
};

#endif