aboutsummaryrefslogtreecommitdiffstats
path: root/mass_profile/fg_2500_500.py
diff options
context:
space:
mode:
authorAaron LI <aaronly.me@outlook.com>2017-02-20 12:11:25 +0800
committerAaron LI <aaronly.me@outlook.com>2017-02-20 12:11:25 +0800
commitcba36462e9e70f45341e432074c726dda5e31092 (patch)
tree0bc9885d0ff86b563e0fb5321c56d51d5687ba70 /mass_profile/fg_2500_500.py
parent2a069ed00d6f1c83153be9174c296e5540f37d30 (diff)
downloadchandra-acis-analysis-cba36462e9e70f45341e432074c726dda5e31092.tar.bz2
Move shell/python scripts from 'mass_profile/' to 'bin/'
Diffstat (limited to 'mass_profile/fg_2500_500.py')
-rwxr-xr-xmass_profile/fg_2500_500.py153
1 files changed, 0 insertions, 153 deletions
diff --git a/mass_profile/fg_2500_500.py b/mass_profile/fg_2500_500.py
deleted file mode 100755
index 67a1a11..0000000
--- a/mass_profile/fg_2500_500.py
+++ /dev/null
@@ -1,153 +0,0 @@
-#!/usr/bin/env python
-
-import sys
-import numpy
-import scipy.interpolate
-
-confidence_level=.68
-def read_file(param):
- delta=float(param[0])
-
- file_mass_center=open("mass_int_center.qdp").readlines();
- file_delta_center=open("overdensity_center.qdp").readlines();
-
- center_r=0
- center_m=0
- center_gm=0
- center_gf=0
-
-
- for i in range(0,len(file_mass_center)):
- lm=file_mass_center[i].strip();
- ld=file_delta_center[i].strip();
- r,m=lm.split()
- r,d=ld.split()
- r=float(r)
- d=float(d)
- m=float(m)
- if m<1e11:
- continue
- if d<delta:
- center_r=r
- center_m=m
- for j in open("gas_mass_int_center.qdp"):
- rgm,gm=j.strip().split()
- rgm=float(rgm)
- gm=float(gm)
- if rgm>r:
-
- center_gm=gm
- center_gf=gm/m
- break
- break
- if len(param)>1 and param[1]=='c':
- #print("%s(<r%d)=%E solar mass"%("mass",delta,center_m))
- #print("%s%d=%E kpc"%("r",delta,center_r))
- #print("%s(<r%d)=%E solar mass"%("gas mass",delta,center_gm))
- #print("%s(<r%d)=%E"%("gas fraction",delta,center_gf))
- return center_m,center_r,center_gm,center_gf,None,None,None,None
-
-
-#print(center_gm,center_gf)
- file_mass=open('summary_mass_profile.qdp').readlines()
- file_delta=open('summary_overdensity.qdp').readlines()
- file_gm=open('summary_gas_mass_profile.qdp')
-
-
- flag=True
- rlist=[]
- mlist=[]
- gmlist=[]
- gflist=[]
- old_m=0
- invalid_count=0
- for i in range(0,len(file_mass)):
- lm=file_mass[i].strip()
- ld=file_delta[i].strip()
- if lm[0]=='n':
- flag=True
- old_m=0
- continue
- if not flag:
- continue
- r,m=lm.split()
- m=float(m)
- if m<1e12:
- continue
- if m<old_m:
- invalid_count+=1
- flag=False
- continue
- r,d=ld.split()
- d=float(d)
- r=float(r)
-
- if d<delta:
- #print("%s %e"%(d,m))
- mlist.append(m)
- rlist.append(r)
- flag1=True
- while True:
- lgm=file_gm.readline().strip()
- if lgm[0]=='n':
- break
- rgm,gm=lgm.split()
- rgm=float(rgm)
- gm=float(gm)
- if rgm>r and flag1:
- gmlist.append(gm)
-
- flag1=False
- gflist.append(gm/mlist[-1])
- #print(gm,gflist[-1])
- flag=False
- old_m=m
- print("%d abnormal data dropped"%(invalid_count))
-
-
- return center_m,center_r,center_gm,center_gf,mlist,rlist,gmlist,gflist
-#center_m=numpy.mean(mlist)
-#center_r=numpy.mean(rlist)
-
-if len(sys.argv)>1:
- center_m2500,center_r2500,center_gm2500,center_gf2500,mlist2500,rlist2500,gmlist2500,gflist2500=read_file([2500,sys.argv[1]])
- center_m500,center_r500,center_gm500,center_gf500,mlist500,rlist500,gmlist500,gflist500=read_file([500,sys.argv[1]])
-else:
- center_m2500,center_r2500,center_gm2500,center_gf2500,mlist2500,rlist2500,gmlist2500,gflist2500=read_file([2500])
- center_m500,center_r500,center_gm500,center_gf500,mlist500,rlist500,gmlist500,gflist500=read_file([500])
-
-if mlist2500!=None and len(mlist2500)!=len(mlist500):
- raise Exception("Something wrong, the number of 2500 and 500 data are different")
-
-
-if mlist2500==None:
- print("gas fraction between r2500 and r500 is %E"%((center_gm500-center_gm2500)/(center_m500-center_m2500)))
- sys.exit(0)
-
-gf_2500_500=[]
-
-for i in range(0,len(mlist500)):
- if mlist500[i]-mlist2500[i]<=0:
- continue
- gf_2500_500.append((gmlist500[i]-gmlist2500[i])/(mlist500[i]-mlist2500[i]))
-
-gf_2500_500.sort();
-
-
-center_gf_2500_500=(center_gm500-center_gm2500)/(center_m500-center_m2500)
-gf_idx=-1
-
-for i in range(len(gf_2500_500)-1):
- if (center_gf_2500_500-gf_2500_500[i])*(center_gf_2500_500-gf_2500_500[i+1])<=0:
- gf_idx=i
- break
-if gf_idx==-1:
- raise Exception("Something wrong!")
-
-gflidx=int(gf_idx*(1-confidence_level))
-gfuidx=gf_idx-1+int((len(gf_2500_500)-gf_idx)*confidence_level)
-
-gferr1=gf_2500_500[gflidx]-center_gf_2500_500
-gferr2=gf_2500_500[gfuidx]-center_gf_2500_500
-
-print("gas_fraction between r2500 and r500=\t%e\t %e/+%e (1 sigma)"%(center_gf_2500_500,gferr1,gferr2))