Skip to content
GitLab
Commit 6722fafe authored by Guillaume Blanchy's avatar Guillaume Blanchy
Browse files

Move sequence generator to sequence.py and refactor the code

parent bd2f2bf9
master v2024.0.35 v2024.0.34 v2024.0.33 v2024.0.32 v2024.0.31 v2024.0.30 v2024.0.29 v2024.0.28 v2024.0.27 v2024.0.26 v2024.0.25 v2024.0.24
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 172 additions and 411 deletions
+172 -411
ohmpi/ohmpi.py
+ 13
− 29
View file @ 6722fafe
...@@ -251,47 +251,31 @@ class OhmPi(object): ...@@ -251,47 +251,31 @@ class OhmPi(object):
Each tuple is the form (<array_name>, param1, param2, ...) Each tuple is the form (<array_name>, param1, param2, ...)
Dipole spacing is specified in terms of "number of electrode spacing". Dipole spacing is specified in terms of "number of electrode spacing".
Dipole spacing is often referred to 'a'. Number of levels is a multiplier Dipole spacing is often referred to 'a'. Number of levels is a multiplier
of 'a', often refere to 'n'. For multigradient array, an additional parameter of 'a', often referred to 'n'. For multigradient array, an additional parameter
's' is needed. 's' is needed.
Types of sequences available are : Types of sequences available are :
- ('wenner', a) - ('wenner', a)
- ('dpdp', a, n) - ('dpdp', a, n)
- ('schlum', a, n) - ('schlum', a, n)
- ('multigrad', a, n, s) - ('multigrad', a, n, s)
ireciprocal : bool, optional By default, if an integer is provided for a, n and s, the parameter
will be considered varying from 1 to this value. For instance, for
('wenner', 3), the sequence will be generated for a = 1, a = 2 and a = 3.
If only some levels are desired, the user can use a list instead of an int.
For instance ('wenner', [3]) will only generate quadrupole for a = 3.
include_reciprocal : bool, optional
If True, will add reciprocal quadrupoles (so MNAB) to the sequence. If True, will add reciprocal quadrupoles (so MNAB) to the sequence.
ip_optimization: bool, optional
If True, will optimize for induced polarization measurement (i.e. will
try to put as much time possible between injection and measurement at
the same electrode).
fpath : str, optional fpath : str, optional
Path where to save the sequence (including filename and extension). By Path where to save the sequence (including filename and extension). By
default, sequence is saved in ohmpi/sequences/sequence.txt. default, sequence is saved in ohmpi/sequences/sequence.txt.
""" """
# dictionnary of function to create sequence dfseq = create_sequence(nelec, params=params, include_reciprocal=include_reciprocal)
fdico = { self.sequence = dfseq.astype(int).values
'dpdp': dpdp1,
'wenner': wenner,
'schlum': schlum1,
'multigrad': multigrad,
}
# check arguments
if fpath is None:
fpath = os.path.join(os.path.dirname(__file__), '../sequences/sequence.txt')
qs = []
# create sequence
for p in params:
pok = [int(p[i]) for i in np.arange(1, len(p))] # make sure all are int
qs.append(fdico[p[0]](nelec, *pok).values.astype(int))
quad = np.vstack(qs)
if len(quad.shape) == 1: # only one quadrupole
quad = quad[None, :]
# add reciprocal
if ireciprocal:
quad = np.r_[quad, quad[:, [2, 3, 0, 1]]]
self.sequence = quad
# save sequence
np.savetxt(fpath, self.sequence, delimiter=' ', fmt='%d') np.savetxt(fpath, self.sequence, delimiter=' ', fmt='%d')
print('{:d} quadrupoles generated.'.format(self.sequence.shape[0]))
@staticmethod @staticmethod
def _find_identical_in_line(quads): def _find_identical_in_line(quads):
......
ohmpi/sequence.py
+ 159
− 382
View file @ 6722fafe
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Nov 14 16:05:59 2018
@author: pmclachlan
originally produced for the ResIPy project (https://gitlab.com/hkex/resipy)
"""
import numpy as np import numpy as np
import pandas as pd import pandas as pd
import os
def create_sequence(nelec, params=[('dpdp', 1, 8)], include_reciprocal=False, ip_optimization=False):
"""Creates a sequence of quadrupole.Several type of
sequence or sequence with different parameters can be combined together.
def dpdp1(elec_num, a, n):
''' Generates quadrupole matrix for dipole-dipole survey.
Parameters Parameters
---------- ----------
elec_num : int nelec : int
Number of electrodes Number of electrodes.
a : int or list of int params : list of tuple, optional
Spacing between C and V pairs in electrode spacing Each tuple is the form (<array_name>, param1, param2, ...)
(a = 1 is the same as skip 0). Dipole spacing is specified in terms of "number of electrode spacing".
n : int or list of int Dipole spacing is often referred to 'a'. Number of levels is a multiplier
Quadrupole separation in electrode spacing. of 'a', often referred to 'n'. For multigradient array, an additional parameter
's' is needed.
Notes Types of sequences available are :
----- - ('wenner', a)
Length of `a` should match `n`. - ('dpdp', a, n)
''' - ('schlum', a, n)
elec_id = np.arange(elec_num)+1 - ('multigrad', a, n, s)
abmn = [] By default, if an integer is provided for a, n and s, the parameter
if isinstance(a, list) is False: will be considered varying from 1 to this value. For instance, for
if a < 1: ('wenner', 3), the sequence will be generated for a = 1, a = 2 and a = 3.
raise ValueError('a must be >= 1, it\'s the electrode spacing between C and V pairs (a = 1 is the same as skip 0)') If only some levels are desired, the user can use a list instead of an int.
n_ = np.array(range(0,n))+1 For instance ('wenner', [3]) will only generate quadrupole for a = 3.
for j in np.array(range(0,len(n_))): include_reciprocal : bool, optional
A = elec_id If True, will add reciprocal quadrupoles (so MNAB) to the sequence.
B = A + a ip_optimization: bool, optional
M = B + n_[j] * a If True, will optimize for induced polarization measurement (i.e. will
N = M + a try to put as much time possible between injection and measurement at
abmn.append(np.vstack([A, B, M, N]).T) the same electrode).
"""
else: # dictionary of function to create sequence
for i in np.array(range(0,len(a))): fdico = {
if a[i] < 1: 'dpdp': dpdp,
raise ValueError('a must be >= 1, it\'s the electrode spacing between C and V pairs (a = 1 is the same as skip 0)') 'wenner': wenner,
n_ = np.array(range(0,n[i]))+1 'schlum': schlum,
for j in np.array(range(0,len(n_))): 'multigrad': multigrad,
A = elec_id }
B = A + a[i]
M = B + n_[j] * a[i] # check parameters (all int and > 0)
N = M + a[i] for param in params:
abmn.append(np.vstack([A, B, M, N]).T) for i, p in enumerate(param[1:]):
if isinstance(p, float) or isinstance(p, int):
proto_mtrx = pd.DataFrame(np.vstack(abmn)) if p <= 0:
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,3] <= elec_num] raise ValueError('parameters of sequence (a, n or s) needs to be > 0')
return proto_mtrx else:
for pp in p:
if pp <= 0:
def dpdp2(elec_num, a, n): raise ValueError('parameters of sequence (a, n or s) needs to be > 0')
# create sequence
dfs = []
for param in params:
dfs.append(fdico[param[0]](nelec, *param[1:]))
dfseq = pd.concat(dfs, axis=0)
# add reciprocal
if include_reciprocal:
dfseqr = dfseq.copy().rename(columns={'a': 'm', 'b': 'n', 'm': 'a', 'n': 'b'})
dfseqr = dfseqr.sort_values(['a', 'b', 'm', 'n'])
dfseq = pd.concat([dfseq, dfseqr]).reset_index(drop=True)
# optimize for IP
if ip_optimization:
print('NOT IMPLEMENTED YET')
# stats
print('{:d} quadrupoles generated.'.format(dfseq.shape[0]))
return dfseq
def dpdp(nelec, a, n):
''' Generates quadrupole matrix for dipole-dipole survey. ''' Generates quadrupole matrix for dipole-dipole survey.
Parameters Parameters
---------- ----------
elec_num : int nelec : int
Number of electrodes Number of electrodes
a : int or list of int a : int or list of int
Spacing between C and V pairs in electrode spacing Spacing between AB (current) and MN (voltage) pairs in electrode spacing
(a = 1 is the same as skip 0). (a = 1 is the same as skip 0).
n : int or list of int n : int or list of int
Quadrupole seperation in electrode spacing. Quadrupole separation in electrode spacing.
Notes
-----
Length of `a` should match `n`.
'''
elec_id = np.arange(elec_num)+1
abmn = []
if isinstance(a, list) is False:
if a < 1:
raise ValueError('a must be >= 1, it\'s the electrode spacing between C and V pairs (a = 1 is the same as skip 0)')
n_ = np.array(range(0,n))+1
for j in np.array(range(0,len(n_))):
A = elec_id
B = A + a
M = B + n_[j]
N = M + a
abmn.append(np.vstack([A, B, M, N]).T)
else:
for i in np.array(range(0,len(a))):
if a[i] < 1:
raise ValueError('a must be >= 1, it\'s the electrode spacing between C and V pairs (a = 1 is the same as skip 0)')
n_ = np.array(range(0,n[i]))+1
for j in np.array(range(0,len(n_))):
A = elec_id
B = A + a[i]
M = B + n_[j]
N = M + a[i]
abmn.append(np.vstack([A, B, M, N]).T)
proto_mtrx = pd.DataFrame(np.vstack(abmn))
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,3] <= elec_num]
return proto_mtrx
def wenner(elec_num, n):
''' Generate quadrupoles matrix for Wenner alpha survey for a = 1 ... n.
''' '''
seq = [] # check parameters
for a in np.arange(n) + 1: if isinstance(a, int):
seq.append(wenner_alpha(elec_num, a)) a = np.arange(a) + 1
return pd.concat(seq).reset_index(drop=True) if isinstance(n, int):
n = np.arange(n) + 1
def wenner_alpha(elec_num, a):
''' Generates quadrupole matrix for Wenner alpha survey.
Parameters
----------
elec_num : int
Number of electrodes
a : int or list of int
Spacing between electrodes (in electrode spacing).
'''
elec_id = np.arange(elec_num)+1 elec = np.arange(nelec) + 1
abmn = [] abmn = []
if isinstance(a, list) is False: for aa in a:
A = elec_id if aa < 1:
M = A + a raise ValueError('a must be >= 1, it is the electrode spacing between AB and MN pairs (a = 1 is the same as skip 0)')
N = M + a for nn in n:
B = N + a if nn < 1:
abmn.append(np.vstack([A, B, M, N]).T) raise ValueError('n must be >= 1, it is the number of level between AB and MN')
A = elec
else: B = A + aa
for i in np.array(range(0,len(a))): M = B + nn * aa
A = elec_id N = M + aa
M = A + a[i]
N = M + a[i]
B = N + a[i]
abmn.append(np.vstack([A, B, M, N]).T) abmn.append(np.vstack([A, B, M, N]).T)
abmn = np.vstack(abmn)
proto_mtrx = pd.DataFrame(np.vstack(abmn)) abmn = abmn[(abmn <= nelec).all(1), :]
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num] df = pd.DataFrame(abmn, columns=['a', 'b', 'm', 'n'])
return(proto_mtrx) df = df.sort_values(['a', 'b', 'm', 'n']).reset_index(drop=True)
return df
def wenner_beta(elec_num, a): def wenner(nelec, a):
''' Generates quadrupole matrix for Wenner beta survey. '''Generates quadrupole matrix for Wenner alpha survey.
Parameters Parameters
---------- ----------
elec_num : int nelec : int
Number of electrodes Number of electrodes
a : int or list of int a : int or list of int
Spacing between electrodes (in electrode spacing). Spacing between electrodes (in electrode spacing).
''' '''
if isinstance(a, int):
elec_id = np.arange(elec_num)+1 a = np.arange(a) + 1
elec = np.arange(nelec) + 11
abmn = [] abmn = []
for aa in a:
if type(a) == int: A = elec
B = elec_id M = A + aa
A = B + a N = M + aa
M = A + a B = N + aa
N = M + a
abmn.append(np.vstack([A, B, M, N]).T) abmn.append(np.vstack([A, B, M, N]).T)
abmn = np.vstack(abmn)
if type(a) == list: abmn = abmn[(abmn <= nelec).all(1), :]
for i in np.array(range(0,len(a))): df = pd.DataFrame(abmn, columns=['a', 'b', 'm', 'n'])
B = elec_id df = df.sort_values(['a', 'b', 'm', 'n']).reset_index(drop=True)
A = B + a[i] return df
M = A + a[i]
N = M + a[i]
abmn.append(np.vstack([A, B, M, N]).T)
proto_mtrx = pd.DataFrame(np.vstack(abmn))
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num]
return(proto_mtrx)
def wenner_gamma(elec_num, a):
''' Generates quadrupole matrix for Wenner gamma survey.
Parameters
----------
elec_num : int
Number of electrodes
a : int or list of int
Spacing between electrodes (in electrode spacing).
'''
elec_id = np.arange(elec_num)+1
abmn = []
if isinstance(a, list) is False:
A = elec_id
M = A + a
B = M + a
N = B + a
abmn.append(np.vstack([A, B, M, N]).T)
else:
for i in np.array(range(0,len(a))):
A = elec_id
M = A + a[i]
B = M + a[i]
N = B + a[i]
abmn.append(np.vstack([A, B, M, N]).T)
proto_mtrx = pd.DataFrame(np.vstack(abmn))
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num]
return(proto_mtrx)
def schlum(nelec, a, n):
def schlum1(elec_num, a, n):
''' Generates quadrupole matrix for Schlumberger survey. ''' Generates quadrupole matrix for Schlumberger survey.
Parameters Parameters
---------- ----------
elec_num : int nelec : int
Number of electrodes Number of electrodes
a : int or list of int a : int or list of int
Spacing between electrodes (in electrode spacing). Spacing between electrodes (in electrode spacing).
n : int or list of int n : int or list of int
Quadrupole seperation in electrode spacing. Quadrupole separation in electrode spacing.
''' '''
elec_id = np.arange(elec_num)+1 # check parameters
abmn = [] if isinstance(a, int):
if isinstance(a, list) is False: a = np.arange(a) + 1
n_ = np.array(range(0,n))+1 if isinstance(n, int):
for j in np.array(range(0,len(n_))): n = np.arange(n) + 1
if n_[j] >= a:
#print(n_[j],a)
A = elec_id
M = A + n_[j]
N = M + a
B = N + n_[j]
abmn.append(np.vstack([A, B, M, N]).T)
else:
for i in np.array(range(0,len(a))):
n_ = np.array(range(0,n[i]))+1
for j in np.array(range(0,len(n_))):
if n_[j] >= a[i]-1: #minus 1 here to permit measurements on level 2 of pseudosection
A = elec_id
M = A + n_[j]
N = M + a[i]
B = N + n_[j]
abmn.append(np.vstack([A, B, M, N]).T)
proto_mtrx = pd.DataFrame(np.vstack(abmn))
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num]
return(proto_mtrx)
#def schlum1(elec_num, a, n):
# ''' Generates quadrupole matrix for Schlumberger survey.
# Parameters
# ----------
# elec_num : int
# Number of electrodes
# a : int or list of int
# Spacing between electrodes (in electrode spacing).
# n : int or list of int
# Quadrupole seperation in electrode spacing.
# '''
# elec_id = np.arange(elec_num)+1
# proto_mtrx = pd.DataFrame([])
# if isinstance(a, list) is False:
# n_ = np.array(range(0,n))+1
# for j in np.array(range(0,len(n_))):
# if n_[j] >= a:
# #print(n_[j],a)
# A = elec_id
# M = A + n_[j]
# N = M + a
# B = N + n_[j]
# # if j == 0:
# # proto_mtrx = pd.DataFrame(np.column_stack((A, B, M, N)))
# # if j > 0:
# new_rows = pd.DataFrame(np.column_stack((A, B, M, N)))
# proto_mtrx = proto_mtrx.append(new_rows)
# else:
# for i in np.array(range(0,len(a))):
# n_ = np.array(range(0,n[i]))+1
# for j in np.array(range(0,len(n_))):
# if n_[j] >= a[i]:
# A = elec_id
# M = A + n_[j]
# N = M + a[i]
# B = N + n_[j]
#
# new_rows = pd.DataFrame(np.column_stack((A, B, M, N)))
# proto_mtrx = proto_mtrx.append(new_rows)
# proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num]
# return(proto_mtrx)
elec = np.arange(nelec) + 1
#def schlum1(elec_num, a, n):
# ''' Generates quadrupole matrix for Schlumberger survey.
# Parameters
# ----------
# elec_num : int
# Number of electrodes
# a : int or list of int
# Spacing between electrodes (in electrode spacing).
# n : int or list of int
# Quadrupole seperation in electrode spacing.
# '''
# elec_id = np.arange(elec_num)+1
# if isinstance(a, list) is False:
# n_ = np.array(range(0,n))+1
# for j in np.array(range(0,len(n_))):
# A = elec_id
# M = A + n_[j] * a
# N = M + a
# B = N + n_[j] * a
# if j == 0:
# proto_mtrx = pd.DataFrame(np.column_stack((A, B, M, N)))
# if j > 0:
# new_rows = pd.DataFrame(np.column_stack((A, B, M, N)))
# proto_mtrx = proto_mtrx.append(new_rows)
#
# else:
# for i in np.array(range(0,len(a))):
# n_ = np.array(range(0,n[i]))+1
# for j in np.array(range(0,len(n_))):
# A = elec_id
# M = A + n_[j] * a[i]
# N = M + a[i]
# B = N + n_[j] * a[i]
# if (i + j) == 0:
# proto_mtrx = pd.DataFrame(np.column_stack((A, B, M, N)))
# if (i + j) > 0:
# new_rows = pd.DataFrame(np.column_stack((A, B, M, N)))
# proto_mtrx = proto_mtrx.append(new_rows)
#
# proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num]
# return(proto_mtrx)
def schlum2(elec_num, a, n):
''' Generates quadrupole matrix for Schlumberger survey.
Parameters
----------
elec_num : int
Number of electrodes
a : int or list of int
Spacing between electrodes (in electrode spacing).
n : int or list of int
Quadrupole seperation in electrode spacing.
'''
elec_id = np.arange(elec_num)+1
abmn = [] abmn = []
if isinstance(a, list) is False: for aa in a:
n_ = np.array(range(0,n))+1 for nn in n:
for j in np.array(range(0,len(n_))): if nn >= aa-1: #minus 1 here to permit measurements on level 2 of pseudosection
A = elec_id A = elec
M = A + n_[j] * a M = A + nn
N = M + a N = M + aa
B = N + n_[j] * a B = N + nn
abmn.append(np.vstack([A, B, M, N]).T)
else:
for i in np.array(range(0,len(a))):
n_ = np.array(range(0,n[i]))+1
for j in np.array(range(0,len(n_))):
A = elec_id
M = A + n_[j]
N = M + a[i]
B = N + n_[j]
abmn.append(np.vstack([A, B, M, N]).T) abmn.append(np.vstack([A, B, M, N]).T)
abmn = np.vstack(abmn)
abmn = abmn[(abmn <= nelec).all(1), :]
df = pd.DataFrame(abmn, columns=['a', 'b', 'm', 'n'])
df = df.sort_values(['a', 'b', 'm', 'n']).reset_index(drop=True)
return df
proto_mtrx = pd.DataFrame(np.vstack(abmn)) def multigrad(nelec, a, n, s):
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num] ''' Generate measurement matrix for multigradient array.
return proto_mtrx
def multigrad(elec_num, a, n, s):
''' Genetrate measurement matrix for multigradient array Torleif Dahlin.
Parameters Parameters
---------- ----------
elec_num : int nelec : int
Number of electrodes Number of electrodes
a : int a : int or list of int
Spacing between potential electrodes (in electrode spacing). Spacing between potential electrodes (in electrode spacing).
n : int n : int or list of int
Multiplier for `a` to determine spacing from A to M. Multiplier for `a` to determine spacing from A to M.
s : int s : int or list of int
Seperation factor for current electrodes, should be the intermediate Separation factor for current electrodes, should be the intermediate
numbers. numbers.
''' '''
elec_id = np.arange(elec_num)+1 # check parameters
if isinstance(a, int):
a = np.arange(a) + 1
if isinstance(n, int):
n = np.arange(n) + 1
if isinstance(s, int):
n = np.arange(s) + 1
elec = np.arange(nelec) + 1
abmn = [] abmn = []
for aa in a:
if isinstance(a, list) is False: for nn in n:
n_ = np.array(range(0, n))+1 for ss in s:
s_ = np.array(range(0, s))+1 A = elec
for j in np.array(range(0, len(n_))): B = A + ss + 2
for k in np.array(range(0, len(s_))): M = A + nn * aa
A = elec_id N = M + aa
B = A + s_[k] + 2
M = A + n_[j] * a
N = M + a
abmn.append(np.vstack([A, B, M, N]).T) abmn.append(np.vstack([A, B, M, N]).T)
abmn = np.vstack(abmn)
else: abmn = abmn[(abmn <= nelec).all(1), :]
for i in np.array(range(0,len(a))): df = pd.DataFrame(abmn, columns=['a', 'b', 'm', 'n'])
n_ = np.array(range(0, n[i]))+1 df = df.sort_values(['a', 'b', 'm', 'n']).reset_index(drop=True)
s_ = np.array(range(0, s[i]))+1 return df
for j in np.array(range(0, len(n_))):
for k in np.array(range(0, len(s_))):
A = elec_id
B = A + s_[k] + 2
M = A + n_[j] * a[i]
N = M + a[i]
abmn.append(np.vstack([A, B, M, N]).T)
proto_mtrx = pd.DataFrame(np.vstack(abmn))
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] > proto_mtrx.iloc[:,3]]
proto_mtrx = proto_mtrx[proto_mtrx.iloc[:,1] <= elec_num]
return proto_mtrx
# test code # test code
# x1 = dpdp1(24, 2, 8) # x1 = dpdp1(24, 2, 8)
# x2 = dpdp2(24, 2, 8)
# x3 = wenner_alpha(24, 1) # x3 = wenner_alpha(24, 1)
# x4 = wenner_beta(24, 1)
# x5 = wenner_gamma(24, 1)
# x6 = schlum1(24, 1, 10) # x6 = schlum1(24, 1, 10)
# x7 = schlum2(24, 1, 10)
# x8 = multigrad(24, 1, 10, 2) # x8 = multigrad(24, 1, 10, 2)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment