本文整理汇总了Python中spectral_cube.SpectralCube.subcube_from_mask方法的典型用法代码示例。如果您正苦于以下问题：Python SpectralCube.subcube_from_mask方法的具体用法？Python SpectralCube.subcube_from_mask怎么用？Python SpectralCube.subcube_from_mask使用的例子？那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类spectral_cube.SpectralCube的用法示例。
在下文中一共展示了SpectralCube.subcube_from_mask方法的1个代码示例，这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞，您的评价将有助于系统推荐出更棒的Python代码示例。
示例1: calc_structure_fcn
# 需要导入模块: from spectral_cube import SpectralCube [as 别名]
# 或者: from spectral_cube.SpectralCube import subcube_from_mask [as 别名]
def calc_structure_fcn(catalog='cohrs_finalcatalog_clumpflux_medaxis.fits',
bootiter=0,doPlot=False):
cat = Table.read(catalog)
if 'sf_offset' not in cat.keys():
keylist = ['sf_offset','sf_index','sf_ngood',
'sf_index_err','sf_offset_err']
for thiskey in keylist:
c = Column(np.zeros(len(cat))+np.nan,name=thiskey)
cat.add_column(c)
current_open_file = ''
for cloud in cat:
orig_file = cloud['orig_file']+'.fits'
asgn_file = cloud['orig_file']+'_fasgn.fits'
if os.path.isfile(datadir+'COHRS/'+orig_file):
if current_open_file != datadir+'COHRS/'+orig_file:
hdu = fits.open(datadir+'COHRS/'+orig_file,memmap=False)
cat.write('cohrs_structurefunc.fits',overwrite=True)
w = wcs.WCS(hdu[0].header)
hdr2 = w.to_header()
hdr2['BMAJ'] = 15./3600
hdr2['BMIN'] = 15./3600
hdr2['BPA'] = 0.
co = SpectralCube(hdu[0].data,w,header=hdr2)
hdu = fits.open(datadir+'ASSIGNMENTS/'+asgn_file,memmap=False)
w = wcs.WCS(hdu[0].header)
hdr2 = w.to_header()
hdr2['BMAJ'] = 15./3600
hdr2['BMIN'] = 15./3600
hdr2['BPA'] = 0.
asgn = SpectralCube(hdu[0].data,w,header=hdr2)
# masked_co = co.with_mask(asgn>0*u.dimensionless_unscaled)
# moment = masked_co.moment(0)
current_open_file = datadir+'COHRS/'+orig_file
cat.write('output_catalog2.fits',overwrite=True)
print(cloud['cloud_id'])
mask = (asgn == cloud['cloud_id']*u.dimensionless_unscaled)
subcube = co.subcube_from_mask(mask)
if subcube.shape[0] > 15:
# zeros = np.zeros_like(subcube) # np.random.randn(*subcube.shape)*0.5
# concatdata = np.vstack([subcube.filled_data[:],zeros])
# hdr2 = subcube.wcs.to_header()
# hdr2['BMAJ'] = 15./3600
# hdr2['BMIN'] = 15./3600
# hdr2['BPA'] = 0.
# newcube = SpectralCube(concatdata,subcube.wcs,header=hdr2)
pcaobj = pca.PCA(subcube,distance=cloud['distance']*u.pc)
try:
pcaobj.run(min_eigval=0.25,verbose=False,mean_sub=True)
cloud['sf_index']= pcaobj.index
cloud['sf_offset'] = pcaobj.intercept.value
cloud['sf_ngood'] = np.min([np.isfinite(pcaobj.spatial_width).sum(),
np.isfinite(pcaobj.spectral_width).sum()])
cloud['sf_index_err'] = (pcaobj.index_error_range[1]-pcaobj.index_error_range[0])*0.5
cloud['sf_offset_err'] = ((pcaobj.intercept_error_range[1]-pcaobj.intercept_error_range[0])*0.5).value
print('{0} +/- {1}'.format(cloud['sf_index'],
cloud['sf_index_err']),
cloud['sf_offset'])
except:
pass
# r, dv = cloudpca.structure_function(subcube,
# meanCorrection=True,
# nScales=nscale,
# noiseScales=nscale/2)
# idx = np.isfinite(r) * np.isfinite(dv)
# n_good = np.sum(idx)
# if n_good >3:
# p = np.polyfit(np.log10(r[idx])+
# np.log10(2.91e-5*cloud['distance']),
# np.log10(dv[idx]),1)
# pboot = np.zeros((2,bootiter))
# if bootiter>0:
# indices= (np.where(idx))[0]
# length = len(indices)
# for ctr in np.arange(bootiter):
# bootidx = np.random.choice(indices,length,True)
# pboot[:,ctr] = np.polyfit(np.log10(r[bootidx])+
# np.log10(2.91e-5*
# cloud['distance']),
# np.log10(dv[bootidx]),1)
# cloud['sf_index_err']=0.5*(\
# np.percentile(pboot[0,:],84.13)-\
# np.percentile(pboot[0,:],15.87))
# cloud['sf_offset_err']=0.5*(\
# np.percentile(pboot[1,:],84.13)-\
# np.percentile(pboot[1,:],15.87))
# if doPlot:
# plt.clf()
# x = np.log10(r[idx])+\
#.........这里部分代码省略.........