Skip to content

Commit

Permalink
add skip plot option to cross_section /plot_cross_section,
Browse files Browse the repository at this point in the history
  • Loading branch information
meihuisu committed Oct 29, 2024
1 parent a5537c4 commit a2e234e
Show file tree
Hide file tree
Showing 8 changed files with 189 additions and 155 deletions.
Binary file modified dist/ucvm_plotting-0.0.6-py3-none-any.whl
View file
Binary file not shown.
Binary file modified dist/ucvm_plotting-0.0.6.tar.gz
View file
Binary file not shown.
36 changes: 6 additions & 30 deletions pycvm/basin_slice.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -106,8 +106,8 @@ def getplotvals(self, mproperty="vs"):
i = i + 1

##
# Plots the basin depth data as a horizontal slice. This code is very similar to the
# HorizontalSlice routine.
# Plots the basin depth data as a horizontal slice.
# Output the basin depth data as a horizontal slice.
#
# @param horizontal_label The horizontal label of the plot. Optional.
def plot(self, horizontal_label = "Depth (km)") :
Expand All @@ -130,33 +130,7 @@ def plot(self, horizontal_label = "Depth (km)") :
self.meta['mproperty']="vs"

HorizontalSlice.plot(self, horizontal_label)

##
# Output the basin depth data as a horizontal slice. This code is very similar to the
# HorizontalSlice routine.
#
# @param horizontal_label The horizontal label of the plot. Optional.
def plot_skip(self, horizontal_label = "Depth (km)") :

if self.upperleftpoint.description == None:
location_text = ""
else:
location_text = self.upperleftpoint.description + " "

# Gets the better CVM description if it exists.
try:
cvmdesc = UCVM_CVMS[self.cvm]
except:
cvmdesc = self.cvm

if 'title' not in self.meta:
self.title = "%sBasin Depth Map For %s" % (location_text, cvmdesc)
self.meta['title'] = self.title;

self.meta['mproperty']="vs"

HorizontalSlice.plot_skip(self, horizontal_label)


##
# @class Z10Slice
# @brief Gets a Z1.0 data map.
Expand Down Expand Up @@ -184,6 +158,7 @@ def getplotvals(self, mproperty):

##
# Plots the Z1.0 slice.
# Output the basin depth data as a horizontal slice.
#
def plot(self) :

Expand All @@ -203,8 +178,9 @@ def plot(self) :
else:
self.title = "%sZ1.0 Map For %s" % (location_text, cvmdesc)
self.meta['title'] = self.title

BasinSlice.plot(self)


##
# @class Z25Slice
Expand Down
182 changes: 128 additions & 54 deletions pycvm/cross_section.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -96,6 +96,11 @@ def __init__(self, startingpoint, endingpoint, meta={}):
self.cvm = self.meta['cvm']
else:
self.cvm = None

if 'data_type' in self.meta :
self.mproperty = self.meta['data_type']
else:
self.mproperty = "vs"

if 'datafile' in self.meta :
self.datafile = self.meta['datafile']
Expand All @@ -107,6 +112,42 @@ def __init__(self, startingpoint, endingpoint, meta={}):
else:
self.filename = None

if 'installdir' in self.meta :
self.installdir = self.meta['installdir']
else:
self.installdir = None

if 'configfile' in self.meta :
self.configfile = self.meta['configfile']
else:
self.configfile = None

# Gets the better CVM description if it exists.
try:
cvmdesc = UCVM_CVMS[self.cvm]
except:
cvmdesc = self.cvm

if self.startingpoint.description == None:
location_text = ""
else:
location_text = self.startingpoint.description + " "

if 'title' in self.meta :
title = self.meta['title']
else:
title = "%s%s Cross Section from (%.2f, %.2f) to (%.2f, %.2f)" % (location_text, cvmdesc, self.startingpoint.longitude, \
self.startingpoint.latitude, self.endingpoint.longitude, self.endingpoint.latitude)
self.meta['title']=title

if 'skip' in self.meta:
self.skip= True;
else
self.skip = None

self.ucvm = UCVM(install_dir=self.installdir, config_file=self.configfile, z_range=self.z_range, floors=self.floors)



##
# Generates the depth profile in a format that is ready to plot.
Expand All @@ -125,7 +166,6 @@ def getplotvals(self, mproperty='vs'):
num_prof = int(math.sqrt((x2-x1)*(x2-x1) + \
(y2-y1)*(y2-y1))/self.hspacing)

# cnt=0
jstart = self.startingdepth
for j in range(int(self.startingdepth), int(self.todepth) + 1, int(self.vspacing)):
depth_list.append( round(j,3))
Expand All @@ -137,9 +177,6 @@ def getplotvals(self, mproperty='vs'):
if ( j == jstart) :
lon_list.append( round(lon,5))
lat_list.append( round(lat,5))
# if(cnt < 10) :
# print("point.. lon "+str(lon)+" lat "+str(lat)+" j "+str(j))
# cnt += 1

self.lon_list=lon_list
self.lat_list=lat_list
Expand All @@ -161,12 +198,12 @@ def getplotvals(self, mproperty='vs'):


if self.datafile.rfind(".binary") != -1 :
data = u.import_binary(self.datafile, self.num_x, self.num_y)
data = ucvm.import_binary(self.datafile, self.num_x, self.num_y)
else :
if self.datafile.rfind(".raw") != -1 :
data = u.import_raw_data(self.datafile, self.num_x, self.num_y)
data = ucvm.import_raw_data(self.datafile, self.num_x, self.num_y)
else: ## with .bin file
data = u.import_np_float_array(self.datafile, self.num_x, self.num_y)
data = ucvm.import_np_float_array(self.datafile, self.num_x, self.num_y)

## this set of data is only for --datatype: either 'vs', 'vp', 'rho', or 'poisson'
## The 2D array of retrieved material properties.
Expand All @@ -187,7 +224,7 @@ def getplotvals(self, mproperty='vs'):

print("\nUsing --> "+self.datafile)
else:
data = u.query(point_list, self.cvm)
data = ucvm.query(point_list, self.cvm)


## Private number of x points.
Expand All @@ -204,18 +241,21 @@ def getplotvals(self, mproperty='vs'):
for x in range(0, self.num_x):
self.materialproperties[y][x] = data[y * self.num_x + x]
##
# Plots the horizontal slice either to an image or a file name.
# Plots the cross section slice either to an image or a file name.
#
def plot(self) :
def plot(self):

self.installdir = None
if 'installdir' in self.meta :
self.installdir = self.meta['installdir']
self.skip :
this._file(self)
else:
this._plot_file(self)

self.configfile = None
if 'configfile' in self.meta :
self.configfile = self.meta['configfile']
##
# Plots the cross section slice to an image and save a data file.
#
def _plot_file(self) :

color_scale = None
if 'color' in self.meta :
color_scale = self.meta['color']

Expand All @@ -226,30 +266,7 @@ def plot(self) :
if color_scale == "b" and scale_gate is None:
scale_gate=2.5

if self.startingpoint.description == None:
location_text = ""
else:
location_text = self.startingpoint.description + " "

if 'data_type' in self.meta :
mproperty = self.meta['data_type']
else:
mproperty = "vs"

# Gets the better CVM description if it exists.
try:
cvmdesc = UCVM_CVMS[self.cvm]
except:
cvmdesc = self.cvm

if 'title' in self.meta :
title = self.meta['title']
else:
title = "%s%s Cross Section from (%.2f, %.2f) to (%.2f, %.2f)" % (location_text, cvmdesc, self.startingpoint.longitude, \
self.startingpoint.latitude, self.endingpoint.longitude, self.endingpoint.latitude)
self.meta['title']=title

self.getplotvals(mproperty)
self.getplotvals(self.mproperty)

# Call the plot object.
p = Plot(None, None, None, None, 10, 10)
Expand Down Expand Up @@ -313,7 +330,7 @@ def plot(self) :
elif mproperty != "poisson" :
datapoints[y][x] = self.materialproperties[y][x].getProperty(mproperty)
else:
datapoints[y][x] = u.poisson(self.materialproperties[y][x].getProperty("vs"), self.materialproperties[y][x].getProperty("vp"))
datapoints[y][x] = ucvm.poisson(self.materialproperties[y][x].getProperty("vs"), self.materialproperties[y][x].getProperty("vp"))


u = UCVM(install_dir=self.installdir, config_file=self.configfile)
Expand All @@ -336,15 +353,15 @@ def plot(self) :
colormap = basemap.cm.GMT_seis

if self.scalemin != None and self.scalemax != None:
BOUNDS= u.makebounds(float(self.scalemin), float(self.scalemax), 5)
TICKS = u.maketicks(float(self.scalemin), float(self.scalemax), 5)
BOUNDS= ucvm.makebounds(float(self.scalemin), float(self.scalemax), 5)
TICKS = ucvm.maketicks(float(self.scalemin), float(self.scalemax), 5)
umax=round(self.scalemax)
umin=round(self.scalemin)
umean=round((umax+umin)/2)
else:
## default BOUNDS are from 0 to 5
BOUNDS = u.makebounds()
TICKS = u.maketicks()
BOUNDS = ucvm.makebounds()
TICKS = ucvm.maketicks()
umax=round(self.max_val)
umin=round(self.min_val)
umean=round(self.mean_val)
Expand All @@ -366,8 +383,8 @@ def plot(self) :
norm = mcolors.Normalize(vmin=BOUNDS[0],vmax=BOUNDS[len(BOUNDS) - 1])
elif color_scale == "sd":
colormap = basemap.cm.GMT_seis
BOUNDS= u.makebounds(umin, umax, 5, umean, substep=5)
TICKS = u.maketicks(umin, umax, 5)
BOUNDS= ucvm.makebounds(umin, umax, 5, umean, substep=5)
TICKS = ucvm.maketicks(umin, umax, 5)
norm = mcolors.Normalize(vmin=BOUNDS[0],vmax=BOUNDS[len(BOUNDS) - 1])
elif color_scale == "b":
C = []
Expand All @@ -385,8 +402,8 @@ def plot(self) :
colormap = plot_cmapDiscretize(basemap.cm.GMT_seis_r, len(BOUNDS) - 1)
norm = mcolors.BoundaryNorm(BOUNDS, colormap.N)
elif color_scale == 'dd':
BOUNDS= u.makebounds(umin, umax, 5, umean, substep=5)
TICKS = u.maketicks(umin, umax, 5)
BOUNDS= ucvm.makebounds(umin, umax, 5, umean, substep=5)
TICKS = ucvm.maketicks(umin, umax, 5)
colormap = plot_cmapDiscretize(basemap.cm.GMT_seis, len(BOUNDS) - 1)
norm = mcolors.BoundaryNorm(BOUNDS, colormap.N)
else:
Expand All @@ -397,6 +414,9 @@ def plot(self) :
colormap = plot_cmapDiscretize(bwr, len(BOUNDS) - 1)


ucvm = self.ucvm


## MEI, TODO this is a temporary way to generate an output of a cross_section input file
if( self.datafile == None ):
self.meta['num_x'] = self.num_x
Expand All @@ -409,14 +429,14 @@ def plot(self) :
self.meta['lat_list'] = self.lat_list
self.meta['depth_list'] = self.depth_list
if self.filename:
u.export_metadata(self.meta,self.filename)
u.export_np_float_array(datapoints,self.filename)
ucvm.export_metadata(self.meta,self.filename)
ucvm.export_np_float_array(datapoints,self.filename)
else:
#https://stackoverflow.com/questions/2257441/random-string-generation-with-upper-case-letters-and-digits-in-python
rnd=''.join(random.SystemRandom().choice(string.ascii_uppercase + string.digits) for _ in range(6))
f = "cross_section"+rnd
u.export_metadata(self.meta,f)
u.export_np_float_array(datapoints,f)
ucvm.export_metadata(self.meta,f)
ucvm.export_np_float_array(datapoints,f)


img = plt.imshow(newdatapoints, cmap=colormap, norm=norm)
Expand All @@ -427,7 +447,7 @@ def plot(self) :
"%.2f" % (self.startingdepth+ ((self.todepth-self.startingdepth)/2)/1000), \
"%.2f" % (self.todepth / 1000)])

plt.title(title)
plt.title(this.title)

cax = plt.axes([0.1, 0.1, 0.8, 0.02])
cbar = plt.colorbar(img, cax=cax, orientation='horizontal',ticks=TICKS,spacing='proportional')
Expand All @@ -446,3 +466,57 @@ def plot(self) :
plt.savefig(self.filename)
else:
plt.show()

##
# Create the cross section slice data file only
#
def _file(self) :

self.getplotvals(self.mproperty)

datapoints = np.arange(self.num_x * self.num_y,dtype=np.float32).reshape(self.num_y, self.num_x)

nancnt=0
zerocnt=0
negcnt=0
for i in range(0, self.num_y):
for j in range(0, self.num_x):
if (self.datafile != None) :
datapoints[i][j] = self.materialproperties[i][j].getProperty(self.mproperty)
elif self.mproperty != "poisson":
datapoints[i][j] = self.materialproperties[i][j].getProperty(self.mproperty)
if (datapoints[i][j] == 0) :
zerocnt=zerocnt+1
if (datapoints[i][j] < 0) :
negcnt=negcnt+1
if(datapoints[i][j] == -1 ) :
datapoints[i][j]=np.nan
nancnt=nancnt+1
else :
datapoints[i][j] = ucvm.poisson(self.materialproperties[i][j].vs, self.materialproperties[i][j].vp)

self.max_val= np.nanmax(datapoints)
self.min_val=np.nanmin(datapoints)
self.mean_val=np.mean(datapoints)

ucvm = self.ucvm

if( self.datafile == None ):
self.meta['num_x'] = self.num_x
self.meta['num_y'] = self.num_y
self.meta['datapoints'] = datapoints.size
self.meta['max'] = self.max_val.item()
self.meta['min'] = self.min_val.item()
self.meta['mean'] = self.mean_val.item()
self.meta['lon_list'] = self.lon_list
self.meta['lat_list'] = self.lat_list
self.meta['depth_list'] = self.depth_list
if self.filename:
ucvm.export_metadata(self.meta,self.filename)
ucvm.export_np_float_array(datapoints,self.filename)
else:
#https://stackoverflow.com/questions/2257441/random-string-generation-with-upper-case-letters-and-digits-in-python
rnd=''.join(random.SystemRandom().choice(string.ascii_uppercase + string.digits) for _ in range(6))
f = "cross_section"+rnd
ucvm.export_metadata(self.meta,f)
ucvm.export_np_float_array(datapoints,f)
1 change: 1 addition & 0 deletions pycvm/cvm_ucvm.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
import math
import struct
import json
import pdb

# Numpy is required.
try:
Expand Down
Loading

0 comments on commit a2e234e

Please sign in to comment.