File: instrumentMultiplyAdd.cpp

package info (click to toggle)
groops 0%2Bgit20250907%2Bds-1
  • links: PTS, VCS
  • area: non-free
  • in suites: forky, sid
  • size: 11,140 kB
  • sloc: cpp: 135,607; fortran: 1,603; makefile: 20
file content (119 lines) | stat: -rw-r--r-- 3,971 bytes parent folder | download | duplicates (2)
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
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
/***********************************************/
/**
* @file instrumentMultiplyAdd.cpp
*
* @brief Multiply instrument data with a factor and add them together.
*
* @author Torsten Mayer-Guerr
* @date 2012-06-24
*
*/
/***********************************************/

// Latex documentation
#define DOCSTRING docstring
static const char *docstring = R"(
This program multiply \file{instrument data}{instrument} with a factor and add them together.
Afterwards the mean of each arc and data column can be removed with \config{removeArcMean}.
The instrument files must be synchronized (\program{InstrumentSynchronize}).

See also \program{InstrumentArcCalculate}.
)";

/***********************************************/

#include "programs/program.h"
#include "files/fileInstrument.h"

/***** CLASS ***********************************/

/** @brief Multiply instrument data with a factor and add them together.
* @ingroup programsGroup */
class InstrumentMultiplyAdd
{
public:
  class Data
  {
    public:
    FileName fileName;
    Double   factor;
  };

  void run(Config &config, Parallel::CommunicatorPtr comm);
};

GROOPS_REGISTER_PROGRAM(InstrumentMultiplyAdd, PARALLEL, "Multiply instrument data with a factor and add them together", Instrument)

/***********************************************/

template<> Bool readConfig(Config &config, const std::string &name, InstrumentMultiplyAdd::Data &var, Config::Appearance mustSet, const std::string &defaultValue, const std::string &annotation)
{
  if(!readConfigSequence(config, name, mustSet, defaultValue, annotation))
    return FALSE;
  readConfig(config, "inputfileInstrument", var.fileName, Config::MUSTSET,  "", "");
  readConfig(config, "factor",              var.factor,   Config::DEFAULT,  "1.0", "");
  endSequence(config);
  return TRUE;
}

/***********************************************/

void InstrumentMultiplyAdd::run(Config &config, Parallel::CommunicatorPtr comm)
{
  try
  {
    FileName fileNameOut;
    std::vector<Data> data;
    Bool removeMean;

    readConfig(config, "outputfileInstrument", fileNameOut, Config::MUSTSET,  "", "");
    readConfig(config, "instrument",           data,        Config::MUSTSET,  "", "");
    readConfig(config, "removeArcMean",        removeMean,  Config::DEFAULT,  "0", "remove mean value of each arc");
    if(isCreateSchema(config)) return;

    // open files and check consistency
    std::vector<InstrumentFilePtr> instrumentFile(data.size());
    for(UInt i=0; i<instrumentFile.size(); i++)
    {
      logStatus<<"read instrument data <"<<data.at(i).fileName<<">"<<Log::endl;
      instrumentFile.at(i) = InstrumentFile::newFile(data.at(i).fileName);
      InstrumentFile::checkArcCount({*instrumentFile.at(0), *instrumentFile.at(i)});
      if(instrumentFile.at(i)->getType() != instrumentFile.at(0)->getType())
        throw(Exception("instruments types are different: "+instrumentFile.at(i)->getTypeName()+", "+instrumentFile.at(0)->getTypeName()));
    }

    logStatus<<"combine instrument data"<<Log::endl;
    std::vector<Arc> arcList(instrumentFile.at(0)->arcCount());
    Parallel::forEach(arcList, [&](UInt arcNo)
    {
      Arc arc = instrumentFile.at(0)->readArc(arcNo);
      Matrix A = data.at(0).factor * arc.matrix();
      for(UInt i=1; i<instrumentFile.size(); i++)
      {
        Arc arc2 = instrumentFile.at(i)->readArc(arcNo);
        Arc::checkSynchronized({arc, arc2});
        axpy(data.at(i).factor, arc2.matrix(), A);
      }

      if(removeMean)
        for(UInt k=0; k<A.columns(); k++)
          A.column(k) -= mean(A.column(k));

      return Arc(arc.times(), A, arc.getType());
    }, comm);


    if(Parallel::isMaster(comm))
    {
      logStatus<<"write instrument data to file <"<<fileNameOut<<">"<<Log::endl;
      InstrumentFile::write(fileNameOut, arcList);
      Arc::printStatistics(arcList);
    }
  }
  catch(std::exception &e)
  {
    GROOPS_RETHROW(e)
  }
}

/***********************************************/