mc_antenna.cc

Go to the documentation of this file.

/* Copyright (C) 2005-2008 Cory Davis <cory@met.ed.ac.uk>
                            
   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. */



/*===========================================================================
  === File description 
  ===========================================================================*/

/*===========================================================================
  === External declarations
  ===========================================================================*/

#include "mc_antenna.h"

#ifdef HAVE_SSTREAM
#include <sstream>
#else
#include "sstream.h"
#endif




Numeric ran_gaussian (
                      Rng& rng, 
                      const Numeric sigma)
{
  Numeric x, y, r2;

  do
    {
      /* choose x,y in uniform square (-1,-1) to (+1,+1) */

      x = -1 + 2 * rng.draw();
      y = -1 + 2 * rng.draw();

      /* see if it is in the unit circle */
      r2 = x * x + y * y;
    }
  while (r2 > 1.0 || r2 == 0);

  /* Box-Muller transform */
  return sigma * y * sqrt (-2.0 * log (r2) / r2);
}

void MCAntenna::set_pencil_beam (void){
  atype=ATYPE_PENCIL_BEAM;
}


void MCAntenna::set_gaussian (const Numeric& za_sigma,
                              const Numeric& aa_sigma)
{
  atype=ATYPE_GAUSSIAN;
  sigma_za=za_sigma;
  sigma_aa=aa_sigma;
}


void MCAntenna::set_gaussian_fwhm (const Numeric& za_fwhm,
                                   const Numeric& aa_fwhm)
{
  atype=ATYPE_GAUSSIAN;
  sigma_za=za_fwhm/2.3548;
  sigma_aa=aa_fwhm/2.3548;
}


void MCAntenna::set_lookup (ConstVectorView& za_grid_,
                            ConstVectorView& aa_grid_,
                            ConstMatrixView& G_lookup_)
{
  atype=ATYPE_LOOKUP;
  za_grid=za_grid_;
  aa_grid=aa_grid_;
  G_lookup=G_lookup_;
}


AType MCAntenna::get_type(void) const
{
  return atype;
}





void MCAntenna::draw_los(VectorView& sampled_rte_los,
                         Rng& rng,
                         ConstVectorView bore_sight_los) const
{

  switch ( atype )
    {
    case ATYPE_PENCIL_BEAM:
      sampled_rte_los=bore_sight_los;
      break;
    case ATYPE_GAUSSIAN:
      sampled_rte_los[0]=bore_sight_los[0]+ran_gaussian(rng,sigma_za);
      sampled_rte_los[1]=bore_sight_los[1]+ran_gaussian(rng,sigma_aa);
      if ( sampled_rte_los[1]>180 )
        {
          sampled_rte_los[1]-=360;
        }
      break;
      //    case ATYPE_LOOKUP:
      //ostringstream os;
      //os << "Antenna type ATYPE_LOOKUP not yet implemented.";
      //throw runtime_error( os.str() );
      //break;
    default:
      ostringstream os;
      os << "invalid Antenna type.";
      throw runtime_error( os.str() );
  }
  
}

ostream& operator<< (ostream& os, const MCAntenna&)
{
  os << "MCAntenna: Output operator not implemented";
  return os;
}

---