subtraction/sub_qplot.py at master · jonnybazookatone/subtraction · 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
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
#!/usr/bin/env python
"""
 ================================
| HOTPANTS Data Analysis Pipeline |
|         v1.0                   |
 ================================
| sub_qplot.py |
 ==============

Summary:
        Quick matlotlib representation of the output magnitudes given by sub_apphot.py and sub_mass_apphot.py.

Usage:
        sub_qplot.py --f file --t time --w y
        f: file with magnitudes given from sub_mass_apphot
	t: time array (MJD)
	   if you include a "MJD.txt" it is used to normalise the light curve to GRB t0

        w: write output; y [Default: False]
"""

import sys
import os
import matplotlib
import numpy
#matplotlib.use("WX")
import matplotlib.pyplot as plt
from optparse import OptionParser

__author__ = "Jonny Elliott"
__copyright__ = "Copyright 2012"
__credits__ =  "Felipe Olivares"
__license__ = "GPL"
__version__ = "1.0"
__maintainer__ = "Jonny Elliott"
__email__ = "jonnyelliott@mpe.mpg.de"
__status__ = "Prototype"

def plot(infile, timefile, write):
	"""
	Plotting routine used.
	It is run for a single file, so can be called multiple times.
	File input format:

		OBName Time Mag MagErr
	"""

	# Load time if there
	try:
		ft0 = open(timefile, "r")
		t0 = float(ft0.readlines()[0].replace("\n", ""))
		ft0.close()
		print "t0 file found, using: %f" % (t0)
	except:
		t0 = 0.0
		print "No t0 file found, using default: %f" % (t0)

	# Load and parse
	infile = open(infile, "r")
	inline = infile.readlines()
	infile.close()

	timeArray = []
	timeErrArray = []
	magArray = []
	magErrArray = []
	for i in inline:
		lsplit = i.replace("\n", "").split(" ")

		timeArray.append(float(lsplit[1]) - t0)
		timeErrArray.append(float(lsplit[2]))
		magArray.append(float(lsplit[3]))

		magErrArray.append(float(lsplit[4]))

	# seconds
	timeArray = numpy.array(timeArray) * 60*60*24.
	timeErrArray = numpy.array(timeErrArray)

	# Plot onto graph
	fig = plt.figure(0)
	ax = fig.add_subplot(111)
	ax.errorbar(timeArray, xerr=numpy.array(timeErrArray), y=magArray, yerr=magErrArray, fmt="o")
	ax.invert_yaxis()
	ax.set_xlabel("Time [s]")
	ax.set_ylabel("Brightness [mag]")
	ax.set_xscale("log")
	plt.draw()

	if write:
		outputf = open("lc.dat", "w")
		for i in range(len(timeArray)):
			outputf.write("%f %f %f %f\n" % (timeArray[i], timeErrArray[i], magArray[i], magErrArray[i]))
		outputf.close()

if __name__ == "__main__":

        parser = OptionParser()
        parser.add_option('--f', dest='filelist', help='input aperture file', default=None)
        parser.add_option('--t', dest='time', help='time of grb t0', default=None)
        parser.add_option('--w', dest='write', help='write output', default=False)
        (options, args) = parser.parse_args()

	if options.filelist:
		plot(options.filelist, options.time, options.write)
		plt.savefig("plot.eps")
		plt.show()
	else:
		print __doc__
# Wed Apr 11 17:06:06 CEST 2012