Code:
#include <stdio.h> // for printf
#include <stdlib.h> // for exit
#include <getopt.h> // for getopt_long
#include <time.h>
#include "parsing.hh"
#include "genfunc.hh"
#include "dynbaseobj.hh"
#include "list.hh"
#include "stack.hh"
#include "tree.hh"
#include "velmod.hh"
#include "map.hh"
#include "ray.hh"
#include "source.hh"
#include "lglevel.hh"
#include "common.hh"
#include "clorfncs.hh"
#include "raytrac_doc.hh"
// Set default values
const double dflt_angli = 0.0; ///< initial shooting angle.
const double dflt_anglf = 360.0; ///< final shooting angle.
const double dflt_dangl = 2.0; ///< increment in the shooting angle.
const Real dflt_tstep = 0.1; ///< time step.
const Real dflt_mdacc = -1.0; ///< mdacc
const Real dflt_mindist = -1.0; ///< mindist
const int dflt_niter = 10; ///< number of iterations.
const char dflt_ofrmt[] = "X T"; ///< Default output format.
const bool dflt_headw = true; ///< Default tracing of headwaves.
const IntgMethod dflt_intgm = euler; ///< Default integration method.
enum RTLevel {
rysht,
hdwvs,
rybrk,
twpnt
};
void write_frmt (
ostream& os,
String& frmt,
int sr,
Vect2 Src,
int ph,
int rc,
Real Rec,
Real A,
Real T,
Vect2 P
);
int main (
int argc,
char* argv[]
) {
for (int i = 0; i < argc;i++) {
printf("***%s***\n", argv[i]);
}
int c;
int digit_optind = 0;
int aopt = 0;
int bopt = 0;
char* val_ifbase = NULL;
char* val_ifmodl = NULL;
char* val_ifrcvs = NULL;
char* val_ifsrcs = NULL;
char* val_lsrcs = NULL;
char* val_lrcvs = NULL;
char* val_ofrays = NULL;
char* val_oftrvt = NULL;
char* val_ofrmt = NULL;
char* val_ifseiphs = NULL;
char* val_lseiphs = NULL;
char* val_intgm = NULL;
char* val_rtlvl = NULL;
char* val_hdwvs = NULL;
char* val_mltrf = NULL;
char* val_tstep = NULL;
char* val_textrp = NULL;
char* val_niter = NULL;
char* val_angli = NULL;
char* val_anglf = NULL;
char* val_dangl = NULL;
char* val_drybrk = NULL;
char* val_dtwpnt = NULL;
char* val_lglvl = NULL;
bool hasargv_ifbase = false;
bool hasargv_ifmodl = false;
bool hasargv_ifrcvs = false;
bool hasargv_ifsrcs = false;
bool hasargv_lsrcs = false;
bool hasargv_lrcvs = false;
bool hasargv_ofrays = false;
bool hasargv_oftrvt = false;
bool hasargv_ofrmt = false;
bool hasargv_ifseiphs = false;
bool hasargv_lseiphs = false;
bool hasargv_intgm = false;
bool hasargv_rtlvl = false;
bool hasargv_hdwvs = false;
bool hasargv_mltrf = false;
bool hasargv_tstep = false;
bool hasargv_textrp = false;
bool hasargv_niter = false;
bool hasargv_angli = false;
bool hasargv_anglf = false;
bool hasargv_dangl = false;
bool hasargv_drybrk = false;
bool hasargv_dtwpnt = false;
bool hasargv_colour = false;
bool hasargv_lglvl = false;
bool hasargv_brief = false;
bool hasargv_detailed = false;
bool hasargv_examples = false;
bool hasargv_help = false;
bool hasargv_invalid = false;
// A string listing valid short options letters.
const char* const short_options = "bdeh";
// An array describing valid long options.
static struct option long_options[] = {
{"ifbase", required_argument, 0, 10},
{"ifmodl", required_argument, 0, 11},
{"ifsrcs", optional_argument, 0, 12},
{"ifrcvs", optional_argument, 0, 13},
{"ifsphs", optional_argument, 0, 14},
{"ofrays", required_argument, 0, 15},
{"oftrvt", required_argument, 0, 16},
{"srcs", optional_argument, 0, 17},
{"rcvs", optional_argument, 0, 18},
{"lseiphs", optional_argument, 0, 19},
{"ofrmt", optional_argument, 0, 20},
{"intgm", optional_argument, 0, 21},
{"rtlvl", optional_argument, 0, 22},
{"hdwvs", no_argument, 0, 23},
{"mltrf", no_argument, 0, 24},
{"textrp", optional_argument, 0, 25},
{"tstep", optional_argument, 0, 26},
{"niter", optional_argument, 0, 27},
{"angli", optional_argument, 0, 28},
{"anglf", optional_argument, 0, 29},
{"dangl", optional_argument, 0, 30},
{"drybrk", optional_argument, 0, 31},
{"dtwpnt", optional_argument, 0, 32},
{"colour", no_argument, 0, 33},
{"lglvl", optional_argument, 0, 34},
{"brief", no_argument, 0, 'b'},
{"detailed", no_argument, 0, 'd'},
{"examples", no_argument, 0, 'e'},
{"help", no_argument, 0, 'h'},
{NULL, 0, NULL, 0} // Required at end of array.
};
int option_index = 0;
while ((c = getopt_long(argc, argv, short_options, long_options, &option_index)) != -1) {
int this_option_optind = optind ? optind : 1;
switch (c) {
case 0:
printf("option h number\n");
printf("option %s", long_options[option_index].name);
if (optarg) printf (" with arg %s", optarg);
printf("\n");
break;
case 10:
printf("--ifbase, option 10 with value '%s'\n\n", optarg );
val_ifbase = optarg;
hasargv_ifbase = true;
break;
case 11:
printf("--ifmodl, option 11, with value '%s'\n\n", optarg );
val_ifmodl = optarg;
hasargv_ifmodl = true;
break;
case 12:
printf("--ifsrcs, option 12, with value '%s'\n\n", optarg );
val_ifsrcs = optarg;
hasargv_ifsrcs = true;
break;
case 13:
printf("--ifrcvs, option 13, with value '%s'\n\n", optarg );
val_ifrcvs = optarg;
hasargv_ifrcvs = true;
break;
case 14:
printf("--ifseiphs, option 19, with value '%s'\n\n", optarg );
val_ifseiphs = optarg;
hasargv_ifseiphs = true;
break;
case 15:
printf("--ofrays, option 16, with value '%s'\n\n", optarg );
val_ofrays = optarg;
hasargv_ofrays = true;
break;
case 16:
printf("--oftrvt, option 17, with value '%s'\n\n", optarg );
val_oftrvt = optarg;
hasargv_oftrvt = true;
break;
case 17:
printf("--lsrcs, option 14, with value '%s'\n\n", optarg );
val_lsrcs = optarg;
hasargv_lsrcs = true;
break;
case 18:
printf("--lrcvs, option 15, with value '%s'\n\n", optarg );
val_lrcvs = optarg;
hasargv_lrcvs = true;
break;
// RAYTRACING ARGUMENTS ----------------------------------------------------
case 19:
printf("--lseiphs, option 20, with value '%s'\n\n", optarg );
val_lseiphs = optarg;
hasargv_lseiphs = true;
break;
case 20:
printf("--trvt-frmt, option 18, with value '%s'\n\n", optarg );
val_ofrmt = optarg;
hasargv_ofrmt = true;
break;
case 21:
printf("--intgm, option 21, with value '%s'\n\n", optarg );
val_intgm = optarg;
hasargv_intgm = true;
break;
case 22:
printf("--rtlvl, option 22, with value '%s'\n\n", optarg );
val_rtlvl = optarg;
hasargv_rtlvl = true;
break;
case 23:
printf("--hdwvs, option 23\n\n" );
hasargv_hdwvs = true;
break;
case 24:
printf("--mltrf, option 24\n\n" );
hasargv_mltrf = true;
break;
case 25:
printf("--tstep, option 25 with value '%s'\n\n", optarg );
val_tstep = optarg;
hasargv_tstep = true;
break;
case 26:
printf("--textrp, option 26\n\n" );
hasargv_textrp = true;
break;
case 27:
printf("--niter, option 27, with value '%s'\n\n", optarg );
val_niter = optarg;
hasargv_niter = true;
break;
case 28:
printf("--angli, option 28, with value '%s'\n\n", optarg );
val_angli = optarg;
hasargv_angli = true;
break;
case 29:
printf("--anglf, option 29, with value '%s'\n\n", optarg );
val_anglf = optarg;
hasargv_anglf = true;
break;
case 30:
printf("--dangl, option 30, with value '%s'\n\n", optarg );
val_dangl = optarg;
hasargv_dangl = true;
break;
case 31:
printf("--drbrck, option 31, with value '%s'\n", optarg );
val_drybrk = optarg;
hasargv_drybrk = true;
break;
case 32:
printf("--dtwpnt, option 32, with value '%s'\n", optarg );
val_dtwpnt = optarg;
hasargv_dtwpnt = true;
break;
// INFORMATION ARGUMENTS ---------------------------------------------------
case 33:
printf("--colour, option 33, with value '%s'\n", optarg );
hasargv_colour = true;
break;
case 34:
printf("--lglvl, option 34, with value '%s'\n", optarg );
val_lglvl = optarg;
hasargv_lglvl = true;
break;
// -b or --brief
case 'b':
hasargv_brief = true;
//print_usage(stdout, 0, clor);
break;
// -d or --detailed
case 'd':
hasargv_detailed = true;
//print_help(stdout, 0, clor);
break;
// -h or --help
case 'h':
hasargv_help = true;
//print_help(stdout, 0, clor);
break;
// -e or --examples
case 'e':
hasargv_examples = true;
//print_examples(stdout, 0, clor);
break;
// invalid option.
case '?':
hasargv_invalid = true;
//print_usage(stderr, 1, clor);
break;
default:
printf("?? getopt returned character code 0%o ??\n", c);
}
}
// Process the remaining non-option arguments. OPTIND points to first non-option argument.
if (optind < argc) {
printf ("Non-option ARGV-elements:\n");
for (int i = optind; i < argc; ++i) {
printf("%s ", argv[i]);
printf("\n");
}
}
//------------------------------------------------------------------------------------------------
bool colourise = false;
if (hasargv_colour) {
colourise = true;
}
if (argc < 2) {
print_usage (stdout, 0, colourise);
}
if (hasargv_brief) {
print_usage (stdout, 0, colourise);
}
if (hasargv_detailed || hasargv_help) {
print_help (stdout, 0, colourise);
}
if (hasargv_examples) {
print_examples (stdout, 0, colourise);
}
//------------------------------------------------------------------------------------------------
Parsing Pc(argc, argv);
//------------------------------------------------------------------------------------------------
// Set logging level
LgLevel lglvl;
String s = "medium";
if (hasargv_lglvl) {
Pc.get_string ("lglvl", s);
}
if ( !get_log_level (s,lglvl) ) {
error ("invalid value for --lglvl");
}
if (lglvl >= medium) {
cerr << "--lglvl = " << s << " " << lglvl << endl;
}
//------------------------------------------------------------------------------------------------
// Program description
String msg = "raytrac 2-D raytracing program";
if (lglvl >= medium) {
cout << "\nPROGRAM\n" << msg << "\n\n" << endl;
}
Real mdacc = dflt_mdacc;
Real mindist = dflt_mindist;
Map M;
Velmod vm;
Source** Src;
List<Phase> Ph;
List<Vect2> Sr;
List<Real> Rec;
List<String> Ls;
String s1;
int nPhases;
int nsrcs;
int nrcvs;
int i;
//------------------------------------------------------------------------------------------------
// ifmodl
// name of file containing the model
String ifmodl;
ifstream ifs_modl;
Parsing Pvm;
if ( !Pc.get_string("--ifmodl", ifmodl) ) {
error("--ifmodl required");
}
if (lglvl >= low) {
cerr << "reading 2-D sound speed model " << ifmodl << endl;
}
ifs_modl.open(ifmodl);
if (ifs_modl.bad()) {
error("--ifl-modl failed to open");
}
Pvm.parse_file(ifs_modl);
vm.read(Pvm);
//------------------------------------------------------------------------------------------------
// ifseiphs
// File containing seismic phases to model
String ifseiphs;
String lseiphs;
Phase phase;
if (Pc.get_string("--ifseiphs", ifseiphs)) {
ifstream pph;
Parsing Pph;
pph.open(ifseiphs);
if (pph.bad()) {
error ("seismis phases file failed to open");
}
Pph.parse_file(pph);
if ( !Pph.get_list("PHASES", Ls) ) {
error ("No phases declared in seismic phases file");
}
nPhases = Ls.size();
for (i = 0; i < nPhases; i++) {
phase.convert(Ls[i]);
Ph += phase;
phase.reset();
}
//------------------------------------------------------------------------------------------------
// lst-sphs
// List of seismic phases to model
} else if (Pc.get_string("--lseiphs", lseiphs)) {
cerr << "lseiphs = " << lseiphs;
nPhases = lseiphs.nfields(':');
cerr << ", nPhases = " << nPhases << endl;
for (i = 0; i < nPhases; i++) {
s1 = lseiphs.get_token(':');
cerr << "Phase " << i+1 << ", s1 = " << s1 << endl;
phase.convert(s1); // Convert phase string to Phase class
Ph += phase; // Append phase to the list of phases
cerr << "phase is " << phase;
phase.reset();
}
} else {
error ("--ifseiphs or --lseiphs argument required");
}
if (lglvl == high) {
for (i = 0; i < nPhases; i++) {
cerr << endl << "Ph[" << i+1 << "] = " << Ph[i] << endl;
}
cerr << endl;
}
// Definitions: List<Phase> Ph, Map M
// In M.Add(Ph[i]), we create Tree<ElMap>* T = M
for (i = 0; i < nPhases; i++) {
M.add(Ph[i]);
} // Calls Map::Add
if (lglvl >= medium) {
M.print(cerr);
} // Calls Map::Print
cerr << endl;
// READ SOURCE LOCATIONS -------------------------------------------------------
Vect2 source;
//------------------------------------------------------------------------------------------------
// ifsrcs File containing list of sources
String ifsrcs;
String lsrcs;
if (Pc.get_string("--ifsrcs", ifsrcs)) {
ifstream psr;
Parsing Psr;
cerr << "reading sources from file" << s;
psr.open(ifsrcs);
if (psr.bad()) {
error ("srfile failed to open");
}
Psr.parse_file(psr);
if ( !Psr.get_list("SOURCES", Ls) ) {
error ("no sources declared in sources file");
}
nsrcs = Ls.size();
Src = new Source *[nsrcs];
for (i = 0; i < nsrcs; i++) {
Ls[i] >> source;
Src[i] = new Source(source, Ph);
Sr += source;
}
//------------------------------------------------------------------------------------------------
// lsrcs List of sources
} else if (Pc.get_string("--lsrcs", lsrcs)) {
nsrcs = lsrcs.nfields(':');
Src = new Source *[nsrcs];
for (i = 0; i < nsrcs; i++) {
s1 = lsrcs.get_token(':');
s1 >> source;
Src[i] = new Source(source, Ph);
Sr += source;
}
} else {
error ("--ifsrcs or --lsrcs argument required");
}
if (lglvl >= medium) {
cerr << ", nsrcs = " << nsrcs << endl;
}
if (lglvl == high) {
for (i = 0; i < nsrcs; i++) {
cerr << " Source " << (i+1) << " Location " << Sr[i] << endl;
}
cerr << endl;
}
// READ RECEIVERS LOCATIONS
Real receiver;
//------------------------------------------------------------------------------------------------
// ifrcvs File containing list of receivers
String ifrcvs;
String lrcvs;
if (Pc.get_string("--ifrcvs", ifrcvs)) {
ifstream prc;
Parsing Prc;
prc.open(ifrcvs);
cerr << "Reading receiver locations from file" << s;
if (prc.bad()) {
error ("receivers file failed to open");
}
Prc.parse_file(prc);
if ( !Prc.get_list("RECEIVERS", Ls) ) {
error ("No receivers declared in receivers file");
}
nrcvs = Ls.size();
cerr << ", nrcvs = " << nrcvs << endl;
for (i = 0; i < nrcvs; i++) {
Ls[i] >> receiver;
Rec += receiver;
}
//------------------------------------------------------------------------------------------------
// lrcvs List of receivers
} else if (Pc.get_string("--lrcvs", lrcvs)) {
if (lrcvs[0] == '@') {
Real xi;
Real xf;
Real dx;
lrcvs >> xi >> xf >> nrcvs;
dx = (xf - xi) / (nrcvs - 1);
for (i = 0; i < nrcvs; i++) {
Rec += xi + (i * dx);
}
} else {
nrcvs = lrcvs.nfields(':');
for (i=0; i<nrcvs; i++) {
s1 = lrcvs.get_token(':');
s1 >> receiver;
Rec += receiver;
}
}
} else {
nrcvs = 0;
}
if (lglvl == high) {
for (i = 0; i < nrcvs; i++) {
cerr << " Receiver " << (i+1) << ". Location " << Rec[i] << endl;
}
cerr << endl;
}
//------------------------------------------------------------------------------------------------
// rtlvl Level of raytracing
s = String("twpnt");
RTLevel rtlvl = twpnt;
if (Pc.get_string("--rtlvl",s)) {
if (s == String("rysht")) { // ray shooting
rtlvl = rysht;
} else if (s == String("hdwvs")) { // trace headwaves
rtlvl = hdwvs;
} else if (s == String("rybrk")) { // receiver bracketing
rtlvl = rybrk;
} else if (s == String("twpnt")) { // two point raytracing
rtlvl = twpnt;
} else { // invalid trace level
error("invalid value for argument --rtlvl");
}
}
if ((nrcvs == 0) && (rtlvl == twpnt)) {
rtlvl = rybrk;
s = String("rybrk");
}
if (lglvl >= medium) {
cerr << "rtlvl = " << s << endl;
}
//------------------------------------------------------------------------------------------------
// ofrays File containing of file containing ray paths
String ofrays;
ofstream ofs_rays;
if (Pc.get_string("--ofrays", ofrays)) {
ofs_rays.open(ofrays);
if (ofs_rays.bad()) {
error ("File containing raypaths failed to open for writing");
}
if (lglvl >= medium) {
cerr << "--ofrays = " << ofrays << endl;
}
}
//------------------------------------------------------------------------------------------------
// oftrvt File containing ray travel-times
String oftrvt;
ofstream ofs_trvt;
if ( !Pc.get_string("--oftrvt", oftrvt) ) {
error ("--oftrvt is a required argument");
}
ofs_trvt.open(oftrvt);
if (ofs_trvt.bad()) {
error ("File containing travel times failed to open for writing");
}
if (lglvl >= medium) {
cerr << "--oftrvt = " << oftrvt << endl;
}
//------------------------------------------------------------------------------------------------
// ofrmt
String ofrmt = "xrcv trvt";
if (hasargv_ofrmt) Pc.get_string("--ofrmt", ofrmt);
// if ( !Pc.set_string("--ofrmt", ofrmt) ) {
// }
if (lglvl >= medium) {
cerr << "--ofrmt = " << ofrmt << endl;
}
//------------------------------------------------------------------------------------------------
// intgm Ray integration method
s = "euler";
IntgMethod intgm;
if (hasargv_intgm) {
Pc.get_string("--intgm", s);
}
if ( !get_integration_method(s, intgm) ) {
error ("--intgm is a required argument");
}
if (lglvl >= medium) {
cerr << "--intgm = " << s << endl;
}
//------------------------------------------------------------------------------------------------
// tstep
// integration time
Real tstep = 0.1;
if (hasargv_tstep) {
Pc.get_real("--tstep", tstep);
}
if (lglvl >= medium) {
cerr << "--tstep = " << tstep << endl;
}
//------------------------------------------------------------------------------------------------
// angli
// initial shooting angle
double angli = 0.0;
if (hasargv_angli) {
Pc.get_double("--angli", angli);
}
if (lglvl >= medium) {
cerr << "--angli = " << angli << endl;
}
//------------------------------------------------------------------------------------------------
// anglf
// final shooting angle
double anglf = 360.0;
if (hasargv_anglf) {
Pc.get_double("--anglf", anglf);
}
// if (Pc.get_double("--anglf", anglf)) {
// }
if (lglvl >= medium) {
cerr << "--anglf = " << anglf << endl;
}
//------------------------------------------------------------------------------------------------
// dangl increment in the shooting angle
double dangl = 1.0;
if (hasargv_dangl) {
Pc.get_double("--dangl", dangl);
}
// if (Pc.get_double("--dangl", dangl)) {
// }
if (lglvl >= medium) {
cerr << "--dangl = " << dangl << endl;
}
//------------------------------------------------------------------------------------------------
// niter maximum number of allowable iterations
int niter = dflt_niter;
if (Pc.get_int("--niter", niter)) {
if (lglvl >= medium) {
cerr << "--niter = " << niter << endl;
}
}
//------------------------------------------------------------------------------------------------
// rybrk mdacc
if (Pc.get_string("--rybrk", s)) {
if (s == String("always")) {
mdacc = -1;
} else {
s >> mdacc;
}
if (lglvl >= medium) {
cerr << "--mdacc = " << mdacc << endl;
}
}
//------------------------------------------------------------------------------------------------
// twpnt
if (Pc.get_string("--twpnt", s)) {
if (s == String("auto")) {
mindist = -1;
} else {
s >> mindist;
}
if (lglvl >= medium) {
cerr << "--mindist = " << mindist << endl;
}
}
//------------------------------------------------------------------------------------------------
// hdwvs Trace head-waves
String hdwvs_sw;
if (Pc.get_string("--hdwvs", hdwvs_sw)) {
s.to_upper();
if (hdwvs_sw == String("on")) {
key_hdwvs = true;
} else if (hdwvs_sw == String("off")) {
key_hdwvs = false;
} else {
error("invalid value for --headwv");
}
if (lglvl >= medium) {
cerr << "--headwv = " << s << endl;
}
}
// multr ---------------------------------------------------------------------
// traces multiple reflections
String mltrf_sw;
if ( Pc.get_string("--mltrf", mltrf_sw) ) {
if (mltrf_sw == String("on")) {
key_mltrf = true;
} else if (mltrf_sw == String("off")) {
key_mltrf = false;
} else {
error("invalid value for --multrf");
}
if (lglvl >= medium) {
cerr << "--multrf = " << s << endl;
}
}
// textrp --------------------------------------------------------------------
// extrapolates traveltimes
String textrp_sw;
if ( Pc.get_string("--textrp", textrp_sw) ) {
if (textrp_sw == String("on")) {
key_textrp = true;
} else if (textrp_sw == String("off")) {
key_textrp = false;
} else {
error("invalid value for argument --textrp");
}
if (lglvl >= medium) {
cerr << "--textrp = " << textrp_sw << endl;
}
}
// Ray Tracing -----------------------------------------------------------------
time_t Begin = time(NULL);
int sr;
int rc;
int ph;
List<Real> A;
List<Real> T;
List<Vect2> P;
List<Vect2> X;
if (lglvl >= low) {
cerr << endl << endl
<< endl << "raytrac.cc:" << endl
<< "Ray shooting"
<< endl << "Store results in the .ry and .tx files" << endl;
}
for (i = 0; i < nsrcs; i++) {
Src[i]->set_time_step(tstep);
Src[i]->set_vm(&vm);
Src[i]->set_integration_method(intgm);
}
// Output
if (hasargv_ofrays) {
ofs_rays << '>' << endl;
}
ofs_trvt << '>' << endl;
// Start for loop < sr >
for (sr = 0; sr < nsrcs; sr++) {
if (lglvl >= medium) {
cerr << "Shot " << sr+1 << ". " << Src[sr]->get_src();
}
if (lglvl >= medium) {
cerr << " rshoot";
}
if ( !Src[sr]->ray_shooting(angli, anglf, dangl) ) {
if (sr < 9) {
cerr << "Shot " << sr+1 << ". " << endl;
} else {
cerr << "Shot " << sr+1 << ". " << endl;
}
error ("Src[sr]->ray_shooting failed");
}
if (hasargv_hdwvs && (rtlvl >= hdwvs)) {
if (lglvl >= medium) cerr << ", hdwvs";
if (!Src[sr]->trace_headwaves()) {
error ("Src[sr]->trace_headwaves failed");
}
}
for (ph = 0; ph < nPhases; ph++) {
if (rtlvl >= rybrk) {
if (lglvl >= medium) {
cerr << ", rybrk";
}
if (!Src[sr]->receiver_bracketing(ph)) {
error ("receiver_bracketing failed");
}
}
if (rtlvl == twpnt) {
if (lglvl >= medium) {
cerr << ", twopnt";
}
if (!Src[sr]->two_point(ph, Rec, mindist, niter)) {
error("Src[sr]->two_point failed");
}
}
}
cerr << endl;
// ----- ----- ----- ----- ----- ----- ----- ----- ----- ----- ----- ---
if (nrcvs > 0) {
for (ph = 0; ph < nPhases; ph++) {
for (rc = 0; rc < Rec.size(); rc++) {
Src[sr]->get(Rec[rc], A, T, P, ph, key_textrp, mdacc);
if (hasargv_mltrf || (A.size() == 1)) {
for (i = 0; i < A.size(); i++) {
write_frmt (ofs_trvt, ofrmt, sr, Sr[sr], ph, rc, Rec[rc], A[i], T[i], P[i]);
}
} else if (A.size() > 0) {
int min = 0;
for (i = 1; i < T.size(); i++) {
if (T[i] < T[min]) { min = i; }
}
write_frmt (ofs_trvt, ofrmt, sr, Sr[sr], ph, rc, Rec[rc], A[min], T[min], P[min]);
}
}
ofs_trvt << '>' << endl;
}
if (hasargv_ofrays) {
for (ph = 0; ph < nPhases; ph++) {
for (rc = 0; rc < Rec.size(); rc++) {
Src[sr]->write_rays (ofs_rays, Rec[rc], ph, key_mltrf, mdacc, key_textrp);
}
}
}
} else {
for (ph = 0; ph < nPhases; ph++) {
Src[sr]->get(A, T, P, X ,ph);
for (i = 0; i < A.size(); i++) {
write_frmt (ofs_trvt, ofrmt, sr, Sr[sr], ph, rc, X[i].X, A[i], T[i], P[i]);
}
ofs_trvt << '>' << endl;
}
if (hasargv_ofrays) {
Src[sr]->write_all_rays(ofs_rays);
}
} // End loop < if nrcvs >
Src[sr]->reset();
} // End for loop < sr >
ofs_rays.close();
if (lglvl >= medium) {
cerr << endl << "computation time = "
<< (time(NULL)-Begin) << "s " << "(" << (time(NULL)-Begin)/60 << "min)" << endl;
}
return (0);
}
void write_frmt (
ostream& os,
String& frmt,
int sr,
Vect2 Src,
int ph,
int rc,
Real Rec,
Real A,
Real T,
Vect2 P
) {
int nf = frmt.nfields();
String F = frmt;
int i;
String s;
for (i = 0; i < nf; i++) {
s = F.get_token();
if (s == String("xrcv")) {
os << Rec << ' ';
} else if (s == String("xsrc")) {
os << Src.X << ' ';
} else if (s == String("zsrc")) {
os << Src.Z << ' ';
} else if (s == String("trvt")) {
os << T << ' ';
} else if (s == String("angl")) {
os << A << ' ';
} else if (s == String("pray")) {
os << P.X << ' ';
} else if (s == String("pxray")) {
os << P.Z << ' ';
} else if (s == String("pzray")) {
os << P.X << ' ';
} else if (s == String("srdist")) {
os << ( Rec - Src.X ) << ' ';
} else if (s == String("tau")) {
os << (T - P.X * (Rec - Src.X)) << ' ';
} else if (s == String("isrc")) {
os << sr << ' ';
} else if (s == String("ircv")) {
os << rc << ' ';
} else if (s == String("iphs")) {
os << ph << ' ';
} else {
error("invalid value for argument --ofrmt");
}
}
os << endl;
}