SHARE
TWEET

Untitled

a guest Oct 21st, 2019 93 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. import numpy as np
  2. from pylab import *
  3. from pylab import loadtxt
  4. from scipy import *
  5. import matplotlib.pyplot as plt
  6. import scipy.constants as const
  7. import scipy.stats as stats
  8. from scipy import optimize
  9.  
  10. #Let's plot TTau in microns and
  11. #First column is wavelength, second is flux density, third is something we ignore.
  12. TTmic = loadtxt('/home/kmiskove/Documents/TTau_merged_microns.txt', usecols=(0,1))
  13. TTw, TTf = (TTmic[:,0],TTmic[:,1])
  14. M2V = loadtxt('/home/kmiskove/Documents/M2V_Gl806.txt', skiprows=45, usecols=(0,1))
  15. M2Vw, M2Vf = (M2V[:,0],M2V[:,1])
  16.  
  17. #Let's normalize the fluxes:
  18. TTfn = TTf/median(TTf)
  19. M2fn = M2Vf/median(M2Vf)
  20.  
  21. #This is where we specify our Av, T, and r values.
  22. Av = 0.3 #veiling
  23. T = 1500 #temp
  24. r = 0.7 #scaling
  25. #And let's define some constants:
  26. h = 6.26e-34
  27. c = 3.0e8
  28. k = 1.38e-23
  29. omega = 5.67e-8
  30.  
  31. #The following defines the planck function
  32. def planck(wave, T):
  33.     """
  34.     Input: Wavelength values
  35.     Returns: B_nu(wave,temp), the Planck distribution function in SI units
  36.     """
  37.     numer = 2.*h*(c)**2/(wave*1.e-6)**5
  38.     denom = exp(h*c/(k*T*(wave*1.e-6)))-1.
  39.     return numer/denom
  40.  
  41. pl = planck(M2Vw,T)
  42. pln = r*planck(M2Vw,T)/median(planck(M2Vw,T))
  43. ext = 0.302*Av/M2Vw**2
  44.  
  45. #our star with the added planck function, normalized
  46. star = M2fn + pln
  47. redd = star/10**(ext/2.5)
  48. reddn = redd/median(redd)
  49.  
  50. #Now we plot:import numpy as np
  51. from pylab import *
  52. from pylab import loadtxt
  53. from scipy import *
  54. import matplotlib.pyplot as plt
  55. import scipy.constants as const
  56. import scipy.stats as stats
  57. from scipy import optimize
  58.  
  59. #Let's plot TTau in microns and
  60. #First column is wavelength, second is flux density, third is something we ignore.
  61. TTmic = loadtxt('/home/kmiskove/Documents/TTau_merged_microns.txt', usecols=(0,1))
  62. TTw, TTf = (TTmic[:,0],TTmic[:,1])
  63. M2V = loadtxt('/home/kmiskove/Documents/M2V_Gl806.txt', skiprows=45, usecols=(0,1))
  64. M2Vw, M2Vf = (M2V[:,0],M2V[:,1])
  65.  
  66. #Let's normalize the fluxes:
  67. TTfn = TTf/median(TTf)
  68. M2fn = M2Vf/median(M2Vf)
  69.  
  70. #This is where we specify our Av, T, and r values.
  71. Av = 0.3 #veiling
  72. T = 1500 #temp
  73. r = 0.7 #scaling
  74. #And let's define some constants:
  75. h = 6.26e-34
  76. c = 3.0e8
  77. k = 1.38e-23
  78. omega = 5.67e-8
  79.  
  80. #The following defines the planck function
  81. def planck(wave, T):
  82.     """
  83.     Input: Wavelength values
  84.     Returns: B_nu(wave,temp), the Planck distribution function in SI units
  85.     """
  86.     numer = 2.*h*(c)**2/(wave*1.e-6)**5
  87.     denom = exp(h*c/(k*T*(wave*1.e-6)))-1.
  88.     return numer/denom
  89.  
  90. pl = planck(M2Vw,T)
  91. pln = r*planck(M2Vw,T)/median(planck(M2Vw,T))
  92. ext = 0.302*Av/M2Vw**2
  93.  
  94. #our star with the added planck function, normalized
  95. star = M2fn + pln
  96. redd = star/10**(ext/2.5)
  97. reddn = redd/median(redd)
  98.  
  99. #Now we plot:
  100. plt.plot(TTw,TTfn)
  101. plt.plot(M2Vw, reddn)
  102. plt.plot(M2Vw, pln)
  103. plt.xlabel('Wavelength (microns)')
  104. plt.ylabel('Flux density')
  105. plt.show(_)
  106.  
  107. #Calculating velocity or smth...c = 3.0e8
  108. yeet = 0.302*.3/2.1660**2
  109. print(yeet)
  110. #^ This is for lambda sub b
  111.  
  112. #Now we do our extinction of 0.0193 = -2.5log(F) and solve for F. I'll just do that in my calculator.
  113. flux = 0.982381
  114. #So we're losing 0.0176 of our light. So we take our flux from IRAF, which was 3.821E-13, and multiply it by 1.0176
  115. #Our new, corrected flux is 3.888E-13! erg/s/cm^2
  116. #When we use figure 3, our  mass accretion rate ends up being E-7, so we'd need E7 years to form the sun. That's 10,000,000 years. Ten million years. That sounds fine.
  117.  
  118. #Now for part 3, we measure our lamba sub b to be 10820 and our lambda not to be 10834.
  119. lb = 10820
  120. ln = 10834
  121. c = 3.0e8
  122.  
  123. v = c*((ln - lb)/ln)
  124. print(v)
  125.  
  126. #So that's 387km/s. The escape velocity of the sun from Earth's distance is about 42km/s, so we can assume this material would escape and be ejected.
  127. plt.plot(TTw,TTfn)
  128. plt.plot(M2Vw, reddn)
  129. plt.plot(M2Vw, pln)
  130. plt.xlabel('Wavelength (microns)')
  131. plt.ylabel('Flux density')
  132. plt.show(_)
  133.  
  134. #Calculating velocity or smth...c = 3.0e8
  135. yeet = 0.302*.3/2.1660**2
  136. print(yeet)
  137. #^ This is for lambda sub b
  138.  
  139. #Now we do our extinction of 0.0193 = -2.5log(F) and solve for F. I'll just do that in my calculator.
  140. flux = 0.982381
  141. #So we're losing 0.0176 of our light. So we take our flux from IRAF, which was 3.821E-13, and multiply it by 1.0176
  142. #Our new, corrected flux is 3.888E-13! erg/s/cm^2
  143. #When we use figure 3, our  mass accretion rate ends up being E-7, so we'd need E7 years to form the sun. That's 10,000,000 years. Ten million years. That sounds fine.
  144.  
  145. #Now for part 3, we measure our lamba sub b to be 10820 and our lambda not to be 10834.
  146. lb = 10820
  147. ln = 10834
  148. c = 3.0e8
  149.  
  150. v = c*((ln - lb)/ln)
  151. print(v)
  152.  
  153. #So that's 387km/s. The escape velocity of the sun from Earth's distance is about 42km/s, so we can assume this material would escape and be ejected.
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top