当前位置: 首页>>代码示例>>Python>>正文


Python MCUtils.area方法代码示例

本文整理汇总了Python中MCUtils.area方法的典型用法代码示例。如果您正苦于以下问题:Python MCUtils.area方法的具体用法?Python MCUtils.area怎么用?Python MCUtils.area使用的例子?那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在MCUtils的用法示例。


在下文中一共展示了MCUtils.area方法的7个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的Python代码示例。

示例1: sigmaclip_bg

# 需要导入模块: import MCUtils [as 别名]
# 或者: from MCUtils import area [as 别名]
def sigmaclip_bg(data,radius,annulus,skypos,maxiter=10,sigmaclip=3.,
				 gausslim=50.,verbose=0,pixsz=0.000416666666666667):
	"""Produce an estimate of background counts within an aperture (radius)
	using a sigma clipping method for extracting the background from an
	annulus.
	This attempts to reproduce the calcuations of the backcalc() function in
	mapaps/poissonbg.c of the mission pipeline. (Probably written by Ted Wyder.)
	"""
	# FIXME: Does not apply response!

	# This cut is now handled by the ugly loop below, which barely dodges a
	# conceptula issue about fractional pixels...
	#ix = np.where((d>annulus[0]) & (d<annulus[1]))

	imsz=gxt.deg2pix(skypos,[annulus[1]*2,annulus[1]*2])
	wcs=define_wcs(skypos,[annulus[1]*2,annulus[1]*2])
	foc_ra,foc_dec=wcs.sip_pix2foc(wcs.wcs_world2pix(data['ra'],data['dec'],1),1)
	H,xedges,yedges=np.histogram2d(foc_ra-0.5,foc_dec-0.5,bins=imsz,
				       range=([ [0,imsz[0]],[0,imsz[1]]]))

	# Convert Gaussian sigma to a probability
	problim = 0.5*scipy.special.erfc(sigmaclip/np.sqrt(2.0))

	# Mask out non-annulus regions... there's probalby a more pythonic way
	bgimg=np.copy(H)
	for i in range(H.shape[0]):
		for j in range(H.shape[1]):
			# Add a little buffer to account for pixel widths?
			# FIXME? including everything within the annulus...
#			if (mc.distance(H.shape[0]/2.,H.shape[1]/2.,i,j)<annulus[0]/pixsz or
			if	mc.distance(H.shape[0]/2.,H.shape[1]/2.,i,j)>annulus[1]/pixsz:#):

				bgimg[i,j]=-1

	ix=np.where(bgimg>=0)
	m,s=bgimg[ix].mean(),bgimg[ix].std()
	d = 1.
	for i in range(maxiter):
		if d<=10e-5 or m<2:
			continue
		if m>=gausslim:
			# Mask anything outside of 3 sigma from the mean (of unmasked data)
			klim=m+sigmaclip*np.sqrt(m)#s
			klo=m-sigmaclip*np.sqrt(m)#s
			if verbose:
				print 'Gaussian cut: {klo} to {klim}'.format(klo=klo,klim=klim)
		else:
			klim = scipy.special.gammainccinv(m,problim)
			klo = -1 # None
			if verbose:
				print 'Poisson cut: {klo} to {klim}'.format(klo=klo,klim=klim)
		ix = np.where((bgimg>=klim) | (bgimg<=klo))
		bgimg[ix]=-1
		ix=np.where(bgimg>=0)
		d = np.abs((bgimg[ix].mean()-m)/m)# - 1)
		m,s=bgimg[ix].mean(),bgimg[ix].std()
	ix = np.where(bgimg>=0)
	return mc.area(radius)*bgimg[ix].mean()/mc.area(pixsz)
开发者ID:jvc2688,项目名称:gPhoton,代码行数:60,代码来源:imagetools.py

示例2: error

# 需要导入模块: import MCUtils [as 别名]
# 或者: from MCUtils import area [as 别名]
def error(data,band,radius,annulus):
    N_a = 1
    N_b0 = (mc.area(annulus[1])-mc.area(annulus[0]))/mc.area(radius)
    N_b = data[band]['bg_eff_area']/mc.area(radius)
    B0 = data[band]['bg']
    B = data[band]['bg_cheese']
    S = gt.mag2counts(data[band]['mag'],band)*data[band]['t_eff']
    s2 = {'bg_cheese_err':(S-B)+(N_a+(N_a**2.)/N_b),
          'bg_err':(S-B0)+(N_a+(N_a**2.)/N_b0)}
    return s2
开发者ID:cmillion,项目名称:gPhoton,代码行数:12,代码来源:script.py

示例3: cheese_bg_area

# 需要导入模块: import MCUtils [as 别名]
# 或者: from MCUtils import area [as 别名]
def cheese_bg_area(band,ra0,dec0,annulus,sources,nsamples=10e5,ntests=10):
    # This is just a really naive Monte Carlo
    ratios = np.zeros(ntests)
    for i in range(ntests):
        ann_events = bg_mask_annulus(band,ra0,dec0,annulus,
             np.random.uniform(ra0-annulus[1],ra0+annulus[1],int(nsamples)),
             np.random.uniform(dec0-annulus[1],dec0+annulus[1],int(nsamples)),
             np.ones(nsamples))
        mask_events= bg_mask_sources(band,ra0,dec0,
                 ann_events[0],ann_events[1],ann_events[2],sources)
        try:
            ratios[i] = float(mask_events[2].sum())/float(ann_events[2].sum())
        except ZeroDivisionError:
            ratios[i] = 0.

    return (mc.area(annulus[1])-mc.area(annulus[0]))*ratios.mean()
开发者ID:jvc2688,项目名称:gPhoton,代码行数:18,代码来源:curvetools.py

示例4: gphot_params

# 需要导入模块: import MCUtils [as 别名]
# 或者: from MCUtils import area [as 别名]
def gphot_params(band,skypos,radius,annulus=None,
                 verbose=0.,detsize=1.25,stepsz=None,
                 trange=None,maskdepth=None,maskradius=None):
    """Populate a dict() with parameters that are constant over all bins."""
    return {'band':band,'ra0':skypos[0],'dec0':skypos[1],'skypos':skypos,
            'trange':trange,'radius':radius,'annulus':annulus,
            'stepsz':stepsz,'verbose':verbose,
            'maskdepth':maskdepth,'maskradius':maskradius,'detsize':detsize,
            'apcorrect1':gxt.apcorrect1(radius,band),
            'apcorrect2':gxt.apcorrect2(radius,band),
            'detbg':gxt.detbg(mc.area(radius),band)}
开发者ID:jvc2688,项目名称:gPhoton,代码行数:13,代码来源:curvetools.py

示例5: cheese_bg

# 需要导入模块: import MCUtils [as 别名]
# 或者: from MCUtils import area [as 别名]
def cheese_bg(band,ra0,dec0,radius,annulus,ras,decs,responses,maskdepth=20.,
              maskradius=1.5,eff_area=False,sources=False):
    """ Returns the estimate number of counts (not count rate) within the
    aperture based upon a masked background annulus.
    """
    if not sources:
        sources = bg_sources(band,ra0,dec0,annulus[1],maskdepth=maskdepth)
    bg_counts = bg_mask(band,ra0,dec0,annulus,ras,decs,responses,
                        sources)[2].sum()
    if not eff_area:
        eff_area = cheese_bg_area(band,ra0,dec0,annulus,sources)
    return mc.area(radius)*bg_counts/eff_area if eff_area else 0.
开发者ID:jvc2688,项目名称:gPhoton,代码行数:14,代码来源:curvetools.py

示例6: quickmag

# 需要导入模块: import MCUtils [as 别名]
# 或者: from MCUtils import area [as 别名]
def quickmag(band, ra0, dec0, tranges, radius, annulus=None, data={},
             stepsz=None, verbose=0, maskdepth=20.0,
             maskradius=1.5,detsize=1.25,coadd=False, photonfile=None):
    if verbose:
        mc.print_inline("Retrieving all of the target events.")
    trange = [np.array(tranges).min(),np.array(tranges).max()]
    try:
        searchradius = annulus[1]
    except TypeError:
        searchradius = radius
    data = pullphotons(band, ra0, dec0, tranges, searchradius,
                       verbose=verbose, photonfile=photonfile)
    if verbose:
        mc.print_inline("Isolating source from background.")
    angSep = mc.angularSeparation(ra0, dec0, data['ra'], data['dec'])
    if verbose:
        mc.print_inline("Binning data according to requested depth.")
    # Multiple ways of defining bins
    if coadd:
        bins = np.array(trange)
    elif stepsz:
        bins = np.append(np.arange(min(trange), max(trange), stepsz),
                                                                max(trange))
    else:
        bins = np.unique(np.array(tranges).flatten())
    # This is equivalent in function to np.digitize(data['t'],bins) except
    # that it's much, much faster. See numpy issue #2656.
    ix = np.searchsorted(bins,data['t'],"right")
    # Initialize histogrammed arrays
    # FIXME: allocate these from a dict of constructors
    lcurve_cols = ['counts', 'sources', 'bg_counts','responses',
                   'detxs', 'detys', 't0_data', 't1_data', 't_mean', 'racent',
                   'deccent']
    lcurve = {'params':gphot_params(band,[ra0,dec0],radius,annulus=annulus,
                                    verbose=verbose,
                                    detsize=detsize,stepsz=stepsz,
                                    trange=trange,maskdepth=maskdepth,
                                    maskradius=maskradius)}
    for col in lcurve_cols:
        lcurve[col] = np.zeros(len(bins)-1)
    # FIXME: Bottleneck. There's probably a way to do this without looping.
    # Don't bother looping through anything with no data.
    lcurve['bg'] = {'simple':np.zeros(len(bins)-1),
                    'cheese':np.zeros(len(bins)-1)}
    if annulus is not None:
        lcurve['bg']['sources'] = bg_sources(band,ra0,dec0,annulus[1],
                                             maskdepth=maskdepth)
        lcurve['bg']['eff_area'] = cheese_bg_area(band,ra0,dec0,annulus,
                                                  lcurve['bg']['sources'])
    else:
        lcurve['bg']['sources'] = None
        lcurve['bg']['eff_area'] = 0.
    if verbose:
        mc.print_inline("Populating histograms.")
    for cnt,i in enumerate(np.unique(ix)):
        # Exclude data outside of the bins in searchsorted.
        if i-1<0 or i==len(bins):
            continue
        if verbose:
            mc.print_inline('Binning {i} of {l}.'.format(
                                                i=cnt,l=len(np.unique(ix))))
        t_ix = np.where(ix==i)
        # TODO: Optionally limit data to specific parts of detector.
        rad_ix = np.where((angSep <= radius) & (ix == i))
        # NOTE: This checks for the dim edge case where you have photons in
        #  the annulus but not in the aperture.
        if not len(rad_ix[0]):
            continue
        lcurve['t0_data'][i-1] = data['t'][rad_ix].min()
        lcurve['t1_data'][i-1] = data['t'][rad_ix].max()
        lcurve['t_mean'][i-1] = data['t'][rad_ix].mean()
        lcurve['counts'][i-1] = len(rad_ix[0])
        lcurve['sources'][i-1] = (1./data['response'][rad_ix]).sum()
        lcurve['responses'][i-1] = data['response'][rad_ix].mean()
        lcurve['detxs'][i-1] = data['col'][rad_ix].mean()
        lcurve['detys'][i-1] = data['row'][rad_ix].mean()
        lcurve['racent'][i-1] = data['ra'][rad_ix].mean()
        lcurve['deccent'][i-1] = data['dec'][rad_ix].mean()
        if annulus is not None:
            ann_ix = np.where((angSep > annulus[0]) &
                              (angSep <= annulus[1]) & (ix == i))
            lcurve['bg_counts'][i-1] = len(ann_ix[0])
            # Background is reported as counts within the aperture
            lcurve['bg']['simple'][i-1] = (mc.area(radius) *
                (1./data['response'][ann_ix]).sum() /
                (mc.area(annulus[1])-mc.area(annulus[0])))
            lcurve['bg']['cheese'][i-1] = cheese_bg(band, ra0, dec0, radius,
                annulus, data['ra'][t_ix], data['dec'][t_ix],
                data['response'][t_ix], maskdepth=maskdepth,
                eff_area=lcurve['bg']['eff_area'],
                sources=lcurve['bg']['sources'])
        else:
            lcurve['bg_counts'][i-1]=0.
            lcurve['bg']['simple'][i-1]=0.
            lcurve['bg']['cheese'][i-1]=0.
    # Only return bins that contain data.
    ix = np.where((np.isfinite(lcurve['sources'])) &
                  (np.array(lcurve['sources']) > 0))
    lcurve['t0'] = bins[ix]
    lcurve['t1'] = bins[ix[0]+1]
#.........这里部分代码省略.........
开发者ID:jvc2688,项目名称:gPhoton,代码行数:103,代码来源:curvetools.py

示例7:

# 需要导入模块: import MCUtils [as 别名]
# 或者: from MCUtils import area [as 别名]
    data[band] = pd.read_csv('{base}{band}.csv'.format(
                                                        base=base,band=band))
    print '{band} sources: {cnt}'.format(
                                band=band,cnt=data[band]['objid'].shape[0])

"""dMag vs. Mag"""
for band in bands:
    dmag = {'gphot_cheese':(lambda band:
                        data[band]['aper4']-data[band]['mag_bgsub_cheese']),
                 'gphot_nomask':(lambda band:
                        data[band]['aper4']-data[band]['mag_bgsub']),
                 'gphot_sigma':(lambda band:
                        data[band]['aper4']-data[band]['mag_bgsub_sigmaclip']),
                 'mcat':lambda band: data[band]['aper4']-
                        gt.counts2mag(gt.mag2counts(data[band]['mag'],band)-
                        data[band]['skybg']*3600**2*mc.area(gt.aper2deg(4)),
                        band)}
bgmodekeys={'gphot_cheese':'mag_bgsub_cheese',
            'gphot_nomask':'mag_bgsub',
            'gphot_sigma':'mag_bgsub_sigmaclip',
            'mcat':'skybg'}
for bgmode in dmag.keys():
    for band in bands:
        fig = plt.figure(figsize=(8*scl,4*scl))
        fig.subplots_adjust(left=0.12,right=0.95,wspace=0.02,
                                                        bottom=0.15,top=0.9)
        dmag_err=gu.dmag_errors(100.,band,sigma=1.41)
        # Make a cut on crazy outliers in the MCAT. Also on det radius and expt.
        ix = ((data[band]['aper4']>0) & (data[band]['aper4']<30) &
              (data[band]['distance']<300) & (data[band]['t_eff']<300) &
              (np.isfinite(np.array(data[band][bgmodekeys[bgmode]]))))
开发者ID:cmillion,项目名称:gPhoton,代码行数:33,代码来源:script.py


注:本文中的MCUtils.area方法示例由纯净天空整理自Github/MSDocs等开源代码及文档管理平台,相关代码片段筛选自各路编程大神贡献的开源项目,源码版权归原作者所有,传播和使用请参考对应项目的License;未经允许,请勿转载。