GyotoScreen.h

Go to the documentation of this file.

/*
    Copyright 2011-2020, 2023-2024 Thibaut Paumard, Frederic Vincent

    This file is part of Gyoto.

    Gyoto 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 3 of the License, or
    (at your option) any later version.

    Gyoto 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 Gyoto.  If not, see <http://www.gnu.org/licenses/>.
 */

#ifndef __GyotoScreen_H_
#define __GyotoScreen_H_ 

#include <iostream>
#include <fstream>
#include <string>
#if defined HAVE_BOOST_ARRAY_HPP
# include <boost/array.hpp>
# define GYOTO_ARRAY boost::array
# if defined HAVE_MPI
#  include <boost/version.hpp>
#  if BOOST_VERSION >= 106400 
#   include <boost/serialization/boost_array.hpp>
#   include <boost/serialization/array_wrapper.hpp>
#  endif
# endif
#else
template <typename T, size_t sz> class GYOTO_ARRAY {
 private:
  T buf[sz];
 public:
  T& operator[](size_t c) { return buf[c] ; }
};
#endif

namespace Gyoto {
  class Screen;
}

#include <GyotoDefs.h>
#include <GyotoUtils.h>
#include <GyotoSmartPointer.h>
#include <GyotoObject.h>
#include <GyotoMetric.h>
#include <GyotoSpectrometer.h>

class Gyoto::Screen
: public Gyoto::SmartPointee,
  public Gyoto::Object
{
  friend class Gyoto::SmartPointer<Gyoto::Screen>;

 private:
  double tobs_; 
  double fov_;  
  double azimuthal_fov_; 
  //  double tmin_;
  size_t npix_; 

  double * mask_;

  std::string mask_filename_;

  double distance_; 
  double dmax_; 

  enum anglekind_e { equatorial_angles=0, rectilinear=1, spherical_angles=2};
  typedef int anglekind_t;

  anglekind_t anglekind_; 
  
  double euler_[3]; 
  double ex_[3]; 
  double ey_[3]; 
  double ez_[3]; 

  double fourvel_[4]; 
  double screen1_[4]; 
  double screen2_[4]; 
  double screen3_[4]; 

  double dangle1_; 
  double dangle2_; 
  SmartPointer<Metric::Generic> gg_; 

  SmartPointer<Spectrometer::Generic> spectro_;

  double freq_obs_;

  obskind_t observerkind_;

 public:
  GYOTO_OBJECT;
  GYOTO_OBJECT_THREAD_SAFETY;

  // Constructors - Destructor
  // -------------------------
  Screen() ; 
  Screen(const Screen& ) ;                
  Screen * clone() const; 

  virtual ~Screen() ;                        
  
  // Mutators / assignment
  // ---------------------

  void setProjection(const double paln,
                     const double inclination,
                     const double argument);
  void setProjection(const double distance,
                     const double paln,
                     const double inclination,
                     const double argument);


  void distance(double dist);


  void dMax(double dist);


  void distance(double dist, const std::string &unit);


  void inclination(double);


  void inclination(double, const std::string &unit);

  void PALN(double);
  void PALN(double, const std::string &unit);
  void argument(double);
  void argument(double, const std::string &unit);
  void spectrometer(SmartPointer<Spectrometer::Generic> spectro);
  SmartPointer<Spectrometer::Generic> spectrometer() const ;

  void freqObs(double fo);


  void freqObs(double fo, const std::string &unit);

  double freqObs() const ;

  double freqObs(const std::string &unit) const;


  void setObserverPos(const double pos[4]);
  void observerKind(const std::string &kind);
  std::string observerKind() const;
  void setFourVel(const double coord[4]);
  void setScreen1(const double coord[4]);
  void setScreen2(const double coord[4]);
  void setScreen3(const double coord[4]);

  // Accessors
  // ---------


  int coordKind() const;


  double distance() const;


  double distance(const std::string&) const;     


  double dMax() const;


  double inclination() const;


  double inclination(const std::string&) const;

  double PALN() const;   
  double PALN(const std::string&) const;         
  double argument() const;       
  double argument(const std::string&) const;     

  SmartPointer<Metric::Generic> metric() const; 
  void metric(SmartPointer<Metric::Generic> gg); 

  double time() const;

  double time(const std::string &) const;

  void time(double, const std::string &);

  void time(double);

  double fieldOfView() const;

  double fieldOfView(std::string const &unit) const;

  void fieldOfView(double);

  void fieldOfView(double, const std::string &unit);

  double azimuthalFieldOfView() const;

  void azimuthalFieldOfView(double ff);

  void dangle1(double);
  void dangle1(double, const std::string &unit);
  double dangle1() const;
  double dangle1(std::string const &unit)const;
  void dangle2(double);
  void dangle2(double, const std::string &unit);
  double dangle2() const;
  double dangle2(std::string const &unit)const;

  void anglekind(int);
  void anglekind(std::string const&);
  std::string anglekind() const;

  size_t resolution() const;
  void resolution(size_t);


  void mask(double const * const mm, size_t resolution=0);

  double const * mask() const ;
  void maskFile(std::string const &fname);
  std::string maskFile() const;
# ifdef GYOTO_USE_CFITSIO

  void fitsReadMask(std::string const &fname);

  void fitsWriteMask(std::string const &fname);
# endif


  bool operator()(size_t, size_t);



  void getObserverPos(double dest[4]) const;


  void getFourVel(double dest[4]) const;

  void fourVel(std::vector<double> const &);
  std::vector<double> fourVel() const;
  void screenVector1(std::vector<double> const &);
  std::vector<double> screenVector1() const;
  void screenVector2(std::vector<double> const &);
  std::vector<double> screenVector2() const;
  void screenVector3(std::vector<double> const &);
  std::vector<double> screenVector3() const;


  void getScreen1(double dest[4]) const;


  void getScreen2(double dest[4]) const;


  void getScreen3(double dest[4]) const;


  void getRayTriad(double x, double y,
                   double dest[8],
                   bool compute_polar_basis=false,
                   double Ephi[4]=NULL, double Etheta[4]=NULL) const;


  void getRayTriad(const size_t i, const size_t j,
                   double dest[8],
                   bool compute_polar_basis=false,
                   double Ephi[4]=NULL, double Etheta[4]=NULL) const;
  
  void coordToSky(const double pos[4], double dest[3], bool geometrical=false) const;

  void skyToCoord(const double sky[3], double dest[4], bool geometrical=false) const;

  void coordToXYZ(const double pos[4], double dest[3]) const;

  void computeBaseVectors() ;

  std::ostream& print(std::ostream&) const ; 
  std::ostream& printBaseVectors(std::ostream&) const ; 

  // UDUNITS
# ifdef HAVE_UDUNITS

  void mapPixUnit();


  void unmapPixUnit();
# endif


#ifdef GYOTO_USE_XERCES
 public:
  void fillProperty(Gyoto::FactoryMessenger *fmp, Property const &p) const;

    static   SmartPointer<Screen> Subcontractor(FactoryMessenger* fmp);
#endif

    enum CoordType_e {angle, pixel};


    class Coord1dSet {
    public:
      const CoordType_e kind;
    public:
      Coord1dSet(CoordType_e k);
      virtual ~Coord1dSet();
      virtual void begin() =0;
      virtual bool valid() =0;
      virtual size_t size()=0;
      virtual size_t operator*() const ;
      virtual double angle() const ;
      virtual Coord1dSet& operator++()=0;

      virtual size_t index() const=0;
    };


    class Coord2dSet {
    public:
      const CoordType_e kind;
      Coord2dSet(CoordType_e k);
      virtual ~Coord2dSet();
      virtual Coord2dSet& operator++()    =0;
      virtual GYOTO_ARRAY<size_t, 2> operator*  () const;
      virtual GYOTO_ARRAY<double, 2> angles() const ;
      virtual void begin() =0;
      virtual bool valid() =0;
      virtual size_t size()=0;
    };

    class Grid: public Coord2dSet {
    protected:
    protected:
      char * prefix_;
      Coord1dSet &iset_;
      Coord1dSet &jset_;
    public:
      Grid(Coord1dSet &iset, Coord1dSet &jset, const char * const p=NULL);
      virtual ~Grid();
      virtual Coord2dSet& operator++();
      virtual GYOTO_ARRAY<size_t, 2> operator*  () const;
      virtual void begin();
      virtual bool valid();
      virtual size_t size();
    };


    class Bucket : public Coord2dSet {
    protected:
      Coord1dSet &alpha_;
      Coord1dSet &delta_;
    public:
      Bucket(Coord1dSet &iset, Coord1dSet &jset);
      virtual Coord2dSet& operator++();
      virtual GYOTO_ARRAY<double, 2> angles() const;
      virtual GYOTO_ARRAY<size_t, 2> operator*() const;
      virtual void begin();
      virtual bool valid();
      virtual size_t size();
    };

    class Empty: public Coord2dSet {
    public:
      Empty();
      virtual Coord2dSet& operator++();
      virtual void begin();
      virtual bool valid();
      virtual size_t size();
    };

    class Range : public Coord1dSet {
    protected:
      const size_t mi_, ma_, d_, sz_;
      size_t cur_;
    public:
      Range(size_t mi, size_t ma, size_t d);
      void begin();
      bool valid();
      size_t size();
      Coord1dSet& operator++();
      size_t operator*() const ;
      virtual size_t index() const ;
    };

    class Indices : public Coord1dSet {
    protected:
      size_t * indices_;
      size_t const sz_;
      size_t i_;
    public:
      Indices (size_t const*const buf, size_t sz);
      ~Indices();
      void begin();
      bool valid();
      size_t size();
      Coord1dSet& operator++();
      size_t operator*() const ;
      virtual size_t index() const ;
    };

    class Angles : public Coord1dSet {
    protected:
      double * buf_;
      size_t const sz_;
      size_t i_;
    public:
      Angles (double const*const buf, size_t sz);
      ~Angles();
      void begin();
      bool valid();
      size_t size();
      Coord1dSet& operator++();
      double angle() const ;
      virtual size_t index() const ;
    };

    class RepeatAngle : public Coord1dSet {
    protected:
      double const val_;
      size_t const sz_;
      size_t i_;
    public:
      RepeatAngle (double val, size_t sz);
      void begin();
      bool valid();
      size_t size();
      Coord1dSet& operator++();
      double angle() const ;
      virtual size_t index() const ;
    };
};

#endif