如何从巨大的（Mbs 级）.txt 文件中制作更平滑的二维地图？
How to make smoother 2d maps out of a huge (Mbs-scale) .txt file?
我有一个 .txt 文件，其中包含五列数据，其中第 1 列和第 2 列对应于沿维度 1 和 2 的投影距离（像素比例为 0.50），然后第 3、4 和 5 列是投影 0 "z"、1 "x" 和 2 "y" 中任意（否则未知）函数的值。 （原始文件长约 1e6 行，但这是我设法上传的一个小版本：m11q.txt）我想将此 .txt table 转换为实际的二维地图，其中每个像素的颜色强度将对应于该像素中函数的大小。这意味着我将拥有三个平面 xy、yz 和 xz，其中每个像素都填充了一种颜色。
我对第 3、4 和 5 列的值大部分为零，其中很少有 $10^{14}$ 到 $10^{21}$ 的顺序。我有兴趣拥有三个轮廓级别。假设如果函数值小于 $10^{14}$，颜色应该是黑色，如果函数值在 $10^{14}$ 和 $10^{15.5}$ 之间，颜色应该是黄色，如果如果函数值在 $10^{15.5}$ 和 $10^{17}$ 之间，颜色应该是橙色，如果函数值在 $10^{17}$ 和 $10^{20.3}$ 之间，颜色应该是红色最后，如果超过 $10^{20.3}$，颜色应该是蓝色。换句话说，我的二维地图实际上是一组由 .txt table 生成的 5 个等高线图。
这是我的原始代码版本：
import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
nInterp = 4000
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
LAF_covering_fraction_XY = np.count_nonzero((np.power(10, 14) <= z0) & (z0 < np.power(10, 17))) * 1.0 / z0.size
LAF_covering_fraction_YZ = np.count_nonzero((np.power(10, 14) <= z1) & (z1 < np.power(10, 17))) * 1.0 / z1.size
LAF_covering_fraction_XZ = np.count_nonzero((np.power(10, 14) <= z2) & (z2 < np.power(10, 17))) * 1.0 / z2.size
LAF = np.sqrt(np.power(LAF_covering_fraction_XY, 2)+np.power(LAF_covering_fraction_YZ, 2)+np.power(LAF_covering_fraction_XZ, 2))/np.sqrt(3)
LAF_sigma = max([abs(LAF-LAF_covering_fraction_XY),abs(LAF-LAF_covering_fraction_YZ),abs(LAF-LAF_covering_fraction_XZ)])
LLS_covering_fraction_XY = np.count_nonzero((np.power(10, 17) <= z0) & (z0 < np.power(10, 20.3))) * 1.0 / z0.size
LLS_covering_fraction_YZ = np.count_nonzero((np.power(10, 17) <= z1) & (z1 < np.power(10, 20.3))) * 1.0 / z1.size
LLS_covering_fraction_XZ = np.count_nonzero((np.power(10, 17) <= z2) & (z2 < np.power(10, 20.3))) * 1.0 / z2.size
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
DLA_covering_fraction_XY = np.count_nonzero(np.power(10, 20.3) <= z0) * 1.0 / z0.size
DLA_covering_fraction_YZ = np.count_nonzero(np.power(10, 20.3) <= z1) * 1.0 / z1.size
DLA_covering_fraction_XZ = np.count_nonzero(np.power(10, 20.3) <= z2) * 1.0 / z2.size
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
z0plot = plt.subplot(221)
plt.imshow(z0, vmin=0.0001+covering_fraction_z.min(), vmax=covering_fraction_z.max(), origin='lower',
extent=[x.min(), x.max(), y.min(), y.max()],
norm=colors.LogNorm())
plt.ylabel("Y [kpc]")
z0plot.set(aspect=1)
axparameters = plt.subplot(222)
frame = plt.gca()
frame.axes.get_xaxis().set_ticks([])
frame.axes.get_yaxis().set_ticks([])
axparameters.set(aspect=1)
plt.text(0.02, 0.2, 'LAF=%.5f$\pm$%.5f \n\n\n LLS=%.5f$\pm$%.5f \n\n\n DLA=%.5f$\pm$%.5f' % (LAF, LAF_sigma, LLS, LLS_sigma, DLA, DLA_sigma))
plt.colorbar()
plt.tight_layout()
z1plot = plt.subplot(223)
plt.imshow(z1, vmin=0.0001+covering_fraction_y.min(), vmax=covering_fraction_y.max(), origin='lower',
extent=[x.min(), x.max(), y.min(), y.max()],
norm=colors.LogNorm())
plt.xlabel("X [kpc]")
plt.ylabel("Z [kpc]")
z1plot.set(aspect=1)
z2plot = plt.subplot(224)
plt.imshow(z2, vmin=0.0001+covering_fraction_x.min(), vmax=covering_fraction_x.max(), origin='lower',
extent=[x.min(), x.max(), y.min(), y.max()],
norm=colors.LogNorm())
plt.xlabel("Y [kpc]")
z2plot.set(aspect=1)
myfig = fig.Figure(figsize=(8, 8), dpi=900, facecolor='w', edgecolor='k', tight_layout=True)
plt.suptitle('HI Column Density Maps', size=20)
plt.tight_layout(rect=[0, 0, 1, 0.9])
plt.savefig('m11q.png', dpi=900)
plt.show()
生成以下地图
这是 Xevaquor 编写的代码，作为对我原来的回答 post:
# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
#import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
import matplotlib.ticker as ticker
from matplotlib import gridspec
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
def fmt(x, pos):
return '10e{}'.format(np.log10(x))
nInterp = 4300
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
LAF_covering_fraction_XY = np.count_nonzero((np.power(10, 14) <= z0) & (z0 < np.power(10, 17))) * 1.0 / z0.size
LAF_covering_fraction_YZ = np.count_nonzero((np.power(10, 14) <= z1) & (z1 < np.power(10, 17))) * 1.0 / z1.size
LAF_covering_fraction_XZ = np.count_nonzero((np.power(10, 14) <= z2) & (z2 < np.power(10, 17))) * 1.0 / z2.size
LAF = np.sqrt(np.power(LAF_covering_fraction_XY, 2)+np.power(LAF_covering_fraction_YZ, 2)+np.power(LAF_covering_fraction_XZ, 2))/np.sqrt(3)
LAF_sigma = max([abs(LAF-LAF_covering_fraction_XY),abs(LAF-LAF_covering_fraction_YZ),abs(LAF-LAF_covering_fraction_XZ)])
LLS_covering_fraction_XY = np.count_nonzero((np.power(10, 17) <= z0) & (z0 < np.power(10, 20.3))) * 1.0 / z0.size
LLS_covering_fraction_YZ = np.count_nonzero((np.power(10, 17) <= z1) & (z1 < np.power(10, 20.3))) * 1.0 / z1.size
LLS_covering_fraction_XZ = np.count_nonzero((np.power(10, 17) <= z2) & (z2 < np.power(10, 20.3))) * 1.0 / z2.size
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
DLA_covering_fraction_XY = np.count_nonzero(np.power(10, 20.3) <= z0) * 1.0 / z0.size
DLA_covering_fraction_YZ = np.count_nonzero(np.power(10, 20.3) <= z1) * 1.0 / z1.size
DLA_covering_fraction_XZ = np.count_nonzero(np.power(10, 20.3) <= z2) * 1.0 / z2.size
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
fig = plt.figure(figsize=(8, 8))
gs = gridspec.GridSpec(2, 3, width_ratios=[5, 5, 1], height_ratios=[1,1,1])
ax0 = plt.subplot(gs[0])
ax1 = plt.subplot(gs[3])
ax2 = plt.subplot(gs[4])
#cmap = colors.ListedColormap(['#ffdddd','#000000', '#ff9999', '#ff7777', '#ff5555', '#ff3333' ,'#ff1111'])
#'gray','green','pink','brown','blue','violet','yellow','orange',
cmap = colors.ListedColormap(['#000000', 'gray','#ffdddd', '#ff9999', '#ff7777', '#ff5555', '#ff3333', '#ff1111', 'blue'])
boundaries = [np.power(10., 10), np.power(10., 11), np.power(10., 12), np.power(10., 13), np.power(10., 14), np.power(10., 15), np.power(10., 16), np.power(10.,17), np.power(10.,20.3), np.power(10.,25)]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
pcmz0 = ax0.pcolormesh(xi, yi, z0, cmap=cmap, norm=norm)
ax0.set_xticks([x.min(), x.max()])
ax0.set_yticks([y.min(), y.max()])
ax0.set_xlim([x.min(), x.max()])
ax0.set_ylim([y.min(), y.max()])
ax0.set(aspect=1)
pcmz1 = ax1.pcolormesh(xi, yi, z1, cmap=cmap, norm=norm)
ax1.set_xticks([x.min(), x.max()])
ax1.set_yticks([y.min(), y.max()])
ax1.set_xlim([x.min(), x.max()])
ax1.set_ylim([y.min(), y.max()])
ax1.set(aspect=1)
pcmz2 = ax2.pcolormesh(xi, yi, z2, cmap=cmap, norm=norm)
ax2.set_xticks([x.min(), x.max()])
ax2.set_yticks([y.min(), y.max()])
ax2.set_xlim([x.min(), x.max()])
ax2.set_ylim([y.min(), y.max()])
ax2.set(aspect=1)
axparameters = plt.subplot(gs[1])
frame = plt.gca()
frame.axes.get_xaxis().set_ticks([])
frame.axes.get_yaxis().set_ticks([])
axparameters.set(aspect=1)
plt.text(0.1, 0.2, 'LAF=%.5f$\pm$%.5f \n\n\n LLS=%.5f$\pm$%.5f \n\n\n DLA=%.5f$\pm$%.5f' % (LAF, LAF_sigma, LLS, LLS_sigma, DLA, DLA_sigma), size=10)
axc = plt.subplot(gs[:,2])
plt.colorbar(pcmz1, cax=axc,format=ticker.FuncFormatter(fmt))
ax0.set_ylabel("Y [kpc]")
ax1.set_xlabel("X [kpc]")
ax1.set_ylabel("Z [kpc]")
ax2.set_xlabel("Y [kpc]")
fig.suptitle('HI Column Density Maps', size=20)
gs.tight_layout(fig, rect=[0, 0, 1, 0.9])
plt.savefig('m11q.png', dpi=900)
plt.show()
生成以下地图：
这里是代码的最终版本，除了在 post:
的答案（由我）中提到的一些警告外，它完成了我的所有目标
import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
#import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
import matplotlib.ticker as ticker
from matplotlib import gridspec
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
def fmt(x, pos):
return '10e{}'.format(np.log10(x))
nInterp = 3000
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
Pixel_Smoothing_Scale = 0.497
Halo_Pixels_No = (x.max()-x.min())*(y.max()-y.min())/ np.power(Pixel_Smoothing_Scale, 2)
DLA_Pixels_No_xy, DLA_Pixels_No_yz, DLA_Pixels_No_xz = 0, 0, 0
LLS_Pixels_No_xy, LLS_Pixels_No_yz, LLS_Pixels_No_xz = 0, 0, 0
LAF_Pixels_No_xy, LAF_Pixels_No_yz, LAF_Pixels_No_xz = 0, 0, 0
for i in np.arange(0, np.size(x)):
if covering_fraction_z[i] >= np.power(10., 20.3):
DLA_Pixels_No_xy+=1
if ((covering_fraction_z[i] >= np.power(10., 17)) & (covering_fraction_z[i] < np.power(10., 20.3))):
LLS_Pixels_No_xy+=1
if ((covering_fraction_z[i] >= np.power(10., 14)) & (covering_fraction_z[i] < np.power(10., 17))):
LAF_Pixels_No_xy+=1
if covering_fraction_x[i] >= np.power(10., 20.3):
DLA_Pixels_No_yz+=1
if ((covering_fraction_x[i] >= np.power(10., 17)) & (covering_fraction_x[i] < np.power(10., 20.3))):
LLS_Pixels_No_yz+=1
if ((covering_fraction_x[i] >= np.power(10., 14)) & (covering_fraction_x[i] < np.power(10., 17))):
LAF_Pixels_No_yz+=1
if covering_fraction_y[i] >= np.power(10., 20.3):
DLA_Pixels_No_xz+=1
if ((covering_fraction_y[i] >= np.power(10., 17)) & (covering_fraction_y[i] < np.power(10., 20.3))):
LLS_Pixels_No_xz+=1
if ((covering_fraction_y[i] >= np.power(10., 14)) & (covering_fraction_y[i] < np.power(10., 17))):
LAF_Pixels_No_xz+=1
DLA_covering_fraction_XY, DLA_covering_fraction_YZ, DLA_covering_fraction_XZ = DLA_Pixels_No_xy/Halo_Pixels_No, DLA_Pixels_No_yz/Halo_Pixels_No, DLA_Pixels_No_xz/Halo_Pixels_No
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
LLS_covering_fraction_XY, LLS_covering_fraction_YZ, LLS_covering_fraction_XZ = LLS_Pixels_No_xy/Halo_Pixels_No, LLS_Pixels_No_yz/Halo_Pixels_No, LLS_Pixels_No_xz/Halo_Pixels_No
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
LAF_covering_fraction_XY, LAF_covering_fraction_YZ, LAF_covering_fraction_XZ = LAF_Pixels_No_xy/Halo_Pixels_No, LAF_Pixels_No_yz/Halo_Pixels_No, LAF_Pixels_No_xz/Halo_Pixels_No
LAF = np.sqrt(np.power(LAF_covering_fraction_XY, 2)+np.power(LAF_covering_fraction_YZ, 2)+np.power(LAF_covering_fraction_XZ, 2))/np.sqrt(3)
LAF_sigma = max([abs(LAF-LAF_covering_fraction_XY),abs(LAF-LAF_covering_fraction_YZ),abs(LAF-LAF_covering_fraction_XZ)])
fig = plt.figure(figsize=(8, 8))
gs = gridspec.GridSpec(2, 3, width_ratios=[5, 5, 1], height_ratios=[1,1,1])
ax0 = plt.subplot(gs[0])
ax1 = plt.subplot(gs[3])
ax2 = plt.subplot(gs[4])
#cmap = colors.ListedColormap(['#ffdddd','#000000', '#ff9999', '#ff7777', '#ff5555', '#ff3333' ,'#ff1111'])
#'gray','green','pink','brown','blue','violet','yellow','orange',
cmap = colors.ListedColormap(['#ffdddd', '#ff7777', '#ff3333', 'yellow', 'orange','green', 'red', 'blue'])
boundaries = [np.power(10., 4), np.power(10., 6), np.power(10., 8), np.power(10., 10), np.power(10., 12), np.power(10., 14), np.power(10.,17), np.power(10.,20.3), np.power(10.,25)]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
pcmz0 = ax0.pcolormesh(xi, yi, z0, cmap=cmap, norm=norm)
ax0.set_xticks([x.min(), x.max()])
ax0.set_yticks([y.min(), y.max()])
ax0.set_xlim([x.min(), x.max()])
ax0.set_ylim([y.min(), y.max()])
ax0.set(aspect=1)
pcmz1 = ax1.pcolormesh(xi, yi, z1, cmap=cmap, norm=norm)
ax1.set_xticks([x.min(), x.max()])
ax1.set_yticks([y.min(), y.max()])
ax1.set_xlim([x.min(), x.max()])
ax1.set_ylim([y.min(), y.max()])
ax1.set(aspect=1)
pcmz2 = ax2.pcolormesh(xi, yi, z2, cmap=cmap, norm=norm)
ax2.set_xticks([x.min(), x.max()])
ax2.set_yticks([y.min(), y.max()])
ax2.set_xlim([x.min(), x.max()])
ax2.set_ylim([y.min(), y.max()])
ax2.set(aspect=1)
axparameters = plt.subplot(gs[1])
frame = plt.gca()
frame.axes.get_xaxis().set_ticks([])
frame.axes.get_yaxis().set_ticks([])
axparameters.set(aspect=1)
plt.text(0.1, 0.2, 'LAF=%.6f$\pm$%.6f \n\n\n LLS=%.6f$\pm$%.6f \n\n\n DLA=%.6f$\pm$%.6f' % (LAF, LAF_sigma, LLS, LLS_sigma, DLA, DLA_sigma), size=9)
axc = plt.subplot(gs[:,2])
plt.colorbar(pcmz1, cax=axc,format=ticker.FuncFormatter(fmt))
ax0.set_ylabel("Y [kpc]")
ax1.set_xlabel("X [kpc]")
ax1.set_ylabel("Z [kpc]")
ax2.set_xlabel("Y [kpc]")
fig.suptitle('HI Column Density Maps', size=20)
gs.tight_layout(fig, rect=[0, 0, 1, 0.9])
plt.savefig('m11q.png', dpi=900)
plt.show()
生成以下地图：
(1) 和(2) 是相关的。您可以定义自定义配色方案。对于任意值，它可以是例如：
cmap = colors.ListedColormap(['red', '#000000','#444444', '#666666', '#ffffff', 'blue', 'orange'])
boundaries = [-1, -0.9, -0.6, -0.3, 0, 0.3, 0.6, 1]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
这将按照以下方式映射值：
- -1.0 x < -0.9 到红色
- -0.9 < x < -0.6 至黑色 (#000000)
- -0.6 < x < -0.3 至 #444444
- 等等
颜色 i 将用于边界 i 和 i+1 之间的值。颜色可以通过名称 ('red'、'green')、HTML 代码 ('#ffaa44'、'#441188') 或 RGB 元组 ((0.2, 0.9, 0.45)) 指定。
完整示例：
# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
#import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
import matplotlib.ticker as ticker
from matplotlib import gridspec
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
def fmt(x, pos):
return '10e{}'.format(np.log10(x))
nInterp = 402
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
fig = plt.figure(figsize=(6, 6))
gs = gridspec.GridSpec(2, 3, width_ratios=[4, 4, 1], height_ratios=[1,1,1])
ax0 = plt.subplot(gs[0])
ax1 = plt.subplot(gs[3])
ax2 = plt.subplot(gs[4])
cmap = colors.ListedColormap(['#ffdddd','#000000', '#ff9999', '#ff7777', '#ff5555', '#ff3333' ,'#ff1111'])
boundaries = [np.power(10., 14), np.power(10., 15.5), np.power(10.,17), np.power(10.,20.3), np.power(10.,25)]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
pcmz0 = ax0.pcolormesh(xi, yi, z0, cmap=cmap, norm=norm)
ax0.set_xticks([x.min(), x.max()]) #arbitrary, also for Y axis ticks
ax0.set_xlim([x.min(), x.max()])
ax0.set_ylim([y.min(), y.max()])
pcmz1 = ax1.pcolormesh(xi, yi, z1, cmap=cmap, norm=norm)
ax1.set_xticks([x.min(), x.max()])
ax1.set_xlim([x.min(), x.max()])
ax1.set_ylim([y.min(), y.max()])
pcmz2 = ax2.pcolormesh(xi, yi, z2, cmap=cmap, norm=norm)
ax2.set_xticks([x.min(), x.max()])
ax2.set_xlim([x.min(), x.max()])
ax2.set_ylim([y.min(), y.max()])
axc = plt.subplot(gs[:,2])
plt.colorbar(pcmz1, cax=axc,format=ticker.FuncFormatter(fmt))
ax0.set_ylabel("Y [kpc]")
ax1.set_xlabel("X [kpc]")
ax1.set_ylabel("Z [kpc]")
ax2.set_xlabel("Y [kpc]")
plt.tight_layout()
plt.savefig('m11q.png', dpi=900)
plt.show()
对于 (3)。有很多方法可以做到这一点。例如：
fraction = np.count_nonzero((np.power(10, 14) <= Z) & (Z < np.power(10, 15.5)) * 1.0 / Z.size
将为您提供 [10e14, 10e15.5)
范围内的部分值
(4) 取决于创建子图的方法。对于 GridSpec，如果您不明确声明它，它将不会出现。
这是问题更新版本的答案。最初我无法为我的地图制作颜色条。 Xevaquor 通过上述回答帮助我解决了这个问题。所以，答案是由 Xevaquor 合法发布的，现在被标记为 "the answer"。然而，当我花时间解决一些覆盖分数的问题时，我终于能够从数学上理解这个问题，然后尝试在 python 代码中用数字来解决它。我的问题是我没有意识到我所追求的实际上是像素数量的比率。这就是为什么我必须计算位于列密度范围内的像素数，然后将其除以进行地图每个投影的像素总数。由于已知平滑比例为 $0.497 kpc$，我们知道像素总数将是地图的物理面积 $kpc^{2}$ 除以每个像素的面积即 $0.497^{2}$。因此，可以将以下代码片段添加到原始 and/or Xevaquor 的代码中以正确说明覆盖分数：
Pixel_Smoothing_Scale = 0.497
Halo_Pixels_No = (x.max()-x.min())*(y.max()-y.min())/ np.power(Pixel_Smoothing_Scale, 2)
DLA_Pixels_No_xy, DLA_Pixels_No_yz, DLA_Pixels_No_xz = 0, 0, 0
LLS_Pixels_No_xy, LLS_Pixels_No_yz, LLS_Pixels_No_xz = 0, 0, 0
for i in np.arange(0, np.size(x)):
if covering_fraction_z[i] >= np.power(10., 20.3):
DLA_Pixels_No_xy+=1
if ((covering_fraction_z[i] >= np.power(10., 17)) & (covering_fraction_z[i] < np.power(10., 20.3))):
LLS_Pixels_No_xy+=1
if covering_fraction_x[i] >= np.power(10., 20.3):
DLA_Pixels_No_yz+=1
if ((covering_fraction_x[i] >= np.power(10., 17)) & (covering_fraction_x[i] < np.power(10., 20.3))):
LLS_Pixels_No_yz+=1
if covering_fraction_y[i] >= np.power(10., 20.3):
DLA_Pixels_No_xz+=1
if ((covering_fraction_y[i] >= np.power(10., 17)) & (covering_fraction_y[i] < np.power(10., 20.3))):
LLS_Pixels_No_xz+=1
DLA_covering_fraction_XY, DLA_covering_fraction_YZ, DLA_covering_fraction_XZ = DLA_Pixels_No_xy/Halo_Pixels_No, DLA_Pixels_No_yz/Halo_Pixels_No, DLA_Pixels_No_xz/Halo_Pixels_No
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
LLS_covering_fraction_XY, LLS_covering_fraction_YZ, LLS_covering_fraction_XZ = LLS_Pixels_No_xy/Halo_Pixels_No, LLS_Pixels_No_yz/Halo_Pixels_No, LLS_Pixels_No_xz/Halo_Pixels_No
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
print(Pixel_Smoothing_Scale)
print(LLS_covering_fraction_XY, LLS_covering_fraction_YZ, LLS_covering_fraction_XZ)
print(DLA_covering_fraction_XY, DLA_covering_fraction_YZ, DLA_covering_fraction_XZ)
话虽这么说，但唯一未解决的问题是某些文件的代码速度很慢，这些文件很大，而且增加 interp 值不能无限期地向更大的值继续，以获得更平滑的地图.所以，我不知道如何解决这个问题，以便在不使代码变得越来越慢的情况下制作更平滑的地图（颜色均匀填充并且地图上的任何单个区域都没有任何颜色）。对此问题的任何帮助表示赞赏。
我有一个 .txt 文件，其中包含五列数据，其中第 1 列和第 2 列对应于沿维度 1 和 2 的投影距离（像素比例为 0.50），然后第 3、4 和 5 列是投影 0 "z"、1 "x" 和 2 "y" 中任意（否则未知）函数的值。 （原始文件长约 1e6 行，但这是我设法上传的一个小版本：m11q.txt）我想将此 .txt table 转换为实际的二维地图，其中每个像素的颜色强度将对应于该像素中函数的大小。这意味着我将拥有三个平面 xy、yz 和 xz，其中每个像素都填充了一种颜色。
我对第 3、4 和 5 列的值大部分为零，其中很少有 $10^{14}$ 到 $10^{21}$ 的顺序。我有兴趣拥有三个轮廓级别。假设如果函数值小于 $10^{14}$，颜色应该是黑色，如果函数值在 $10^{14}$ 和 $10^{15.5}$ 之间，颜色应该是黄色，如果如果函数值在 $10^{15.5}$ 和 $10^{17}$ 之间，颜色应该是橙色，如果函数值在 $10^{17}$ 和 $10^{20.3}$ 之间，颜色应该是红色最后，如果超过 $10^{20.3}$，颜色应该是蓝色。换句话说，我的二维地图实际上是一组由 .txt table 生成的 5 个等高线图。
这是我的原始代码版本：
import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
nInterp = 4000
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
LAF_covering_fraction_XY = np.count_nonzero((np.power(10, 14) <= z0) & (z0 < np.power(10, 17))) * 1.0 / z0.size
LAF_covering_fraction_YZ = np.count_nonzero((np.power(10, 14) <= z1) & (z1 < np.power(10, 17))) * 1.0 / z1.size
LAF_covering_fraction_XZ = np.count_nonzero((np.power(10, 14) <= z2) & (z2 < np.power(10, 17))) * 1.0 / z2.size
LAF = np.sqrt(np.power(LAF_covering_fraction_XY, 2)+np.power(LAF_covering_fraction_YZ, 2)+np.power(LAF_covering_fraction_XZ, 2))/np.sqrt(3)
LAF_sigma = max([abs(LAF-LAF_covering_fraction_XY),abs(LAF-LAF_covering_fraction_YZ),abs(LAF-LAF_covering_fraction_XZ)])
LLS_covering_fraction_XY = np.count_nonzero((np.power(10, 17) <= z0) & (z0 < np.power(10, 20.3))) * 1.0 / z0.size
LLS_covering_fraction_YZ = np.count_nonzero((np.power(10, 17) <= z1) & (z1 < np.power(10, 20.3))) * 1.0 / z1.size
LLS_covering_fraction_XZ = np.count_nonzero((np.power(10, 17) <= z2) & (z2 < np.power(10, 20.3))) * 1.0 / z2.size
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
DLA_covering_fraction_XY = np.count_nonzero(np.power(10, 20.3) <= z0) * 1.0 / z0.size
DLA_covering_fraction_YZ = np.count_nonzero(np.power(10, 20.3) <= z1) * 1.0 / z1.size
DLA_covering_fraction_XZ = np.count_nonzero(np.power(10, 20.3) <= z2) * 1.0 / z2.size
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
z0plot = plt.subplot(221)
plt.imshow(z0, vmin=0.0001+covering_fraction_z.min(), vmax=covering_fraction_z.max(), origin='lower',
extent=[x.min(), x.max(), y.min(), y.max()],
norm=colors.LogNorm())
plt.ylabel("Y [kpc]")
z0plot.set(aspect=1)
axparameters = plt.subplot(222)
frame = plt.gca()
frame.axes.get_xaxis().set_ticks([])
frame.axes.get_yaxis().set_ticks([])
axparameters.set(aspect=1)
plt.text(0.02, 0.2, 'LAF=%.5f$\pm$%.5f \n\n\n LLS=%.5f$\pm$%.5f \n\n\n DLA=%.5f$\pm$%.5f' % (LAF, LAF_sigma, LLS, LLS_sigma, DLA, DLA_sigma))
plt.colorbar()
plt.tight_layout()
z1plot = plt.subplot(223)
plt.imshow(z1, vmin=0.0001+covering_fraction_y.min(), vmax=covering_fraction_y.max(), origin='lower',
extent=[x.min(), x.max(), y.min(), y.max()],
norm=colors.LogNorm())
plt.xlabel("X [kpc]")
plt.ylabel("Z [kpc]")
z1plot.set(aspect=1)
z2plot = plt.subplot(224)
plt.imshow(z2, vmin=0.0001+covering_fraction_x.min(), vmax=covering_fraction_x.max(), origin='lower',
extent=[x.min(), x.max(), y.min(), y.max()],
norm=colors.LogNorm())
plt.xlabel("Y [kpc]")
z2plot.set(aspect=1)
myfig = fig.Figure(figsize=(8, 8), dpi=900, facecolor='w', edgecolor='k', tight_layout=True)
plt.suptitle('HI Column Density Maps', size=20)
plt.tight_layout(rect=[0, 0, 1, 0.9])
plt.savefig('m11q.png', dpi=900)
plt.show()
生成以下地图
这是 Xevaquor 编写的代码，作为对我原来的回答 post:
# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
#import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
import matplotlib.ticker as ticker
from matplotlib import gridspec
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
def fmt(x, pos):
return '10e{}'.format(np.log10(x))
nInterp = 4300
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
LAF_covering_fraction_XY = np.count_nonzero((np.power(10, 14) <= z0) & (z0 < np.power(10, 17))) * 1.0 / z0.size
LAF_covering_fraction_YZ = np.count_nonzero((np.power(10, 14) <= z1) & (z1 < np.power(10, 17))) * 1.0 / z1.size
LAF_covering_fraction_XZ = np.count_nonzero((np.power(10, 14) <= z2) & (z2 < np.power(10, 17))) * 1.0 / z2.size
LAF = np.sqrt(np.power(LAF_covering_fraction_XY, 2)+np.power(LAF_covering_fraction_YZ, 2)+np.power(LAF_covering_fraction_XZ, 2))/np.sqrt(3)
LAF_sigma = max([abs(LAF-LAF_covering_fraction_XY),abs(LAF-LAF_covering_fraction_YZ),abs(LAF-LAF_covering_fraction_XZ)])
LLS_covering_fraction_XY = np.count_nonzero((np.power(10, 17) <= z0) & (z0 < np.power(10, 20.3))) * 1.0 / z0.size
LLS_covering_fraction_YZ = np.count_nonzero((np.power(10, 17) <= z1) & (z1 < np.power(10, 20.3))) * 1.0 / z1.size
LLS_covering_fraction_XZ = np.count_nonzero((np.power(10, 17) <= z2) & (z2 < np.power(10, 20.3))) * 1.0 / z2.size
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
DLA_covering_fraction_XY = np.count_nonzero(np.power(10, 20.3) <= z0) * 1.0 / z0.size
DLA_covering_fraction_YZ = np.count_nonzero(np.power(10, 20.3) <= z1) * 1.0 / z1.size
DLA_covering_fraction_XZ = np.count_nonzero(np.power(10, 20.3) <= z2) * 1.0 / z2.size
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
fig = plt.figure(figsize=(8, 8))
gs = gridspec.GridSpec(2, 3, width_ratios=[5, 5, 1], height_ratios=[1,1,1])
ax0 = plt.subplot(gs[0])
ax1 = plt.subplot(gs[3])
ax2 = plt.subplot(gs[4])
#cmap = colors.ListedColormap(['#ffdddd','#000000', '#ff9999', '#ff7777', '#ff5555', '#ff3333' ,'#ff1111'])
#'gray','green','pink','brown','blue','violet','yellow','orange',
cmap = colors.ListedColormap(['#000000', 'gray','#ffdddd', '#ff9999', '#ff7777', '#ff5555', '#ff3333', '#ff1111', 'blue'])
boundaries = [np.power(10., 10), np.power(10., 11), np.power(10., 12), np.power(10., 13), np.power(10., 14), np.power(10., 15), np.power(10., 16), np.power(10.,17), np.power(10.,20.3), np.power(10.,25)]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
pcmz0 = ax0.pcolormesh(xi, yi, z0, cmap=cmap, norm=norm)
ax0.set_xticks([x.min(), x.max()])
ax0.set_yticks([y.min(), y.max()])
ax0.set_xlim([x.min(), x.max()])
ax0.set_ylim([y.min(), y.max()])
ax0.set(aspect=1)
pcmz1 = ax1.pcolormesh(xi, yi, z1, cmap=cmap, norm=norm)
ax1.set_xticks([x.min(), x.max()])
ax1.set_yticks([y.min(), y.max()])
ax1.set_xlim([x.min(), x.max()])
ax1.set_ylim([y.min(), y.max()])
ax1.set(aspect=1)
pcmz2 = ax2.pcolormesh(xi, yi, z2, cmap=cmap, norm=norm)
ax2.set_xticks([x.min(), x.max()])
ax2.set_yticks([y.min(), y.max()])
ax2.set_xlim([x.min(), x.max()])
ax2.set_ylim([y.min(), y.max()])
ax2.set(aspect=1)
axparameters = plt.subplot(gs[1])
frame = plt.gca()
frame.axes.get_xaxis().set_ticks([])
frame.axes.get_yaxis().set_ticks([])
axparameters.set(aspect=1)
plt.text(0.1, 0.2, 'LAF=%.5f$\pm$%.5f \n\n\n LLS=%.5f$\pm$%.5f \n\n\n DLA=%.5f$\pm$%.5f' % (LAF, LAF_sigma, LLS, LLS_sigma, DLA, DLA_sigma), size=10)
axc = plt.subplot(gs[:,2])
plt.colorbar(pcmz1, cax=axc,format=ticker.FuncFormatter(fmt))
ax0.set_ylabel("Y [kpc]")
ax1.set_xlabel("X [kpc]")
ax1.set_ylabel("Z [kpc]")
ax2.set_xlabel("Y [kpc]")
fig.suptitle('HI Column Density Maps', size=20)
gs.tight_layout(fig, rect=[0, 0, 1, 0.9])
plt.savefig('m11q.png', dpi=900)
plt.show()
生成以下地图：
这里是代码的最终版本，除了在 post:
的答案（由我）中提到的一些警告外，它完成了我的所有目标import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
#import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
import matplotlib.ticker as ticker
from matplotlib import gridspec
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
def fmt(x, pos):
return '10e{}'.format(np.log10(x))
nInterp = 3000
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
Pixel_Smoothing_Scale = 0.497
Halo_Pixels_No = (x.max()-x.min())*(y.max()-y.min())/ np.power(Pixel_Smoothing_Scale, 2)
DLA_Pixels_No_xy, DLA_Pixels_No_yz, DLA_Pixels_No_xz = 0, 0, 0
LLS_Pixels_No_xy, LLS_Pixels_No_yz, LLS_Pixels_No_xz = 0, 0, 0
LAF_Pixels_No_xy, LAF_Pixels_No_yz, LAF_Pixels_No_xz = 0, 0, 0
for i in np.arange(0, np.size(x)):
if covering_fraction_z[i] >= np.power(10., 20.3):
DLA_Pixels_No_xy+=1
if ((covering_fraction_z[i] >= np.power(10., 17)) & (covering_fraction_z[i] < np.power(10., 20.3))):
LLS_Pixels_No_xy+=1
if ((covering_fraction_z[i] >= np.power(10., 14)) & (covering_fraction_z[i] < np.power(10., 17))):
LAF_Pixels_No_xy+=1
if covering_fraction_x[i] >= np.power(10., 20.3):
DLA_Pixels_No_yz+=1
if ((covering_fraction_x[i] >= np.power(10., 17)) & (covering_fraction_x[i] < np.power(10., 20.3))):
LLS_Pixels_No_yz+=1
if ((covering_fraction_x[i] >= np.power(10., 14)) & (covering_fraction_x[i] < np.power(10., 17))):
LAF_Pixels_No_yz+=1
if covering_fraction_y[i] >= np.power(10., 20.3):
DLA_Pixels_No_xz+=1
if ((covering_fraction_y[i] >= np.power(10., 17)) & (covering_fraction_y[i] < np.power(10., 20.3))):
LLS_Pixels_No_xz+=1
if ((covering_fraction_y[i] >= np.power(10., 14)) & (covering_fraction_y[i] < np.power(10., 17))):
LAF_Pixels_No_xz+=1
DLA_covering_fraction_XY, DLA_covering_fraction_YZ, DLA_covering_fraction_XZ = DLA_Pixels_No_xy/Halo_Pixels_No, DLA_Pixels_No_yz/Halo_Pixels_No, DLA_Pixels_No_xz/Halo_Pixels_No
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
LLS_covering_fraction_XY, LLS_covering_fraction_YZ, LLS_covering_fraction_XZ = LLS_Pixels_No_xy/Halo_Pixels_No, LLS_Pixels_No_yz/Halo_Pixels_No, LLS_Pixels_No_xz/Halo_Pixels_No
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
LAF_covering_fraction_XY, LAF_covering_fraction_YZ, LAF_covering_fraction_XZ = LAF_Pixels_No_xy/Halo_Pixels_No, LAF_Pixels_No_yz/Halo_Pixels_No, LAF_Pixels_No_xz/Halo_Pixels_No
LAF = np.sqrt(np.power(LAF_covering_fraction_XY, 2)+np.power(LAF_covering_fraction_YZ, 2)+np.power(LAF_covering_fraction_XZ, 2))/np.sqrt(3)
LAF_sigma = max([abs(LAF-LAF_covering_fraction_XY),abs(LAF-LAF_covering_fraction_YZ),abs(LAF-LAF_covering_fraction_XZ)])
fig = plt.figure(figsize=(8, 8))
gs = gridspec.GridSpec(2, 3, width_ratios=[5, 5, 1], height_ratios=[1,1,1])
ax0 = plt.subplot(gs[0])
ax1 = plt.subplot(gs[3])
ax2 = plt.subplot(gs[4])
#cmap = colors.ListedColormap(['#ffdddd','#000000', '#ff9999', '#ff7777', '#ff5555', '#ff3333' ,'#ff1111'])
#'gray','green','pink','brown','blue','violet','yellow','orange',
cmap = colors.ListedColormap(['#ffdddd', '#ff7777', '#ff3333', 'yellow', 'orange','green', 'red', 'blue'])
boundaries = [np.power(10., 4), np.power(10., 6), np.power(10., 8), np.power(10., 10), np.power(10., 12), np.power(10., 14), np.power(10.,17), np.power(10.,20.3), np.power(10.,25)]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
pcmz0 = ax0.pcolormesh(xi, yi, z0, cmap=cmap, norm=norm)
ax0.set_xticks([x.min(), x.max()])
ax0.set_yticks([y.min(), y.max()])
ax0.set_xlim([x.min(), x.max()])
ax0.set_ylim([y.min(), y.max()])
ax0.set(aspect=1)
pcmz1 = ax1.pcolormesh(xi, yi, z1, cmap=cmap, norm=norm)
ax1.set_xticks([x.min(), x.max()])
ax1.set_yticks([y.min(), y.max()])
ax1.set_xlim([x.min(), x.max()])
ax1.set_ylim([y.min(), y.max()])
ax1.set(aspect=1)
pcmz2 = ax2.pcolormesh(xi, yi, z2, cmap=cmap, norm=norm)
ax2.set_xticks([x.min(), x.max()])
ax2.set_yticks([y.min(), y.max()])
ax2.set_xlim([x.min(), x.max()])
ax2.set_ylim([y.min(), y.max()])
ax2.set(aspect=1)
axparameters = plt.subplot(gs[1])
frame = plt.gca()
frame.axes.get_xaxis().set_ticks([])
frame.axes.get_yaxis().set_ticks([])
axparameters.set(aspect=1)
plt.text(0.1, 0.2, 'LAF=%.6f$\pm$%.6f \n\n\n LLS=%.6f$\pm$%.6f \n\n\n DLA=%.6f$\pm$%.6f' % (LAF, LAF_sigma, LLS, LLS_sigma, DLA, DLA_sigma), size=9)
axc = plt.subplot(gs[:,2])
plt.colorbar(pcmz1, cax=axc,format=ticker.FuncFormatter(fmt))
ax0.set_ylabel("Y [kpc]")
ax1.set_xlabel("X [kpc]")
ax1.set_ylabel("Z [kpc]")
ax2.set_xlabel("Y [kpc]")
fig.suptitle('HI Column Density Maps', size=20)
gs.tight_layout(fig, rect=[0, 0, 1, 0.9])
plt.savefig('m11q.png', dpi=900)
plt.show()
生成以下地图：
(1) 和(2) 是相关的。您可以定义自定义配色方案。对于任意值，它可以是例如：
cmap = colors.ListedColormap(['red', '#000000','#444444', '#666666', '#ffffff', 'blue', 'orange'])
boundaries = [-1, -0.9, -0.6, -0.3, 0, 0.3, 0.6, 1]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
这将按照以下方式映射值：
- -1.0 x < -0.9 到红色
- -0.9 < x < -0.6 至黑色 (#000000)
- -0.6 < x < -0.3 至 #444444
- 等等
颜色 i 将用于边界 i 和 i+1 之间的值。颜色可以通过名称 ('red'、'green')、HTML 代码 ('#ffaa44'、'#441188') 或 RGB 元组 ((0.2, 0.9, 0.45)) 指定。
完整示例：
# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
import scipy.interpolate
#import matplotlib.figure as fig
import matplotlib.cm as cm
import matplotlib.colors as colors
import matplotlib.colorbar as cb
import matplotlib.ticker as tk
import pylab as pl
import matplotlib.ticker as ticker
from matplotlib import gridspec
x, y, covering_fraction_z, covering_fraction_x, covering_fraction_y = np.genfromtxt(r'm11q.txt', unpack=True)
def fmt(x, pos):
return '10e{}'.format(np.log10(x))
nInterp = 402
xi, yi = np.linspace(x.min(), x.max(), nInterp), np.linspace(y.min(), y.max(), nInterp)
xi, yi = np.meshgrid(xi, yi)
z0 = scipy.interpolate.griddata((x, y), covering_fraction_z, (xi, yi), method='nearest')
z1 = scipy.interpolate.griddata((x, y), covering_fraction_x, (xi, yi), method='nearest')
z2 = scipy.interpolate.griddata((x, y), covering_fraction_y, (xi, yi), method='nearest')
fig = plt.figure(figsize=(6, 6))
gs = gridspec.GridSpec(2, 3, width_ratios=[4, 4, 1], height_ratios=[1,1,1])
ax0 = plt.subplot(gs[0])
ax1 = plt.subplot(gs[3])
ax2 = plt.subplot(gs[4])
cmap = colors.ListedColormap(['#ffdddd','#000000', '#ff9999', '#ff7777', '#ff5555', '#ff3333' ,'#ff1111'])
boundaries = [np.power(10., 14), np.power(10., 15.5), np.power(10.,17), np.power(10.,20.3), np.power(10.,25)]
norm = colors.BoundaryNorm(boundaries, cmap.N, clip=True)
pcmz0 = ax0.pcolormesh(xi, yi, z0, cmap=cmap, norm=norm)
ax0.set_xticks([x.min(), x.max()]) #arbitrary, also for Y axis ticks
ax0.set_xlim([x.min(), x.max()])
ax0.set_ylim([y.min(), y.max()])
pcmz1 = ax1.pcolormesh(xi, yi, z1, cmap=cmap, norm=norm)
ax1.set_xticks([x.min(), x.max()])
ax1.set_xlim([x.min(), x.max()])
ax1.set_ylim([y.min(), y.max()])
pcmz2 = ax2.pcolormesh(xi, yi, z2, cmap=cmap, norm=norm)
ax2.set_xticks([x.min(), x.max()])
ax2.set_xlim([x.min(), x.max()])
ax2.set_ylim([y.min(), y.max()])
axc = plt.subplot(gs[:,2])
plt.colorbar(pcmz1, cax=axc,format=ticker.FuncFormatter(fmt))
ax0.set_ylabel("Y [kpc]")
ax1.set_xlabel("X [kpc]")
ax1.set_ylabel("Z [kpc]")
ax2.set_xlabel("Y [kpc]")
plt.tight_layout()
plt.savefig('m11q.png', dpi=900)
plt.show()
对于 (3)。有很多方法可以做到这一点。例如：
fraction = np.count_nonzero((np.power(10, 14) <= Z) & (Z < np.power(10, 15.5)) * 1.0 / Z.size
将为您提供 [10e14, 10e15.5)
范围内的部分值(4) 取决于创建子图的方法。对于 GridSpec，如果您不明确声明它，它将不会出现。
这是问题更新版本的答案。最初我无法为我的地图制作颜色条。 Xevaquor 通过上述回答帮助我解决了这个问题。所以，答案是由 Xevaquor 合法发布的，现在被标记为 "the answer"。然而，当我花时间解决一些覆盖分数的问题时，我终于能够从数学上理解这个问题，然后尝试在 python 代码中用数字来解决它。我的问题是我没有意识到我所追求的实际上是像素数量的比率。这就是为什么我必须计算位于列密度范围内的像素数，然后将其除以进行地图每个投影的像素总数。由于已知平滑比例为 $0.497 kpc$，我们知道像素总数将是地图的物理面积 $kpc^{2}$ 除以每个像素的面积即 $0.497^{2}$。因此，可以将以下代码片段添加到原始 and/or Xevaquor 的代码中以正确说明覆盖分数：
Pixel_Smoothing_Scale = 0.497
Halo_Pixels_No = (x.max()-x.min())*(y.max()-y.min())/ np.power(Pixel_Smoothing_Scale, 2)
DLA_Pixels_No_xy, DLA_Pixels_No_yz, DLA_Pixels_No_xz = 0, 0, 0
LLS_Pixels_No_xy, LLS_Pixels_No_yz, LLS_Pixels_No_xz = 0, 0, 0
for i in np.arange(0, np.size(x)):
if covering_fraction_z[i] >= np.power(10., 20.3):
DLA_Pixels_No_xy+=1
if ((covering_fraction_z[i] >= np.power(10., 17)) & (covering_fraction_z[i] < np.power(10., 20.3))):
LLS_Pixels_No_xy+=1
if covering_fraction_x[i] >= np.power(10., 20.3):
DLA_Pixels_No_yz+=1
if ((covering_fraction_x[i] >= np.power(10., 17)) & (covering_fraction_x[i] < np.power(10., 20.3))):
LLS_Pixels_No_yz+=1
if covering_fraction_y[i] >= np.power(10., 20.3):
DLA_Pixels_No_xz+=1
if ((covering_fraction_y[i] >= np.power(10., 17)) & (covering_fraction_y[i] < np.power(10., 20.3))):
LLS_Pixels_No_xz+=1
DLA_covering_fraction_XY, DLA_covering_fraction_YZ, DLA_covering_fraction_XZ = DLA_Pixels_No_xy/Halo_Pixels_No, DLA_Pixels_No_yz/Halo_Pixels_No, DLA_Pixels_No_xz/Halo_Pixels_No
DLA = np.sqrt(np.power(DLA_covering_fraction_XY, 2)+np.power(DLA_covering_fraction_YZ, 2)+np.power(DLA_covering_fraction_XZ, 2))/np.sqrt(3)
DLA_sigma = max([abs(DLA-DLA_covering_fraction_XY),abs(DLA-DLA_covering_fraction_YZ),abs(DLA-DLA_covering_fraction_XZ)])
LLS_covering_fraction_XY, LLS_covering_fraction_YZ, LLS_covering_fraction_XZ = LLS_Pixels_No_xy/Halo_Pixels_No, LLS_Pixels_No_yz/Halo_Pixels_No, LLS_Pixels_No_xz/Halo_Pixels_No
LLS = np.sqrt(np.power(LLS_covering_fraction_XY, 2)+np.power(LLS_covering_fraction_YZ, 2)+np.power(LLS_covering_fraction_XZ, 2))/np.sqrt(3)
LLS_sigma = max([abs(LLS-LLS_covering_fraction_XY),abs(LLS-LLS_covering_fraction_YZ),abs(LLS-LLS_covering_fraction_XZ)])
print(Pixel_Smoothing_Scale)
print(LLS_covering_fraction_XY, LLS_covering_fraction_YZ, LLS_covering_fraction_XZ)
print(DLA_covering_fraction_XY, DLA_covering_fraction_YZ, DLA_covering_fraction_XZ)
话虽这么说，但唯一未解决的问题是某些文件的代码速度很慢，这些文件很大，而且增加 interp 值不能无限期地向更大的值继续，以获得更平滑的地图.所以，我不知道如何解决这个问题，以便在不使代码变得越来越慢的情况下制作更平滑的地图（颜色均匀填充并且地图上的任何单个区域都没有任何颜色）。对此问题的任何帮助表示赞赏。