|
@@ -0,0 +1,95 @@
|
|
1
|
+/* This file is part of Lemma, a geophysical modelling and inversion API.
|
|
2
|
+ * More information is available at http://lemmasoftware.org
|
|
3
|
+ */
|
|
4
|
+
|
|
5
|
+/* This Source Code Form is subject to the terms of the Mozilla Public
|
|
6
|
+ * License, v. 2.0. If a copy of the MPL was not distributed with this
|
|
7
|
+ * file, You can obtain one at http://mozilla.org/MPL/2.0/.
|
|
8
|
+ */
|
|
9
|
+
|
|
10
|
+/**
|
|
11
|
+ * @file
|
|
12
|
+ * @date 11/11/2016 02:44:37 PM
|
|
13
|
+ * @version $Id$
|
|
14
|
+ * @author Trevor Irons (ti)
|
|
15
|
+ * @email tirons@egi.utah.edu
|
|
16
|
+ * @copyright Copyright (c) 2016, University of Utah
|
|
17
|
+ * @copyright Copyright (c) 2016, Lemma Software, LLC
|
|
18
|
+ */
|
|
19
|
+
|
|
20
|
+#include <Merlin>
|
|
21
|
+using namespace Lemma;
|
|
22
|
+
|
|
23
|
+
|
|
24
|
+int main(int argc, char** argv) {
|
|
25
|
+
|
|
26
|
+ if (argc<4) {
|
|
27
|
+ std::cout << "./KernelV0-2 earth.yaml tx.yaml rx.yaml \n";
|
|
28
|
+ return(EXIT_SUCCESS);
|
|
29
|
+ }
|
|
30
|
+ std::cout << "Using earth model: " << argv[1] << std::endl;
|
|
31
|
+ auto earth = LayeredEarthEM::DeSerialize( YAML::LoadFile(argv[1]) );
|
|
32
|
+
|
|
33
|
+ std::cout << "Using transmitter: " << argv[2] << std::endl;
|
|
34
|
+ auto Tx = PolygonalWireAntenna::DeSerialize( YAML::LoadFile(argv[2]) );
|
|
35
|
+
|
|
36
|
+ std::cout << "Using receivers: " << argv[3] << std::endl;
|
|
37
|
+ auto Rx1 = PolygonalWireAntenna::DeSerialize( YAML::LoadFile(argv[3]) );
|
|
38
|
+
|
|
39
|
+ auto Kern = KernelV0::NewSP();
|
|
40
|
+ Kern->PushCoil( "Coil 1", Tx );
|
|
41
|
+ Kern->PushCoil( "Coil 2", Rx1 );
|
|
42
|
+ Kern->SetLayeredEarthEM( earth );
|
|
43
|
+
|
|
44
|
+ Kern->SetIntegrationSize( (Vector3r() << 200,200,200).finished() );
|
|
45
|
+ Kern->SetIntegrationOrigin( (Vector3r() << 0,0,0).finished() );
|
|
46
|
+ Kern->SetTolerance( 1e-10 ); // 1e-12
|
|
47
|
+
|
|
48
|
+// Kern->AlignWithAkvoDataset( YAML::LoadFile(argv[2]) );
|
|
49
|
+
|
|
50
|
+ Kern->SetPulseDuration(0.020);
|
|
51
|
+ VectorXr I(36);
|
|
52
|
+ // off from VC by 1.075926340216996
|
|
53
|
+ // Pulses from Wyoming Red Buttes exp 0
|
|
54
|
+ I << 397.4208916184016, 352.364477036168, 313.0112765842783, 278.37896394065376, 247.81424224324982,
|
|
55
|
+ 220.77925043190442, 196.76493264105017, 175.31662279234038, 156.0044839325404, 138.73983004230124,
|
|
56
|
+ 123.42064612625474, 109.82713394836259, 97.76534468972267, 87.06061858367781, 77.56000002944572, 69.1280780096311,
|
|
57
|
+ 61.64250263640252, 54.99473044877554, 49.091182970515476, 43.84634004556388, 39.184136917167976, 35.03619319797924,
|
|
58
|
+ 31.347205894128976, 28.06346770557137, 25.139117042424758, 22.53420773366429, 20.214205433283347,
|
|
59
|
+ 18.144318026099942, 16.299965972298878, 14.652633628829891, 13.184271405688083, 11.870540177313893,
|
|
60
|
+ 10.697057141915716, 9.64778948429609, 8.709338689612677, 7.871268012862094;
|
|
61
|
+ //Kern->SetPulseCurrent( VectorXr::LinSpaced( 1, 10, 200 ) ); // nbins, low, high
|
|
62
|
+ Kern->SetPulseCurrent( I ); // nbins, low, high
|
|
63
|
+
|
|
64
|
+ //VectorXr interfaces = VectorXr::LinSpaced( 41, .5, 45.5 ); // nlay, low, high
|
|
65
|
+ VectorXr interfaces = VectorXr::LinSpaced( 51, .5, 45.5 ); // nlay, low, high
|
|
66
|
+ Real thick = .5;
|
|
67
|
+ for (int ilay=1; ilay<interfaces.size(); ++ilay) {
|
|
68
|
+ interfaces(ilay) = interfaces(ilay-1) + thick;
|
|
69
|
+ thick *= 1.05;
|
|
70
|
+ }
|
|
71
|
+ Kern->SetDepthLayerInterfaces( interfaces ); // nlay, low, high
|
|
72
|
+
|
|
73
|
+ // We could, I suppose, take the earth model in here? For non-linear that
|
|
74
|
+ // may be more natural to work with?
|
|
75
|
+ std::vector<std::string> tx = {std::string("Coil 1")};
|
|
76
|
+ std::vector<std::string> rx = {std::string("Coil 2")};
|
|
77
|
+ Kern->CalculateK0( tx, rx, false );
|
|
78
|
+
|
|
79
|
+ std::ofstream dout = std::ofstream(std::string("test-")+ std::string(argv[1])+ std::string(".dat"));
|
|
80
|
+ //std::ofstream dout = std::ofstream(std::string("k-coincident.dat"));
|
|
81
|
+ dout << interfaces.transpose() << std::endl;
|
|
82
|
+ dout << I.transpose() << std::endl;
|
|
83
|
+ dout << "#real\n";
|
|
84
|
+ dout << Kern->GetKernel().real() << std::endl;
|
|
85
|
+ dout << "#imag\n";
|
|
86
|
+ dout << Kern->GetKernel().imag() << std::endl;
|
|
87
|
+ dout.close();
|
|
88
|
+
|
|
89
|
+ std::ofstream out = std::ofstream(std::string("test-")+std::string(argv[1])+std::string(".yaml"));
|
|
90
|
+ //std::ofstream out = std::ofstream(std::string("k-coincident.yaml"));
|
|
91
|
+ out << *Kern;
|
|
92
|
+ out.close();
|
|
93
|
+}
|
|
94
|
+
|
|
95
|
+
|