MPPNPy/abund_plot_z.py at master · Bigash94/MPPNPy · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
def prod_frac(data,z,label1,label2):
    import numpy as np
    import pandas as pd
    import matplotlib.pyplot as plt

    df = pd.read_csv('mass_fractions/'+data,skiprows=[0,1],usecols=[1,2,3],sep='\&',engine='python')
    iso1=df[' specie '].values
    pf1=df[' '+label1+'     '].values
    pf2=df[' '+label2+'     '].values

    iso=[]
    for i in iso1:
        iso.append(i[:3])

    A = np.array(iso)
    Y = np.array(pf1)
    Y2 = np.array(pf2)
    Auq = pd.unique(A)

    xsum1 = []
    xsum2 = []
    for j in Auq:
        dummy11 = sum(Y[np.where(A==j)])
        xsum1.append(dummy11)
        dummy1 = sum(Y2[np.where(A==j)])
        xsum2.append(dummy1)

    df1 = pd.read_csv(z+'E-02GN93.ppn',skiprows=[4,6,7,283,284,285],sep='\s*',engine='python',header=None)
    ele=df1[0].values
    pf=df1[3].values

    A1 = np.array(ele)
    Y1 = np.array(pf)
    Auq1 = np.unique(A1)

    xsum = []
    for j in Auq1:
        dummy = sum(Y1[np.where(A1==j)])
        xsum.append(dummy)

    x = Auq1
    y1 = [a/b for a,b in zip(xsum1,xsum)]
    y2 = [a/b for a,b in zip(xsum2,xsum)]

    df2 = pd.read_csv('isotopi_m3p0z2m2_000_20190118_16311.txt',sep='\s+',engine='python')
    iso3=df2['Z'].values
    pf3=df2['TDU_14'].values
    df3 = pd.read_csv('isotopi_m3p0zsun_000_20190118_16409.txt',sep='\s+',engine='python')
    ini = df3['INI'].values

    A3 = np.array(iso3)
    Y3 = np.array(pf3)
    Auq3 = np.unique(A3)

    xsum3 = []
    xsum4 = []
    for j in Auq3:
        dummy3 = sum(Y3[np.where(A3==j)])
        xsum3.append(dummy3)

        dummy4 = sum(ini[np.where(A3==j)])
        xsum4.append(dummy4)

    x1 = Auq3
    y3 = [a/b for a,b in zip(xsum3,xsum4)]

    Z_list=[]
    for i in Auq:
        Z_list.append(i[1]+i[2].lower())
    y77 = 28*[26,22,18]
    for label3, x7, y7 in zip(Z_list, x, y77):
        plt.annotate(label3,xy=(x7, y7), xytext=(-5, 10),
                    textcoords='offset points', fontsize=15)
    plt.axhline(y=1.,linewidth=1, color='k')
    plt.semilogy(x,y1,'rs',markersize=8.,label=label1)
    plt.semilogy(x,y2,'g^',markersize=8.,label=label2)
    if data=='z2m2_Ritter18.dat':
        plt.semilogy(x1,y3,'bo',label='FRUITY',markersize=8.)
    plt.tick_params(axis='both', which='major', length=8, width=2, labelsize=20)
    plt.tick_params(axis='both', which='minor', length=6, width=2, labelsize=20)
    plt.xlim(0,84)
    plt.ylim(0.1,35)
    plt.xlabel('Atomic Number', fontsize=30)
    plt.ylabel('Production Factors', fontsize=30)
    plt.legend(prop={'size': 30},loc='lower right')

    return(plt.show(),x,y1,y2)

#prod_frac(data='z1m2.dat',z='1.0',label1='2M$\odot$',label2='3M$\odot$')
#prod_frac(data='z1m2_he07.dat',z='1.0',label1='2M$\odot$',label2='3M$\odot$')
#prod_frac(data='z2m2_he07.dat',z='2.0',label1='2M$\odot$',label2='3M$\odot$')
#prod_frac(data='z3m2.dat',z='3.0',label1='2M$\odot$',label2='3M$\odot$')
#prod_frac(data='z3p5m2.dat',z='3.5',label1='2M$\odot$',label2='3M$\odot$')
prod_frac(data='z2m2_Ritter18.dat',z='2.0',label1='m3z2m2',label2='RI18')