Files
TomoATT/test/old_tests/update_input_file/make_test_model.ipynb
2025-12-17 10:53:43 +08:00

348 lines
26 KiB
Plaintext
Raw Permalink Blame History

This file contains invisible Unicode characters
This file contains invisible Unicode characters that are indistinguishable to humans but may be processed differently by a computer. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# notebook for create init and true test model"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"depminmax 10.0 -10.0\n",
"17640\n"
]
}
],
"source": [
"import numpy as np\n",
"import math\n",
"\n",
"# grid\n",
"R_earth = 6371.0\n",
"\n",
"rr1=6361 \n",
"rr2=6381\n",
"tt1=(38.0-0.3)/180*math.pi\n",
"tt2=(42.0+0.3)/180*math.pi\n",
"pp1=(23.0-0.3)/180*math.pi\n",
"pp2=(27.0+0.3)/180*math.pi\n",
"\n",
"n_rtp = [10,50,50]\n",
"dr = (rr2-rr1)/(n_rtp[0]-1)\n",
"dt = (tt2-tt1)/(n_rtp[1]-1)\n",
"dp = (pp2-pp1)/(n_rtp[2]-1)\n",
"rr = np.array([rr1 + x*dr for x in range(n_rtp[0])])\n",
"tt = np.array([tt1 + x*dt for x in range(n_rtp[1])])\n",
"pp = np.array([pp1 + x*dp for x in range(n_rtp[2])])\n",
"\n",
"# initial model\n",
"gamma = 0.0\n",
"s0 = 1.0/6.0\n",
"slow_p=0.06\n",
"ani_p=0.04\n",
"\n",
"eta_init = np.zeros(n_rtp)\n",
"xi_init = np.zeros(n_rtp)\n",
"zeta_init = np.zeros(n_rtp)\n",
"fun_init = np.zeros(n_rtp)\n",
"vel_init = np.zeros(n_rtp)\n",
"\n",
"# true model\n",
"eta_true = np.zeros(n_rtp)\n",
"xi_true = np.zeros(n_rtp)\n",
"zeta_true = np.zeros(n_rtp)\n",
"fun_true = np.zeros(n_rtp)\n",
"vel_true = np.zeros(n_rtp)\n",
"\n",
"c=0\n",
"for ir in range(n_rtp[0]):\n",
" for it in range(n_rtp[1]):\n",
" for ip in range(n_rtp[2]):\n",
" # already initialized above\n",
" #eta_init[ir,it,ip] = 0.0\n",
" #xi_init[ir,it,ip] = 0.0\n",
" zeta_init[ir,it,ip] = gamma*math.sqrt(eta_init[ir,it,ip]**2 + xi_init[ir,it,ip]**2)\n",
" fun_init[ir,it,ip] = s0\n",
" vel_init[ir,it,ip] = 1.0/s0\n",
"\n",
" # true model\n",
" if (tt[it] >= 38.0/180.0*math.pi and tt[it] <= 42.0/180.0*math.pi \\\n",
" and pp[ip] >= 23.0/180.0*math.pi and pp[ip] <= 27.0/180.0*math.pi):\n",
" c+=1\n",
" sigma = math.sin(2.0*math.pi*(tt[it]-38.0/180.0*math.pi)/(4.0/180.0*math.pi))*math.sin(2.0*math.pi*(pp[ip]-23.0/180.0*math.pi)/(4.0/180.0*math.pi))\n",
" else:\n",
" sigma = 0.0\n",
"\n",
" if sigma < 0:\n",
" psi = 60.0/180.0*math.pi\n",
" elif sigma > 0:\n",
" psi = 120.0/180.0*math.pi\n",
" else:\n",
" psi = 0.0\n",
"\n",
" eta_true[ir,it,ip] = ani_p*abs(sigma)*math.sin(2.0*psi)\n",
" xi_true[ir,it,ip] = ani_p*abs(sigma)*math.cos(2.0*psi)\n",
" zeta_true[ir,it,ip] = gamma*math.sqrt(eta_true[ir,it,ip]**2 + xi_true[ir,it,ip]**2)\n",
" fun_true[ir,it,ip] = s0/(1.0+sigma*slow_p)\n",
" vel_true[ir,it,ip] = 1.0/fun_true[ir,it,ip] \n",
"\n",
"\n",
"#r_earth = 6378.1370\n",
"print(\"depminmax {} {}\".format(R_earth-rr1,R_earth-rr2))\n",
"print(c)\n"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"# write out in hdf5 format\n",
"import h5py\n",
"\n",
"fout_init = h5py.File('test_model_init.h5', 'w')\n",
"fout_true = h5py.File('test_model_true.h5', 'w')\n",
"\n",
"# write out the arrays eta_init, xi_init, zeta_init, fun_init, a_init, b_init, c_init, f_init\n",
"fout_init.create_dataset('eta', data=eta_init)\n",
"fout_init.create_dataset('xi', data=xi_init)\n",
"fout_init.create_dataset('zeta', data=zeta_init)\n",
"fout_init.create_dataset('vel', data=vel_init)\n",
"\n",
"# writeout the arrays eta_true, xi_true, zeta_true, fun_true, a_true, b_true, c_true, f_true\n",
"fout_true.create_dataset('eta', data=eta_true)\n",
"fout_true.create_dataset('xi', data=xi_true)\n",
"fout_true.create_dataset('zeta', data=zeta_true)\n",
"fout_true.create_dataset('vel', data=vel_true)\n",
"\n",
"fout_init.close()\n",
"fout_true.close()\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# prepare src station file\n",
"\n",
"The following code creates a src_rec_file for the inversion, which describes the source and receiver positions and arrival times.\n",
"Format is as follows:\n",
"\n",
"```\n",
" 26 1992 1 1 2 43 56.900 1.8000 98.9000 137.00 2.80 8 305644 <- src  : id_src year month day hour min sec lat lon dep_km mag num_recs id_event\n",
" 26 1 PCBI 1.8900 98.9253 1000.0000 P 10.40 18.000 <- arrival : id_src id_rec name_rec lat lon elevation_m phase epicentral_distance_km arrival_time_sec\n",
" 26 2 MRPI 1.6125 99.3172 1100.0000 P 50.84 19.400\n",
" 26 3 HUTI 2.3153 98.9711 1600.0000 P 57.84 19.200\n",
" ....\n",
"\n",
"```"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"import random\n",
"random.seed(1145141919810)\n",
"\n",
"# dummys\n",
"year_dummy = 1998\n",
"month_dummy = 1\n",
"day_dummy = 1\n",
"hour_dummy = 0\n",
"minute_dummy = 0\n",
"second_dummy = 0\n",
"mag_dummy = 3.0\n",
"id_dummy = 1000\n",
"st_name_dummy = 'AAAA'\n",
"phase_dummy = 'P'\n",
"dist_dummy = 100.0\n",
"arriv_t_dummy = 0.0\n",
"\n",
"tt1deg = tt1 * 180.0/math.pi\n",
"tt2deg = tt2 * 180.0/math.pi\n",
"pp1deg = pp1 * 180.0/math.pi\n",
"pp2deg = pp2 * 180.0/math.pi\n",
"\n",
"\n",
"n_srcs = [20,20,20]\n",
"n_src = n_srcs[0]*n_srcs[1]*n_srcs[2]\n",
"n_rec = [30 for x in range(n_src)]\n",
"\n",
"\n",
"lines = []\n",
"\n",
"nij_src = math.sqrt(n_src)\n",
"nij_rec = math.sqrt(n_rec[0])\n",
"\n",
"pos_src=[]\n",
"pos_rec=[]\n",
"\n",
"\n",
"# create receiver coordinates\n",
"elev_recs=[]\n",
"lon_recs=[]\n",
"lat_recs=[]\n",
"rec_names=[]\n",
"for i in range(n_rec[0]):\n",
" #elev_recs.append(random.uniform(-100.0,-100.0)) # elevation in m\n",
" #elev_recs.append(0) # elevation in m\n",
" #lon_recs .append(random.uniform(pp1deg*1.1,pp2deg*0.9))\n",
" #lat_recs .append(random.uniform(tt1deg*1.1,tt2deg*0.9))\n",
" rec_names.append(i)\n",
" # regularly\n",
" elev_recs.append(0.0)\n",
" tmp_ilon = i%nij_rec\n",
" tmp_ilat = int(i/nij_rec)\n",
" lon_recs.append(pp1deg + tmp_ilon*(pp2deg-pp1deg)/nij_rec)\n",
" lat_recs.append(tt1deg + tmp_ilat*(tt2deg-tt1deg)/nij_rec)\n",
"\n",
"\n",
"\n",
"# create source coordinates\n",
"for ir in range(n_srcs[0]):\n",
" for it in range(n_srcs[1]):\n",
" for ip in range(n_srcs[2]):\n",
" i_src = ir*n_srcs[1]*n_srcs[2] + it*n_srcs[2] + ip\n",
" # define one point in the domain (rr1 bottom, rr2 top)\n",
" # random\n",
" #dep = random.uniform((R_earth-rr1)*0.5,(R_earth-rr1)*0.98)\n",
" #lon = random.uniform(pp1deg,pp2deg)\n",
" #lat = random.uniform(tt1deg,tt2deg)\n",
"\n",
" # regular\n",
" dep = (R_earth-rr1)/n_srcs[0]*ir\n",
" lon = pp1deg + ip*(pp2deg-pp1deg)/n_srcs[2]\n",
" lat = tt1deg + it*(tt2deg-tt1deg)/n_srcs[1]\n",
"\n",
" # put independent name for each source\n",
" id_dummy = \"src_\"+str(i_src)\n",
" \n",
" src = [i_src, year_dummy, month_dummy, day_dummy, hour_dummy, minute_dummy, second_dummy, lat, lon, dep, mag_dummy, n_rec[i_src], id_dummy]\n",
" lines.append(src)\n",
"\n",
" pos_src.append([lon,lat,dep])\n",
"\n",
"\n",
" # create dummy station\n",
" for i_rec in range(n_rec[i_src]):\n",
" elev_rec = elev_recs[i_rec]\n",
" lon_rec = lon_recs[i_rec]\n",
" lat_rec = lat_recs[i_rec]\n",
" st_name_dummy = \"rec_\"+str(rec_names[i_rec])\n",
"\n",
" rec = [i_src, i_rec, st_name_dummy, lat_rec, lon_rec, elev_rec, phase_dummy, dist_dummy, arriv_t_dummy]\n",
" lines.append(rec)\n",
" \n",
" pos_rec.append([lon_rec,lat_rec,elev_rec])\n",
"\n",
"\n",
"# write out ev_arrivals file\n",
"fname = 'src_rec_test.dat'\n",
"\n",
"with open(fname, 'w') as f:\n",
" for line in lines:\n",
" for elem in line:\n",
" f.write('{} '.format(elem))\n",
" f.write('\\n')\n",
" "
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAbm0lEQVR4nO2dX4xuV1nGf6u0/GnagrEFj5zOGRMgpkIgpRKj1YRKFY9GTfUCU4kX4pELoaGIhGBQNE3UqBxiIsQgSu0hYuCiQggVwwn+iSnOgcOxLdhgbA9gIpJIojca8PVi9mam37enM7P2+tbez9rPL/lyZr6Z9a2ne737We96116dFBEYY4zR44qpBRhjjMnDBm6MMaLYwI0xRhQbuDHGiGIDN8YYUa6s2dn1118f29vbNbs0xhh5Lly48NWIuGH1/aoGvr29zc7OTs0ujTFGnpTS40Pvu4RijDGi2MCNMUYUG7gxxohiAzfGGFFs4MYYI4oN3JianDsH29twxRW7/547N7UiI4wNfBWFG0xBo1nn3Dk4cwYefxwidv89c2ae4+cY0yAiqr1e+tKXxqy5776Iq6+O2L29dl9XX737/lxQ0Nhz330Rp05FpLT77xw11uTUqSeOW/86dWpqZU9EJcYWFF/ATgx46vwNvOYgKdxgChojbAJDpDQ8dinNS6dCjKnEVyE0Dbz2ICncYAoaI2wCQ+ReE5X7wPG1MTQNvPYgKdxgChojNCaa2vGVOwYK94Hja6NoGnjtQVK4wRQ0jumvphFMYQI57cbozCFnDBxfG0XTwKcYpNo3WO5E06IJRNQ1AhUTGHNNak00TrY2iqaBqwySghEomEBEXSNQia9cnQrx1XqyVQhNA4/QGCQFI1AwgYj6RqAQX7ntFOLLydaR0DXwHFrPBlo1gV7j3I1AxQRc2hvWOPf4GmBZBq4ySN7tL9euphG0HF8RGqWXXmeL8TWAroG3PEje7Z+HToX4cmmvTLtcVDNw4EbgPPAI8DBwV/f+B4CL3esx4OJhnzX7gzwK2YBNoJxOlfhyaW8dhfgqyBgDPwHc3H19LfAocNPK7/we8LbDPmv2B3kUsoHWTSC3v5x2KvGlsOpyaW+jFCuhAPcDt+/7PgFfBJ5/WNvZH+TJbaew5FQwgYi6RqASXwqrLsfXRili4MA2cBm4bt97P3DQh3c/PwPsADtbW1vHU60ySN7tH9Y494lGJb4UVl0Kk0zE5LXsXEYbOHANcAG4Y+X9dwFvPMpn+CDPCgqll15ni9mmSnwp6FSYZCKWeZAHuAp4ALh75f0rgX8HTh7lc3yQZwWF0kvr2aZCfE2h87go3ANj+hPexEzAvcDZgZ+9EvjkYZ/Rv3yQZwXv9g9rnLsRqJhA7VWXwmSoEF8DjDHwW4EALu17bPB097M/BV572Gf0Lx/kKdBf6/X93HY1J0OV+FJYdbWebBXCB3layQYUJpkIjWxTyQRaXXWpTIaqGXjJV7MHebzbv45CtqliAiqrLpVkKwfVGnjJV7MHeVqvv7aabaqYgMJEo5Js9VprTTSF0DRwBRMYq/O4tGwCERr1116nQny1WtpTqO8XRNPAFUxgjM5eq3f791Cov6rEV8ulvdx2E9eyc9E0cAUTGKNTIdtsvfSiYAIK94FKsrXEgzylXs0e5Mltp5BtuvSyjkp8ubRXrj9vYjZ+kCcH7/YPa5y7EajEl0t7w23mHl8DLMvAVQbJu/3l2vkgTzmdLu2t44M8Psgz2Ma7/WWoPRkqxFduO5f2pu9vBU0DVzCBsTq9218GhclQpfTSemmv5qqrEJoGrmACtXW2bgK91hpGoGICChONk62NomngCiYwVudxadkEIlx/HUJhonGytVE0DVzBBMbo7LV6t38P11/L9efSXlmdE6Jp4Co3mMJGn4IJRGgc5PGjleX6aj3ZKoSmgUdomEBuO2eb5fpTmAx7nQrZ5nH7U4kvhWRrAF0Dz0HBBCK8239Qm7kbgYoJ1L4PFOIrt93EtXNdA2/VBHL7a90E+ra1jKBlE1C4D1pPtgqhaeAqJpDbn3f7y6EwGapkm36qavr+VtA0cJVB8m5/WZ05KEyGKtnmGJ0K8VVzVV8ITQNXMIEILzlL6uy11jACFRNQWHWpxFftVX0hNA1cwQQivNtfUqfC0yStl/Zy2qnE18SlkFw0DVzBBCK821+yncLTJC69rKMSXz7IU9HAI+ZvAr3GuRuBgglE+CDPQRodX9Pq7LW6hFIBZwPD+uZuAmP6Uyi99DoV4sulvTLtCrEsA1cwgbE6WzSBvs3cjUDFBFzaK9du4tq5roG3agK5Ols3gb5tLSNo2QQU7gOVZMsHeXyQp0g7m0C5/pRMQGHV1XKy5QzcB3mK0LoJ5PaXMwaOr3I6W0+2XAMX+IMOCtlAyyYwVudxx9zxVU7nEiZDP4Ui8Acd5p4NtGwCtXU6vsq186OVG0XTwFUGybv989B5XBxf02vsdSrElzNwgT/okIOCESiYwBidvdacssbcTUAhvsYkIwrx5Rq4wHPgERrZQKsmMEanQrbZeuklp41KfPkplIoGrmACEd7tL9lOIdt06WUdlfhSfQ4cuBE4DzwCPAzcte9nrwM+373/O4d9VrUSSm47hWzAJjCMwkEel/bK9aWQbBVkjIGfAG7uvr4WeBS4CXg58NfA07qfPfuwz6p6kCcHhWzAJlCuPyUTUJhoXNrbGMVKKMD9wO3AXwCvOE5bH+SZUGOvc+4mkNufwmrGpb2D2849vsa0K0ARAwe2gcvAdcBF4O3Ag8Ange8+oM0ZYAfY2draOp5qBROI0Ki/KpjAWJ01VzMu7U2rUSHZKshoAweuAS4Ad3TfPwT8AZCAlwH/CqQn+wwf5CnUX6smUFunigkoTDQu7W2UUQYOXAU8ANy9772PAS/f9/2/ADc82ef4IM8K3u0vq/O4OL6m19jrVIgvxRJKl2HfC5xdef+1wG90X78A+GLxDDxi/ibQa5y7ESiYwBidvdYaqxmX9sppVIkv1U1M4FYggEtd3fsicBp4KnBfV0r5NHDbYZ/lgzwreLe/nE6FbNOlvXVU4mvi0osP8szNBCK821+ynUK26dLLOirxpXqQp+TLB3kGNHq3vwytH+RpddWlEl/OwH2Qp0g7m0C5/lRMQGXV1XJ8qdbAS758kEdUY6vZpooJKNwHrcfXmHYF0DRwBROIcP21ZH81jUDFBBSeqlpCfE2IpoErmECExm6/gglEuP5aUmev1aW9PSauZeeiaeAqg6RgBAomEFHXCFTiK1enQny1nmwVQtPAIzQGScEIFEwgor4RKMRXbjuF+HKydSR0DTyH1rOBVk2g1zh3I1AxAZf2hjXOPb4GWJaBqwxS7d3+HFSyTf8/VMr1p1B66XW2GF8D6Bp4y4Pk3f556FSIL5f2yrTLxRm4wEEehWzAJlBOp0p8ubS3jkJ8FUTTwFUGSWHJqWICuf3ltFOJL4VVl0t7G0XTwFUGSWHJqWACEXWNQCW+XNor19/EmXQumgauMkje7R/WOPeJRiW+XNobbjP3+CqIpoGrDJKCESiYQIQP8pTUqRBfSyjtFUDTwCM0BknBCBRMIKK+ESjEV247hfhysnUkdA08h9azgVZNoNc4dyNQMYHWS3utJlsDLMvAVQbJu/3l2tWcDFuOrwiXXobwQZ7GD/LUMgKbQDmdSibQ6kSjoHGK/lbQNHAFExir07v9ZfBkeHDbOa9mamvMxTXwhg/ytF5/Vcg2c9qpmIDCRKOSbPVaa000hdA0cAUTGKvzuLRsAhEa9ddep0J8tVraG/PUkQ/yzDwDbz0baNUEIjTqryrx1XJpL7fdxLXsXDQNXMEExuhUyDZbL70omIDCfaCSbPkgT0UDj9Awgdx2CtmmSy/rqMSXS3vl+vMmZuMHeXLwbv+wxrkbgUp8ubQ33Gbu8TXAsgxcZZC821+unQ/ylNPp0t46PsjjgzyDbbzbX4bak6FCfOW2c2lv+v5W0DRwBRMYq9O7/WVQmAxVSi+tl/ZqrroKoWngCiZQW2frJtBrrWEEru+X68/J1kbRNHAFExir87i0bAIRGvXXXqdCfLm0N63OQmgauIIJjNHZa/Vu/x4K9VeV+HJpr6zOCdE0cAUTGKNTIdtsfaOv1fp+bZ1OtjaKpoFHaJhAbjuFbNOll3VU4qumTpX4Uki2Bsg2cOBG4DzwCPAwcFf3/q8DXwYudq/Th32WD/Ks4N3+4TZzNwKV+KqtUyG+cttNXDsfY+AngJu7r68FHgVu6gz8lw9rv/9VLQNXMIHc/lo3gb5tLSNoOb4UdKpMhq0c5AHuB26vYuAqJpDbn3f7y1F7MlSIryl0HpeW46sgRQwc2AYuA9d1Bv4YcAl4L/AtB7Q5A+wAO1tbW8dTrTJI3u0vqzMHhclQJdsco1Mhvmquugox2sCBa4ALwB3d988BngJcAdwDvPewz6j6HHgOCtlA6ybQa61hBComoLDqUomv2quuQowycOAq4AHg7gN+vg08dNjnVD2J2Wo20LIJjGmXg4oJKKy6VOJr4lJILmM2MRNwL3B25f0T+75+A/Dnh31WtRp469lAqyYQUfcGUzEBhVWXSnwt7SAPcCsQXa37m48MAn8G/FP3/l/uN/SDXtWeQnE2UK4vpd3+3LLG3E3A8TW9zl6rYgml1Kvac+DOBob1zd0ExvSnUHrpdSrEl0t7ZdoVYlkGrmACY3W2aAJ9m7kbgYoJuLRXrt3EtXNdA2/VBHJ1tm4CfdtaRtCyCSjcByrJVisHeca8fJCnQDubQLn+lExAYdXVcrLlDFzgDzooZAOtm0Bufzlj4Pgqp7P1ZMs1cB/kKULLJjBW53FvTBUTUFh1LWEy9FMoPsgzur+WTaC2ThUTUFh1OdnaKJoG7t3+g9u1aAJjdR4XFRNQWHU52doomgYe4d3+UiiYwBidvdacssbcTUAhvpxsbRRdA8+h9WyzVRMYo1Mh22y99OJka2Msy8BVZlnv9pdrp5BtuvSyjkp8+TlwH+QpotMmMIzCQR5v9JXrSyHZKoimgTvbLKsxBwUTyO1PyQQUJhonWxtD08AVTCBCo/7a65y7CeT2p7CaUajvR7i0V7pdATQNXMEEIjTqrwomMFZnzdVMq/X9XuPcJ0OFZKsgmgauYAIR3u0/SOPcJ0MVE1CYaFza2yiaBq4ySApGoGACY3UeF8fX9Bp7nQrx5RJKpadQnA2U60ul9NJrrbGacWmvnEaV+PImpsBz4BEa2UCrJjBGp0K26dLeOirxNXHpZVkGrmACEd7tL9lOIdt06WUdlfjyQZ6KJZTcdgrZgE1gmNYP8rS66lKJL2fglQ/y5KCQDdgEyvWnYgIqq66W48s1cB/kkdXYarapYgIK90Hr8TWmXQE0DVzBBCJcfy3ZX00jUDEBhaeqlhBfE6Jp4AomEKGx269gAhGuv5bU2Wt1aW+PiWvZuWgauMogKRiBgglE1DUClfjK1akQX60nW4XQNPAIjUFSMAIFE4iobwQK8ZXbTiG+nGwdCV0Dz6H1bKBVE+g1zt0IVEzApb1hjXOPrwGWZeAqg1R7tz8HlWzT/w+Vcv0plF56nQrxVaD0omvgKoOUg3f756FTIb5c2ivTLpeJ7wNNA1cZpF6rd/v3UCi9qMSXS3vrKMRXQZ2aBq4ySApLThUTyO0vp51KfCmsulzaK69zH5oGrjJICktOBROIqGsEKvHl0l65/lRWXStoGrjKIHm3f1jj3CcalfhyaW+4zdzja4zOFTQNXGWQFIxAwQQifJCnpE6F+HJp70hkGzhwI3AeeAR4GLhr5edvBAK4/rDP8kGeFRRKLyrZ5pj/vrnHV247hfhysnUkxhj4CeDm7utrgUeBm2LP3B8AHt+YgefQejbQqgn0GuduBCom0Hppr9Vka4BiJRTgfuD27usPAi8GHpuVgasMknf7y7WrORm2HF8RLr0M0cJBHmAbuAxcB/wE8M7u/QMNHDgD7AA7W1tbxxYuYQK5/Xm3f1qdrZuAwkSjoHFMf3PZxASuAS4AdwBXAw8Cz4xDDHz/q9mDPDWNwCZQTucSTGDuq5naGnOZ+D4YZeDAVV2t++7u+xcBX+mM+zHg611m/m1P9jnNHuRpvf6qkG3mtFuYCRwZhdLemGsiONGM2cRMwL3A2Sf5nc1k4AomMFbncWnZBCI06q+9ToX4arW0lxsnShPNPsYY+K3dY4KXgIvd6/TK72zGwBVMYIzOXqt3+/dQqL+qxFfLpb3cdioTzQo+yDPHQVLINlsvvbRsAq2X9nJQmWhW0DTwCA0TyG2nkG269LKOSny5tFeuv9oTzQq6Bp6DyiB5t39Y49yNQCW+XNobbjP3+BpgWQauMkje7S/Xzgd5yul0aW+dFg7yjH01/Rd5vNtfpl0utSdDhfjKbefSXrn+vIk5cxMYq9O7/WVQmAxVSi8u7Q1rVD3IU+rlgzwFaN0Eeq01jGBhJnBkXNor127qgzwlXz7IU4CWTSBCo/7a61SIL5f2yrTLxRn4zE1gjM5eq3f791Cov6rEl0t766hMNCtoGriCCYzRqZBttr7R17IJuLS3jspEs4KmgUdomEBuO4Vs06WXdVTiq6ZOx9dG0TXwHFQGqfXd/lZLLyrxVVun42tj6Bp4y4PU8m6/Suml5fhS0Nl6fOX2t4Kmgbc+SC3v9rc+GSrE1xQ6j0vL8RWx8E3MJQxSq7v9Cs9YLyG+autUiC/BVZemgSuYQET7S86aOnutNYxgYSZwZHJ0qsSX0qprH5oGrmACEd7tL6mzphEswQRq3Qcq8aWy6lpB08AVTCDCu/0l29W8wVo3gZr3gUp8qazqV9A08P4CzNkEeo1zNwIFE4jQOMijYgKOr3I6e621JpoVdA08B2cDw/rmbgJj+lMovfQ6FeLLpb0y7QqxLANXMIGxOls0gb7N3I1AxQRc2ivXrvZEs4KugbdqArk6WzeBvm0tI2jZBBTuA5Vkywd5fJCnSDubQLn+lmACLu2VYeKJRtPAFzZIR8ImcLBGH+SZXqdCfAmu6jUNXMEEIrzkLKmz11rDCBZmAkfGpb2D206UyGgauIIJRHi3v6ROhadJlEzApb1pdS46A1cwgQjv9pdsp/A0iYoJuLQ3vc5F18D7CzBnE+g1zt0IFEwgov2DPK2WXlTiS2XVtYKugefQejbQqgmM6c+ll3Vc2tPVucKyDFzBBHJ1tmwCfZu532AKGsf059JeWZ0F0DXwVk0gV2frJtC3nfNqprbGXBTug5aTrf1aF1lCUTGB3P5y2tkEyvW3BBOY+2SoEl8TP1ChaeAKJhChUX/tdc7dBHL7U1jNTGwCR8alvXLtCsWYpoErmECERv1VwQTG6qy5mhE0gSOjMBkqJFsRPsgzexOI0Hj0TcEEautcggm4tPdEVCaaFbINHLgROA88AjwM3NW9/5vAJeAi8FfAtx/2WdVq4CqDVNMIFExgrM7j4viaXmOvUyG+aq7qVxhj4CeAm7uvrwUeBW4Crtv3O68H3n3YZ1V7CsXZQLm+VEovvdYaqxkVE3Bpbx2VVf0KxUoowP3A7SvvvQV412FtZ/8ceIRGNtCqCYzRqZBturS3jkp81da5QhEDB7aBy332DdwDfBF4CLjhgDZngB1gZ2trq8p/rIQJRHi3v2Q7hWzTpZd1VOJrCp37GG3gwDXABeCOgZ+9BXj7YZ9RrYSS204hG7AJDNP6QZ5WV10q8TWxzlEGDlwFPADcfcDPt4CHDvucqgd5clDIBlzfL9ffEkzApb0yTKxzzCZmAu4Fzq68//x9X78O+OBhn+WDPBNq7HU629xjYSZwZFzam4fOfYwx8FuB2PfI4EXgNPChrvZ9Cfgw8NzDPssHeQrQcn1/rM6aqxlBEzgyCpOhQrI1RucKPsjTUjbQan2/ts4lmIBCaa/VZGuMzhU0DVxlkBSMYAnZ5nH7U4kvhVWXk63y7fahaeARGoOkYAQKk8wUOhXiK7edQmnPydaR0DXwHFrPBlqdZFR0qpiAQmlPob4fUf8+WGFZBq4ySN7tn4fO49JyfEVolF56nQrx5RLKqfkPUg7e7Z+HTsfXE1EovSjFlzcxZz5IY6i525+rb+4mkKvT8TXMFDE25/p+r1HxIE/J1+wP8ky8UXEkbALldDq+hpm43nskFOKroE5NA1cZpJrYBMrh+BpGIcYU4ivCGbjEINXGJlAGx9fBzD3GFOIrwjVwiUEyw9gEzCaZe3z1bPAplLT7szrccsstsbOzc7xG587BW98Kly/D1hbccw/ceedmBJrl4fgyAqSULkTELWvvz97AjTFm4Rxk4FdMIcYYY8x4bODGGCOKDdwYY0SxgRtjjCg2cGOMEaXqUygppf8AHq/W4ZNzPfDVqUXMCF+PPXwt9vC12GPKa3EqIm5YfbOqgc+JlNLO0GM5S8XXYw9fiz18LfaY47VwCcUYY0SxgRtjjChLNvA/mlrAzPD12MPXYg9fiz1mdy0WWwM3xhh1lpyBG2OMNDZwY4wRZREGnlJ6ekrpUymlz6aUHk4pvb17/wdTSp9OKV1MKf1dSul5U2vdNE9yLW7rrsVDKaX3pZSunFprLVJKT0kpfSal9JHu++9IKT2YUvpCSukDKaWnTq2xFgPX4pe66xAppeun1leTgWtxLqX0z9098t6U0lVTa1yEgQP/A9wWES8GXgK8MqX0PcC7gDsj4iXA+4FfnUxhPYauxfcC7wNeFREvZPew1c9NJ7E6dwGf2/f9bwPviIjnAf8J/PwkqqZh9Vr8PfAK5nMAryar1+Ic8J3Ai4BnAK+ZQtR+FmHg3R+1+O/u26u6V3Sv67r3nwn82wTyqnLAtfgG8L8R8Wj3/seBn5pCX21SSieBHwXe032fgNuAD3a/8j7gJycRV5nVawEQEZ+JiMcmEzURB1yLj+77CzmfAk5Opa9nEQYO31wOXQS+Anw8Ih5kdwb9aErpS8Crgd+aUGI1Vq8Fu8F4ZUqpP2X208CNE8mrzVngV4D/677/VuBrEfH17vsvAc+dQNcUnOWJ12LJnOWAa9GVTl4NfKyypjUWY+AR8Y2uVHISeFlK6YXAG4DTEXES+BPg9yeUWI3VawF8F/Aq4B0ppU8B/8VuVt40KaUfA74SERem1jI1vhZ7HOFa/CHwNxHxtxVlDbKYjaqeiPhaSuk88CPAi7tMHOADzGBGrcm+a/HKiPhd4PsBUko/BLxgUnF1+D7gx1NKp4Gns1tOeyfwrJTSlV0WfhL48oQaa7F2LVJK90XEz06sawoOvBYppV8DbgB+cVKFHYvIwFNKN6SUntV9/QzgdnY3J56ZUuqNqn+vaQ64Fp9PKT27e+9pwJuBd08mshIR8ZaIOBkR2+yuQD4REXcC59ktI8HuZu79E0msxgHXYonmfeC1SCm9Bvhh4GciYhZlpkUYOHACOJ9SugT8I7s18I8AvwB8KKX0WXZrWm+aUGMtDroWb0opfQ64BHw4Ij4xpciJeTNwd0rpC+zWxP94Yj2TkVJ6fbdHdBK4lFJ6z2FtGubdwHOAf+gePX7b1IJ8lN4YY0RZSgZujDHNYQM3xhhRbODGGCOKDdwYY0SxgRtjjCg2cGOMEcUGbowxovw/g47aIBo60iUAAAAASUVORK5CYII=",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"# draw src and rec positions\n",
"import matplotlib.pyplot as plt\n",
"\n",
"for i_src in range(n_src):\n",
" plt.scatter(pos_src[i_src][1],pos_src[i_src][0],c='r',marker='o')"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAQlUlEQVR4nO3db2zd113H8c8nNGMrLSsQt5T8qaeNDoWhtugyVUR90DDKCIgOhqCTF/ZgI4BWlo6qm1gQf4QijWmk3RM2maWoCIOKSARiqiYqLRoUoRQ7S9OmLrAHyZQSVleia1ClTm2+PLg/r67r6/uzfa5/59zzfklWfH/3j49Oez755vje83VECABQni1dDwAAsD4EOAAUigAHgEIR4ABQKAIcAAp1xWb+sG3btsXk5ORm/kgAKN7c3NzzETGx/PqmBvjk5KRmZ2c380cCQPFsn1/pOlsoAFAoAhwACkWAA0ChCHAAKNTQX2La3inpLyVdJykkTUfE52w/LOmdzcOukfRCRNw8onECAJZpU4G/IuneiNgt6VZJH7W9OyJ+NSJubkL7mKTjIxwngDEwMyNNTkpbtvT/nJnpekRlG1qBR8RFSReb7y/Znpe0XdLTkmTbkn5F0t4RjhNA4WZmpAMHpJde6t8+f75/W5KmprobV8nWtAdue1LSLZJOLrl8m6RvRsR/DXjOAduztmcXFhbWPVCMNyqz8Xfo0Gvhveill/rXsT6tA9z2VepvldwTES8uuesDkv5m0PMiYjoiehHRm5h4wweJgO9UZufPSxGvVWaE+Hj5xjfWdh3DtQpw21vVD++ZiDi+5PoVkn5J0sOjGR6VWQ2ozOqwa9farmO4oQHe7HEflTQfEUeW3f0eSc9ExIVRDI7KrA5UZnU4fFi68srXX7vyyv51rE+bCnyPpP2S9to+3Xzta+67S6tsn2wUlVkdqMzqMDUlTU9LN9wg2f0/p6f5BeZGeDN7YvZ6vVjLYVZbtvQr7+Vs6fLlhANDp5a/O0HqV2YsbqDP9lxE9JZfz/qTmFRmdaAyA9Yn6wBnz6weU1PSuXP9f1mdO0d4A21kHeBUZgAw2KY2dFiPqSkCGwBWknUFDgAYjAAHgEIR4ABQKAIcAApFgANAoQhwACgUAQ4AhSLAAaBQBDgAFIoABzpCsxJsFAGeIRb2+KNZCVIgwDPDwq4DzUqQAgGeGRZ2HWgjhxQI8MywsOtAsxKkQIBnhoVdB5qVIAUCPDMs7DrQrAQpZN/QoTaLC/jQof62ya5d/fBmYY8fmpVgowjwDLGwAbTBFgoAFIoAB4BCEeAAUCgCHAAKRYADQKEIcAAoFAEOAIUiwAGgUAQ4ABSKAAewaWhWkhYfpQewKRablSyed7/YrETi6Ij1ogJHFqjMxh/NStKjAkfnqMzqQLOS9LKvwKnMxh+VWR1oVpJe1gFOg986UJnVgWYl6Q0NcNs7bZ+w/bTts7YPLrnvt20/01z/TOrBUZnVgcqsDnQhSq/NHvgrku6NiFO2r5Y0Z/tRSddJulPSTRHxsu1rUw+OyqwOhw+/fg9cojIbVzQrSWtoBR4RFyPiVPP9JUnzkrZL+i1Jn46Il5v7nks9OCqzOlCZAeuzpj1w25OSbpF0UtKNkm6zfdL2V23/xIDnHLA9a3t2YWFhTYNjz6weU1PSuXPS5cv9PwlvYLjWAW77KknHJN0TES+qv/3y/ZJulXSfpL+17eXPi4jpiOhFRG9iYmJNg6MyA4DBWr0P3PZW9cN7JiKON5cvSDoeESHpcduXJW2TtLYyewj2zABgZW3ehWJJRyXNR8SRJXf9vaTbm8fcKOlNkp4fwRgBACtoU4HvkbRf0pO2TzfXPiXpQUkP2n5K0rclfaipxgEAm2BogEfEY5LesLfd+GDa4QAA2sr6k5gAgMEIcAAoFAEOAIUiwAGgUAQ4ABSKAAc6wln32Cg68gAdoAsRUqACzxCV2fjjrHukQAWeGSqzOnDWPVKgAs8MlVkdOOseKRDgmaEyqwNn3SMFAjwzVGZ14Kx7pECAZ4bKrB50IcJGEeCZoTID0BbvQskQXYgAtEEFDgCFIsABoFAEOAAUigAHgEIR4ABQKAIcAApFgANAoQhwACgUAQ4AhSLAAWwampWkxUfpAWwKmpWkRwWOLFCZjT+alaSXfYCzsMffYmV2/rwU8Vplxn/r8UKzkvSyDnAWdh2ozOpAs5L0sg5wFnYdqMzqQLOS9LIOcBZ2HajM6kCzkvSyDnAWdh2ozOpBG7m0sg5wFnYdqMyA9cn6feCLC/jQof62ya5d/fBmYY8f2sgBa5d1gEssbAAYJOstFADAYAQ4ABRqaIDb3mn7hO2nbZ+1fbC5/oe2n7V9uvnaN/rhAgAWtdkDf0XSvRFxyvbVkuZsP9rcd39EfHZ0wwMADDI0wCPioqSLzfeXbM9L2j7qgQEAVremPXDbk5JukXSyuXS37TO2H7T9fakHBwAYrHWA275K0jFJ90TEi5I+L+ntkm5Wv0L/0wHPO2B71vbswsLCxkcMAJDUMsBtb1U/vGci4rgkRcQ3I+LViLgs6c8lvXul50bEdET0IqI3MTGRatwAUL0270KxpKOS5iPiyJLr1y952C9Keir98AAAg7SpwPdI2i9p77K3DH7G9pO2z0i6XdLHRzlQYNzQrAQb1eZdKI9J8gp3PZJ+OEAd6A+JFPgkZoaozMYfzUqQQvaHWdWGyqwONCtBClTgmaEyqwPNSpACAZ4ZKrM60KwEKRDgmaEyqwNdiJACAZ4ZKrN60B8SG0WAZ4bKDEBbvAslQ7SRA9AGFTgAFIoAB4BCEeAAUCgCHAAKRYADQKEIcAAoFAEOAIUiwAGgUAQ4ABSKAAewaWhWkhYBjiywsMffYrOS8+eliNealfDfev0IcHSOhV0HmpWkl32AU5mNPxZ2HWhWkl7WAU5lVgcWdh1oVpJe1gFOZVYHFnYdaFaSXtYBTmVWBxZ2HWhWkl7WAU5lVgcWdj1oI5dW1gFOZVYPFjawdlkHOJUZAAyWfU9M+kMCwMqyrsABAIMR4ABQKAIcAApFgANAoQhwACgUAQ4AhSLAAaBQBDgAjNAoj8TO/oM8AFCqxSOxF09VXTwSW0rzAUUqcKAjNCsZf6M+EntogNveafuE7adtn7V9cNn999oO29vSDAkYfzQrqcOoj8RuU4G/IuneiNgt6VZJH7W9W+qHu6Q7JHFCd0JUZuOPZiV1GPWR2EMDPCIuRsSp5vtLkuYlbW/uvl/SJyRFmuGAyqwONCupw6iPxF7THrjtSUm3SDpp+05Jz0bEE0Oec8D2rO3ZhYWF9Y+0ElRmdaBZSR1GfSS2I9oVz7avkvRVSYclfVnSCUl3RMS3bJ+T1IuI51d7jV6vF7Ozsxsb8ZjbsqVfeS9n95sdYDwsf3eC1K/MOO8eK7E9FxG95ddbVeC2t0o6JmkmIo5Lerukt0l6ognvHZJO2f7BdEOuE5VZHWhWghTavAvFko5Kmo+II5IUEU9GxLURMRkRk5IuSPrxiPifkY62ArSRqwdt5LBRbSrwPZL2S9pr+3TztW/E46oWlRmAtlrvgafAHjgArN2G9sABAPkhwAGgUAQ4ABSKAAeAQhHgAFAoAhwACkWAA0ChCHAAKBQBDgCFIsABbBqalaRFgCMLLOzxR7OS9AhwdI6FXQealaSXfYBTmY0/FnYdaCOXXtYBTmVWBxZ2HWhWkl7WAU5lVgcWdh1oVpJe1gFOZVYHFnYdaFaSXtYBTmVWBxZ2PWgjl1bWAU5lVg8WNrB2WQc4lRkADHZF1wMYZmqKwAaAlWRdgQMABiPAAaBQBDgAFIoAB4BCEeAAUCgCHAAKRYADQKEIcAAYoVEeiZ39B3kAoFSLR2Ivnqq6eCS2lOYDilTgQEdoVjL+Rn0kNhU40IFRV2bIw6iPxKYCzxCV2fijWUkdRn0kNgGeGdrI1YFmJXUY9ZHYBHhmqMzqQLOSOoz6SGwCPDNUZnWgWUk9RtmshADPDJVZHWhWghQI8MxQmdWDNnLYqKEBbnun7RO2n7Z91vbB5vof2z5j+7Ttf7L9Q6Mf7vijMgPQliNi9QfY10u6PiJO2b5a0pyk90m6EBEvNo/5mKTdEfGbq71Wr9eL2dnZJAMHgFrYnouI3vLrQyvwiLgYEaea7y9Jmpe0fTG8G98jafW/CQAASa3pk5i2JyXdIulkc/uwpF+T9C1Jtw94zgFJByRpF7+JA4BkWv8S0/ZVko5Jumex+o6IQxGxU9KMpLtXel5ETEdELyJ6ExMTKcYMAFDLALe9Vf3wnomI4ys8ZEbS+1MODACwujbvQrGko5LmI+LIkus/vORhd0p6Jv3wAACDtNkD3yNpv6QnbZ9urn1K0odtv1PSZUnnJa36DhQAQFpDAzwiHpPkFe56JP1wAABt8UlMACgUAQ5g03DWfVp05AGwKehClB4VOLJAZTb+OOs+PSpwdI7KrA6cdZ9e9hU4ldn4ozKrA2fdp5d1gNMfsg5UZnXgrPv0sg5wKrM6UJnVgbPu08s6wKnM6kBlVg+6EKWVdYBTmdWBygxYn6wDnMqsHlRmwNplHeBUZgAwWPbvA5+aIrABYCVZV+AAgMEIcAAoFAEOAIUiwAGgUAQ4ABTKEbF5P8xeUL9/Zk62SXq+60FkjPkZjjlaHfMz3LA5uiEiJpZf3NQAz5Ht2YjodT2OXDE/wzFHq2N+hlvvHLGFAgCFIsABoFAEuDTd9QAyx/wMxxytjvkZbl1zVP0eOACUigocAApFgANAoaoJcNtvtv247Sdsn7X9R831n7J9yvZp24/ZfkfXY+3KKnO0t5mjp2w/ZDv7UyxHyfZ32f6a7S81t99m+6Ttr9t+2Pabuh5jl1aYn7ubuQnb27oeXw5WmKMZ2//RrLEHbW9t8zrVBLiklyXtjYibJN0s6b22b5X0eUlTEXGzpL+W9HudjbB7K83RT0p6SNJdEfEu9T+I9aHuhpiFg5Lml9z+E0n3R8Q7JP2vpA93Mqp8LJ+ff5X0HuX3Ib4uLZ+jGUk/IunHJL1F0kfavEg1AR59/9fc3Np8RfP1vc31t0r67w6Gl4UBc/SqpG9HxH821x+V9P4uxpcD2zsk/ZykLza3LWmvpL9rHvKQpPd1MrgMLJ8fSYqIr0XEuc4GlZkBc/RIs/5C0uOSdrR5rWoCXPrOP1tOS3pO0qMRcVL9v+kesX1B0n5Jn+5wiJ1bPkfq/890he3FT4n9sqSdHQ0vBw9I+oSky83tH5D0QkS80ty+IGl7B+PKxQN6/fzgjR7QgDlqtk72S/pymxeqKsAj4tVmq2SHpHfbfpekj0vaFxE7JP2FpCMdDrFzy+dI0o9KukvS/bYfl3RJ/aq8OrZ/XtJzETHX9VhyxPwM12KO/kzSP0fEv7R5vSp/GRURL9g+IelnJd3UVOKS9LBa/s037pbM0Xsj4rOSbpMk23dIurHTwXVnj6RfsL1P0pvV33r7nKRrbF/RVOE7JD3b4Ri79Ib5sf1XEfHBjseVk4FzZPsPJE1I+o22L1ZNBW57wvY1zfdvkfTT6v8S4a22FwNp8VqVBszRM7avba59t6RPSvpCZ4PsUET8bkTsiIhJ9f9V8pWImJJ0Qv2tJan/C95/6GiInRowP4T3EoPmyPZHJP2MpA9EROvtp2oCXNL1kk7YPiPp39XfA/+SpF+XdMz2E+rvPd3X4Ri7NmiO7rM9L+mMpH+MiK90OcgMfVLS79j+uvp74kc7Hk9WbH+s+R3TDklnbH9x2HMq9AVJ10n6t+Ytzb/f5kl8lB4AClVTBQ4AY4UAB4BCEeAAUCgCHAAKRYADQKEIcAAoFAEOAIX6fyDuGKLLtQk6AAAAAElFTkSuQmCC",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"# plot receivers\n",
"for i_rec in range(n_rec[0]):\n",
" plt.scatter(pos_rec[i_rec][1],pos_rec[i_rec][0],c='b',marker='o')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.8 (main, Dec 18 2022, 17:23:16) [GCC 12.2.0]"
},
"orig_nbformat": 4,
"vscode": {
"interpreter": {
"hash": "c8112a80bcc57e66f52aa6c921b20667122cfa5e8b1ad09e9dfbd3a3ba336bd6"
}
}
},
"nbformat": 4,
"nbformat_minor": 2
}