linefunctions.h

Go to the documentation of this file.

/* Copyright (C) 2017
 * Richard Larsson <ric.larsson@gmail.com>
 * 
 * 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, 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 linefunctions_h
#define linefunctions_h

#include "complex.h"
#include "energylevelmap.h"
#include "jacobian.h"
#include "absorptionlines.h"

namespace Linefunctions {

constexpr Index ExpectedDataSize() { return 2; }

void set_lineshape(Eigen::Ref<Eigen::VectorXcd> F,
                   const Eigen::Ref<const Eigen::VectorXd> f_grid,
                   const Absorption::SingleLine& line,
                   const Numeric& temperature,
                   const Numeric& zeeman_df,
                   const Numeric& magnetic_magnitude,
                   const Numeric& doppler_constant,
                   const LineShape::Output& lso,
                   const LineShape::Type lineshape_type,
                   const Absorption::MirroringType mirroring_type,
                   const Absorption::NormalizationType norm_type);

void set_lorentz(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    Eigen::Ref<Eigen::Matrix<Complex, Eigen::Dynamic, ExpectedDataSize()>> data,
    const Eigen::Ref<const Eigen::VectorXd> f_grid,
    const Numeric& zeeman_df,
    const Numeric& magnetic_magnitude,
    const Numeric& F0_noshift,
    const LineShape::Output& lso,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex(),
    const LineShape::Output& dT = {0, 0, 0, 0, 0, 0, 0, 0, 0},
    const LineShape::Output& dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0});

void set_htp(Eigen::Ref<Eigen::VectorXcd> F,
             Eigen::Ref<Eigen::MatrixXcd> dF,
             const Eigen::Ref<const Eigen::VectorXd> f_grid,
             const Numeric& zeeman_df,
             const Numeric& magnetic_magnitude,
             const Numeric& F0_noshift,
             const Numeric& GD_div_F0,
             const LineShape::Output& lso,
             const AbsorptionLines& band=AbsorptionLines(),
             const Index& line_ind=0,
             const ArrayOfRetrievalQuantity& derivatives_data =
                 ArrayOfRetrievalQuantity(),
             const ArrayOfIndex& derivatives_data_position = ArrayOfIndex(),
             const Numeric& dGD_div_F0_dT = 0.0,
             const LineShape::Output& dT = {0, 0, 0, 0, 0, 0, 0, 0, 0},
             const LineShape::Output& dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0});

void set_voigt(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    Eigen::Ref<Eigen::Matrix<Complex, Eigen::Dynamic, ExpectedDataSize()>> data,
    const Eigen::Ref<const Eigen::VectorXd> f_grid,
    const Numeric& zeeman_df,
    const Numeric& magnetic_magnitude,
    const Numeric& F0_noshift,
    const Numeric& GD_div_F0,
    const LineShape::Output& lso,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex(),
    const Numeric& dGD_div_F0_dT = 0.0,
    const LineShape::Output& dT = {0, 0, 0, 0, 0, 0, 0, 0, 0},
    const LineShape::Output& dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0});

void set_doppler(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    Eigen::Ref<Eigen::Matrix<Complex, Eigen::Dynamic, ExpectedDataSize()>> data,
    const Eigen::Ref<const Eigen::VectorXd> f_grid,
    const Numeric& zeeman_df,
    const Numeric& magnetic_magnitude,
    const Numeric& F0_noshift,
    const Numeric& GD_div_F0,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex(),
    const Numeric& dGD_div_F0_dT = 0.0);

void apply_linemixing_scaling_and_mirroring(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    const Eigen::Ref<Eigen::VectorXcd> Fm,
    const Eigen::Ref<Eigen::MatrixXcd> dFm,
    const LineShape::Output& lso,
    const bool with_mirroring,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex(),
    const LineShape::Output& dT = {0, 0, 0, 0, 0, 0, 0, 0, 0},
    const LineShape::Output& dVMR = {0, 0, 0, 0, 0, 0, 0, 0, 0});

void apply_rosenkranz_quadratic_scaling(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    const Eigen::Ref<const Eigen::VectorXd> f_grid,
    const Numeric& F0,
    const Numeric& T,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex());

void apply_VVH_scaling(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    Eigen::Ref<Eigen::Matrix<Complex, Eigen::Dynamic, ExpectedDataSize()>> data,
    const Eigen::Ref<const Eigen::VectorXd> f_grid,
    const Numeric& F0,
    const Numeric& T,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex());

void apply_VVW_scaling(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    const Eigen::Ref<const Eigen::VectorXd> f_grid,
    const Numeric& F0,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex());

Numeric lte_linestrength(Numeric S0,
                         Numeric E0,
                         Numeric F0,
                         Numeric QT0,
                         Numeric T0,
                         Numeric QT,
                         Numeric T);

void apply_linestrength_scaling_by_lte(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    Eigen::Ref<Eigen::VectorXcd> N,
    Eigen::Ref<Eigen::MatrixXcd> dN,
    const Absorption::SingleLine& line,
    const Numeric& T,
    const Numeric& T0,
    const Numeric& isotopic_ratio,
    const Numeric& QT,
    const Numeric& QT0,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex(),
    const Numeric& dQT_dT = 0.0);

void apply_linestrength_scaling_by_vibrational_nlte(
  Eigen::Ref<Eigen::VectorXcd> F,
  Eigen::Ref<Eigen::MatrixXcd> dF,
  Eigen::Ref<Eigen::VectorXcd> N,
  Eigen::Ref<Eigen::MatrixXcd> dN,
  const Absorption::SingleLine& line,
  const Numeric& T,
  const Numeric& T0,
  const Numeric& Tu,
  const Numeric& Tl,
  const Numeric& Evu,
  const Numeric& Evl,
  const Numeric& isotopic_ratio,
  const Numeric& QT,
  const Numeric& QT0,
  const AbsorptionLines& band=AbsorptionLines(),
  const Index& line_ind=0,
  const ArrayOfRetrievalQuantity& derivatives_data =
      ArrayOfRetrievalQuantity(),
  const ArrayOfIndex& derivatives_data_position = ArrayOfIndex(),
  const Numeric& dQT_dT = 0.0);

void apply_lineshapemodel_jacobian_scaling(
  Eigen::Ref<Eigen::MatrixXcd> dF,
  const AbsorptionLines& band,
  const Index& line_ind,
  const ArrayOfRetrievalQuantity& derivatives_data,
  const ArrayOfIndex& derivatives_data_position,
  const Numeric& T,
  const Numeric& P,
  const Vector& vmrs);

Numeric DopplerConstant(Numeric T, Numeric mass);

Numeric dDopplerConstant_dT(const Numeric& T, const Numeric& dc);

void find_cutoff_ranges(Index& start_cutoff,
                        Index& nelem_cutoff,
                        const Eigen::Ref<const Eigen::VectorXd> f_grid,
                        const Numeric& fmin,
                        const Numeric& fmax);

void apply_linestrength_from_nlte_level_distributions(
    Eigen::Ref<Eigen::VectorXcd> F,
    Eigen::Ref<Eigen::MatrixXcd> dF,
    Eigen::Ref<Eigen::VectorXcd> N,
    Eigen::Ref<Eigen::MatrixXcd> dN,
    const Numeric& r1,
    const Numeric& r2,
    const Numeric& g1,
    const Numeric& g2,
    const Numeric& A21,
    const Numeric& F0,
    const Numeric& T,
    const AbsorptionLines& band=AbsorptionLines(),
    const Index& line_ind=0,
    const ArrayOfRetrievalQuantity& derivatives_data =
        ArrayOfRetrievalQuantity(),
    const ArrayOfIndex& derivatives_data_position = ArrayOfIndex());

class InternalData {
public:
  Eigen::VectorXcd F;
  Eigen::VectorXcd N;
  Eigen::MatrixXcd dF;
  Eigen::MatrixXcd dN;
  
  Eigen::Matrix<Complex, 1, 1> Fc;
  Eigen::Matrix<Complex, 1, 1> Nc;
  Eigen::Matrix<Complex, 1, Eigen::Dynamic> dFc;
  Eigen::Matrix<Complex, 1, Eigen::Dynamic> dNc;
  
  Eigen::Matrix<Complex, Eigen::Dynamic, Linefunctions::ExpectedDataSize()> data;
  Eigen::Matrix<Complex, 1, Linefunctions::ExpectedDataSize()> datac;
  
  InternalData(Index nf, Index nj) {
    F.setZero(nf);
    N.setZero(nf);
    dF.setZero(nf, nj);
    dN.setZero(nf, nj);
    data.setZero(nf, Linefunctions::ExpectedDataSize());
    
    Fc.setZero(1);
    Nc.setZero(1);
    dFc.setZero(1, nj);
    dNc.setZero(1, nj);
    datac.setZero(1, Linefunctions::ExpectedDataSize());
  }
  
  void SetZero() {
    F.setZero(F.rows());
    N.setZero(N.rows());
    dF.setZero(dF.rows(), dF.cols());
    dN.setZero(dN.rows(), dN.cols());
  }
};  // InternalData

void set_cross_section_of_band(
  InternalData& scratch,
  InternalData& sum,
  const ConstVectorView f_grid,
  const AbsorptionLines& band,
  const ArrayOfRetrievalQuantity& derivatives_data,
  const ArrayOfIndex& derivatives_data_active,
  const Vector& vmrs,
  const EnergyLevelMap& nlte,  // This must be turned into a map of some kind...
  const Numeric& P,
  const Numeric& T,
  const Numeric& isot_ratio,
  const Numeric& H,
  const Numeric& DC,
  const Numeric& dDCdT,
  const Numeric& QT,
  const Numeric& dQTdT,
  const Numeric& QT0,
  const bool no_negatives=false,
  const bool zeeman=false,
  const Zeeman::Polarization zeeman_polarization=Zeeman::Polarization::Pi);
};  // namespace Linefunctions

#endif  //linefunctions_h