stfinvidentity.cc 4.95 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30
/*! \file stfinvidentity.cc
 * \brief just find a scaling factor (implementation)
 * 
 * ----------------------------------------------------------------------------
 * 
 * \author Thomas Forbriger
 * \date 07/05/2011
 * 
 * just find a scaling factor (implementation)
 * 
 * Copyright (c) 2011 by Thomas Forbriger (BFO Schiltach) 
 *
 * ----
 * 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., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
 * ----
 *
 * 
 * REVISIONS and CHANGES 
31 32 33 34 35
 *  - 07/05/2011   V1.0   Thomas Forbriger (thof)
 *  - 21/02/2014   V1.1   implemented scaling to average weighted energy
 *  - 14/10/2015   V1.2   new end-user usage functions
 *  - 19/10/2015   V1.3   weight factor is defined in terms amplitude; apply
 *                        square of weight to signal energy
36 37 38 39
 * 
 * ============================================================================
 */
#define STFINV_STFINVIDENTITY_CC_VERSION \
40
  "STFINV_STFINVIDENTITY_CC   V1.3"
41 42

#include <stfinv/stfinvidentity.h>
43 44 45 46 47 48
#include <stfinv/stfinvidentity_summary_usage.h>
#include <stfinv/stfinvidentity_description_usage.h>
#include <aff/functions/sqrsum.h>
#include <aff/seriesoperators.h>
#include <cmath>
#include <stfinv/tools.h>
49 50 51 52 53 54

namespace stfinv {

  const char* const STFEngineIdentity::ID="ident";

  const char* const STFEngineIdentity::description
55
    ="scale with amplitude factor";
56 57 58 59 60 61 62 63 64 65

  /*----------------------------------------------------------------------*/

  void STFEngineIdentity::help(std::ostream& os) const
  {
    STFEngineIdentity::classhelp(os);
  } // void STFEngineIdentity::help(std::ostream& os) const

  /*----------------------------------------------------------------------*/

66 67 68 69 70 71 72
  void STFEngineIdentity::usage(std::ostream& os) const
  {
    STFEngineIdentity::classusage(os);
  } // void STFEngineIdentity::usage(std::ostream& os) const

  /*----------------------------------------------------------------------*/

73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89
  const char* STFEngineIdentity::name() const
  {
    return("STFEngineIdentity");
  } //  const char const* STFEngineIdentity::name() const

  /*----------------------------------------------------------------------*/

  void STFEngineIdentity::initialize() 
  {
    // scale energy
    Mscaleenergy=(this->parameter("scaleenergy","false")=="true");
  } // void STFEngineIdentity::initialize()

  /*----------------------------------------------------------------------*/

  void STFEngineIdentity::exec() 
  {
90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109
    // effective amplitude factor
    double fac=1.;

    // scale to reproduce average energy if requested
    if (Mscaleenergy)
    {
      double recording_sqrsum=0.;;
      double synthetic_sqrsum=0.;;
      for (unsigned int i=0; i<this->nreceivers(); ++i)
      {
        synthetic_sqrsum 
          += aff::func::sqrsum(this->synthetic(i)) 
             * this->weight(i) * this->weight(i);
        recording_sqrsum 
          += aff::func::sqrsum(this->recording(i)) 
             * this->weight(i) * this->weight(i);
        fac = std::sqrt(recording_sqrsum/synthetic_sqrsum);
      }
    } // if (Mscaleenergy)

110 111
    Tseries stf=this->stf();
    stf=0.;
112
    stf(0)=fac/this->dt();
113 114 115 116 117
    for (unsigned int i=0; i<this->nreceivers(); ++i)
    {
      Tseries::Tcoc synthetic=this->synthetic(i);
      Tseries convolvedsynthetic=this->convolvedsynthetic(i);
      convolvedsynthetic.copyin(synthetic);
118 119 120 121 122 123 124 125
      if (Mscaleenergy) { convolvedsynthetic *= fac; }
    }
    for (unsigned int i=0; i<this->npairs(); ++i)
    {
      Tseries::Tcoc series=this->series(i);
      Tseries convolvedseries=this->convolvedseries(i);
      convolvedseries.copyin(series);
      if (Mscaleenergy) { convolvedseries *= fac; }
126 127 128 129 130 131 132
    }
  } // void STFEngineIdentity::exec()

  /*----------------------------------------------------------------------*/

  void STFEngineIdentity::classhelp(std::ostream& os)
  {
133 134 135
    os << stfinvidentity_summary_usage;
    os << std::endl;
    stfinv::tools::report_engine_ID<STFEngineIdentity>(os);
136 137
  } // void STFEngineIdentity::classhelp(std::ostream& os)

138 139 140 141 142 143 144 145 146 147 148
  /*----------------------------------------------------------------------*/

  void STFEngineIdentity::classusage(std::ostream& os)
  {
    os << stfinvidentity_description_usage;
    os << std::endl;
    Tbase::classusage(os);
    os << std::endl;
    stfinv::tools::report_engine_ID<STFEngineIdentity>(os);
  } // void STFEngineIdentity::classusage(std::ostream& os)

149 150 151
} // namespace stfinv

/* ----- END OF stfinvidentity.cc ----- */