prasopes.py 11.1 KB
Newer Older
Yan's avatar
Yan committed
1
2
3
4
5
6
7
#!/usr/bin/env python3

from matplotlib.backends.backend_qt5agg import FigureCanvasQTAgg as FigureCanvas
from matplotlib.figure import Figure
from matplotlib.widgets import SpanSelector
from PyQt5 import QtCore
from PyQt5 import QtWidgets
8
9
#from rawFin2 import load_raw
from rawParse import load_raw
Yan's avatar
Yan committed
10
11
12
13
14
15
16
import sys
import matplotlib
import numpy as np
import logging
matplotlib.use("Qt5Agg")


17
def zoom_factory(ax, base_scale, plot=None, an=None):
Yan's avatar
Yan committed
18
    """returns zooming functionality to axis"""
Yan's avatar
Yan committed
19
    def zoom_fun(event, chart, ann):
Yan's avatar
Yan committed
20
        """zoom when scrolling"""
Yan's avatar
Yan committed
21
22
        if event.inaxes == ax:
            if event.button == 'up':
23
                #zoom in
Yan's avatar
Yan committed
24
25
                scale_factor = 1/base_scale
            elif event.button == 'down':
26
                #zoom out
Yan's avatar
Yan committed
27
28
                scale_factor = base_scale
            else:
29
                #should not happen
Yan's avatar
Yan committed
30
31
                scale_factor = 1
                print(event.button)
32
33
34

            if event.key == 'shift':
                data = event.ydata
Yan's avatar
Yan committed
35
36
                new_top = data + (ax.get_ylim()[1] - data)*scale_factor
                ax.set_ylim([new_top*-0.01,new_top])
37
38
            else:
                data = event.xdata
Yan's avatar
Yan committed
39
40
41
42
43
44
45
                x_left = data - ax.get_xlim()[0]
                x_right = ax.get_xlim()[1] - data
                ax.set_xlim([data - x_left*scale_factor,
                            data + x_right*scale_factor])
            if ann is not None:
                if type(chart['x']) != type(None):
                    ann_spec(event.inaxes, chart['x'], chart['y'], ann)
46
            ax.figure.canvas.draw()
Yan's avatar
Yan committed
47
48
49

    fig = ax.get_figure()  # get the figure of interest
    # attach the call back
Yan's avatar
Yan committed
50
51
    fig.canvas.mpl_connect('scroll_event',
                           lambda event: zoom_fun(event, plot, an))
Yan's avatar
Yan committed
52
53
54
55

    return zoom_fun


56
def pan_factory(axis, plot=None, annotations=None):
Yan's avatar
Yan committed
57
    """pan spectrum when you press a button"""
58
    def pan_fun(event, ax):
Yan's avatar
Yan committed
59
        #rescale to origin if doubleclicked
Yan's avatar
Yan committed
60
61
62
63
        if event.dblclick and event.inaxes == ax:
            ax.get_figure()
            ax.autoscale(True)
            ax.set_ylim(ax.get_ylim()[1]*-0.01,ax.get_ylim()[1]*1.1)
64
65
            if annotations is not None:
                ann_spec(ax, plot['x'], plot['y'], annotations)
Yan's avatar
Yan committed
66
            ax.figure.canvas.draw()
Yan's avatar
Yan committed
67
        #otherwise pan
Yan's avatar
Yan committed
68
        elif event.button == 1 and event.inaxes == ax:
Yan's avatar
Yan committed
69
            ax.start_pan(event.x, event.y, event.button)
Yan's avatar
Yan committed
70
71
72
            id_drag = fig.canvas.mpl_connect(
                'motion_notify_event',
                lambda action: drag_fun(action, ax))
Yan's avatar
Yan committed
73
            id_release = fig.canvas.mpl_connect(
Yan's avatar
Yan committed
74
75
76
                'button_release_event',
                lambda action: drag_end(
                    action, id_drag, id_release, plot, annotations, ax))
Yan's avatar
Yan committed
77

78
    def drag_fun(event, ax):
Yan's avatar
Yan committed
79
80
81
        ax.drag_pan(1, 'x', event.x, event.y)
        ax.figure.canvas.draw()

Yan's avatar
Yan committed
82
    def drag_end(event, id_drag, id_release, chart, ann, ax):
Yan's avatar
Yan committed
83
84
85
        if event.button == 1:
            fig.canvas.mpl_disconnect(id_drag)
            fig.canvas.mpl_disconnect(id_release)
Yan's avatar
Yan committed
86
87
88
            if ann is not None:
                if type(chart['x']) != type(None):
                    ann_spec(ax, chart['x'], chart['y'], ann)
89
            ax.figure.canvas.draw()
Yan's avatar
Yan committed
90
    fig = axis.get_figure()
Yan's avatar
Yan committed
91
92
    fig.canvas.mpl_connect('button_press_event',
                           lambda action: pan_fun(action, axis))
Yan's avatar
Yan committed
93

Yan's avatar
Yan committed
94
95
def pick_masses(x_min, x_max, ax, annotations, plot):
    """zoom the spectrum in x axis by mass range"""
Yan's avatar
Yan committed
96
    ax.set_xlim(x_min, x_max)
Yan's avatar
Yan committed
97
    ann_spec(spectrum, plot['x'], plot['y'], annotations)
Yan's avatar
Yan committed
98
    ax.figure.canvas.draw()
Yan's avatar
Yan committed
99

Yan's avatar
Yan committed
100
101
102
def pick_times(x_min, x_max, mpl_spectrum, data_set, mpl_chromatogram,
               annotations, mass_spect):
    """plot averaged spectrum of subselected part of the chromatogram"""
Yan's avatar
Yan committed
103
104
    start_scan = None
    end_scan = None 
Yan's avatar
Yan committed
105
106
    for i, times in enumerate(data_set['chrom_dat'][0, :]):
        if times > x_min and start_scan is None:
Yan's avatar
Yan committed
107
            start_scan = i
Yan's avatar
Yan committed
108
        if times > x_max and end_scan is None:
Yan's avatar
Yan committed
109
            end_scan = i
Yan's avatar
Yan committed
110
    mpl_spectrum.clear()
Yan's avatar
Yan committed
111
112
113
114
115
    annotations.clear()
    mass_spect['x'] = data_set['masses']
    mass_spect['y'] = np.mean(data_set['matrix'][start_scan:end_scan]
                              , axis=0)
    spectrum_plot(spectrum, mass_spect['x'], mass_spect['y'], annotations)
Yan's avatar
Yan committed
116
    mpl_chromatogram.clear()
Yan's avatar
Yan committed
117
118
    chrom_plot(mpl_chromatogram, data_set['chrom_dat'][0, :],
               data_set['chrom_dat'][1, :])
119
120
121
    mpl_chromatogram.plot(data_set['chrom_dat'][0, start_scan:end_scan],
                          data_set['chrom_dat'][1, start_scan:end_scan],
                          'b.')
Yan's avatar
Yan committed
122

Yan's avatar
Yan committed
123
124
125
def populate(mpl_chromatogram, mpl_spectrum, data_set, time_sel, an,
             mass_spect):
    """populate the GUI plots with desired dataset"""
126
    an.clear()
Yan's avatar
Yan committed
127
    def update_span_selector(chrom, spec, d_set, ms_spec):
Yan's avatar
Yan committed
128
129
130
131
132
133
        time_selector = SpanSelector(
            chrom, lambda x_min, x_max: pick_times(
                       x_min, x_max,spec, d_set, chrom, an, ms_spec),
            'horizontal',
            useblit=True,
            rectprops=dict(alpha=0.15, facecolor='purple'), button=3)
Yan's avatar
Yan committed
134
135
136
        return time_selector
    mpl_spectrum.clear()
    mpl_chromatogram.clear()
Yan's avatar
Yan committed
137
138
139
140
141
142
143
    mass_spect['x'] = data_set['masses']
    mass_spect['y'] = np.mean(data_set['matrix'], axis=0)
    spectrum_plot(mpl_spectrum, mass_spect['x'], mass_spect['y'], an)
    chrom_plot(mpl_chromatogram, data_set['chrom_dat'][0, :],
               data_set['chrom_dat'][1, :])
    time_sel[0] = update_span_selector(mpl_chromatogram, mpl_spectrum,
                                       data_set, mass_spect)
Yan's avatar
Yan committed
144
145
146
    mpl_chromatogram.figure.canvas.draw()
    mpl_spectrum.figure.canvas.draw()

Yan's avatar
Yan committed
147
148
149
150
151
152
153
154
155
156
157
158
159
160
def ann_spec(ms_spec, mass, intensities, annotations, ann_limit=0.01):
    """annotate spectrum
    
    First define borders where the annotation should occur. Then the
    function finds local maxima of the data set by forcycling through
    the dataset and then substracts found local minimas so it contains
    only subset of the most intense one which are in reasonable distance
    from each other. Then the function remove old annotation and
    annotate the spectrum by abovementioned subsection of local maximas.
    """
    def sub_peaks(peakz, x_value, y_value, xrange, yrange, coef_x=10,
                  coef_y=10):
        """Returns reasonable subselection of local maximas"""
        sort_peaks = sorted(peakz,
161
162
                key=lambda h: y_value[h],
                reverse=True)
Yan's avatar
Yan committed
163
164
165
166
        subtracted_peaks=[]
        red_x = xrange / coef_x
        red_y = yrange / coef_y
        for peak in sort_peaks:
167
            add = True
Yan's avatar
Yan committed
168
169
170
171
172
            yv = y_value[peak]
            xv = x_value[peak]
            for u in subtracted_peaks:
                if abs(yv - y_value[u]) < red_y \
                        and abs(xv - x_value[u]) < red_x:
173
174
                    add = False
            if add:
Yan's avatar
Yan committed
175
176
                subtracted_peaks.append(peak)
        return subtracted_peaks
177

Yan's avatar
Yan committed
178
    borders=ms_spec.get_xlim()
Yan's avatar
Yan committed
179
180
181
182
    start = np.argmax(mass>borders[0])
    end = np.argmax(mass>borders[1])
    if end == 0:
        end = np.argmax(mass)
183
184
    peaks = []
    l = ms_spec.get_ylim()[1] * ann_limit
Yan's avatar
Yan committed
185
186
    for i,intensity in enumerate(intensities[start:end], start=start):
        #end selects the first element which not to include.
187
        #Thus q[i+1] always exists as it points to that last element.
Yan's avatar
Yan committed
188
        if intensity > l and intensity > q[i-1] and intensity > q[i+1]:
189
            peaks.append(i)
Yan's avatar
Yan committed
190
    s_peaks = sub_peaks(peaks, mass, intensities,
191
192
                        ms_spec.get_xlim()[1] - ms_spec.get_xlim()[0],
                        ms_spec.get_ylim()[1] - ms_spec.get_ylim()[0])
193
    
Yan's avatar
Yan committed
194
195
    for intensity in annotations:
        intensity.remove()
196
    annotations.clear()
Yan's avatar
Yan committed
197

Yan's avatar
Yan committed
198
    for i in s_peaks:
Yan's avatar
Yan committed
199
200
201
        annotations.append(ms_spec.annotate(
            '{0:.2f}'.format(mass[i]), xy=(mass[i], intensities[i]),
            textcoords='data'))
Yan's avatar
Yan committed
202

Yan's avatar
Yan committed
203
204
def spectrum_plot(spc, mass, intensity,an):
    """Define and populate spectrum"""
Yan's avatar
Yan committed
205
206
207
208
209
210
211
212
213
    spc.plot(mass, intensity)
    spc.set_title("Spectrum:", loc ="right")
    spc.set_xlabel("m/z")
    spc.set_ylabel("ion count")
    spc.set_ylim(spc.get_ylim()[1] * -0.01, )
    spc.ticklabel_format(scilimits=(0, 0), axis='y')
    ann_spec(spc, mass, intensity, an)

def chrom_plot(chrom, times, tic):
Yan's avatar
Yan committed
214
    """Define and populate chromatogram"""
Yan's avatar
Yan committed
215
    chrom.plot(times, tic)
Yan's avatar
Yan committed
216
217
    chrom.set_ylim(chrom.get_ylim()[1] * -0.011,
                   chrom.get_ylim()[1] * 1.1)
Yan's avatar
Yan committed
218
219
220
221
222
223
    chrom.set_title("Chromatogram:", loc ="right")
    chrom.set_xlabel("time (min)")
    chrom.set_ylabel("total ion count")
    chrom.ticklabel_format(scilimits=(0, 0), axis='y')

def open_file(chrom, spc, d_set, time_sel, an, ms_s):
Yan's avatar
Yan committed
224
    """Returns chromatogram, masses and ion intensities"""
Yan's avatar
Yan committed
225
226
227
    filename=QtWidgets.QFileDialog.getOpenFileName(
            caption = "Open spectrum",
            filter="Finnigan RAW files (*.raw, *.RAW)")[0]
Yan's avatar
Yan committed
228
229
    if filename is '':
        return
Yan's avatar
Yan committed
230
231
    d_set['chrom_dat'], d_set['masses'], d_set['matrix'] \
        = load_raw(filename)
Yan's avatar
Yan committed
232
    populate(chrom, spc, d_set, time_sel, an, ms_s)
Yan's avatar
Yan committed
233
234
235
236


if __name__=="__main__":
    #ds for data_set
Yan's avatar
Yan committed
237
    ds = dict(chrom_dat=None,masses=None,matrix=None)
238
239
    #mass spectrometry set
    ms = dict(x=None,y=None)
Yan's avatar
Yan committed
240

Yan's avatar
Yan committed
241
    p_logger = logging.getLogger('parseLogger')
Yan's avatar
Yan committed
242
    logging.basicConfig()
Yan's avatar
Yan committed
243
    p_logger.setLevel("WARN")
Yan's avatar
Yan committed
244
245
246
247
248
249
250
251
252
253
254

    graph = Figure(figsize=(5,4),dpi=100)
    graph.patch.set_facecolor("None")

    chromatogram = graph.add_subplot(211,facecolor=(1,1,1,0.8))
    spectrum = graph.add_subplot(212,facecolor=(1,1,1,0.8))
    graph.tight_layout()

    mpl_canvas = FigureCanvas(graph)
    mpl_canvas.setStyleSheet("background-color:transparent;")
    mpl_canvas.setAutoFillBackground(False)
255
    mpl_canvas.setFocusPolicy( QtCore.Qt.ClickFocus )
Yan's avatar
Yan committed
256
257
258
259
    
    timeSelector=[None]
    annotation=[]

260
    pan_factory(chromatogram)
Yan's avatar
Yan committed
261
    zoom_factory(chromatogram, 1.15)
262
263
    pan_factory(spectrum, ms, annotation)
    zoom_factory(spectrum, 1.15, ms, annotation)
Yan's avatar
Yan committed
264
    mass_selector = SpanSelector(
Yan's avatar
Yan committed
265
266
267
268
        spectrum, lambda x_min, x_max: pick_masses(
            x_min, x_max, spectrum, annotation, ms),
        'horizontal', useblit=True, rectprops=dict(alpha=0.15,
                                                   facecolor='purple'),
Yan's avatar
Yan committed
269
        button=3)
Yan's avatar
Yan committed
270
271
272
273

    app = QtWidgets.QApplication(sys.argv)
    main_window = QtWidgets.QMainWindow()

Yan's avatar
Yan committed
274
    file_menu = QtWidgets.QMenu('&File', main_window)
Yan's avatar
Yan committed
275
    main_window.menuBar().addMenu(file_menu)
Yan's avatar
Yan committed
276
    file_menu.addAction('&Open..', 
Yan's avatar
Yan committed
277
278
        lambda: open_file(chromatogram, spectrum, ds, timeSelector,
                          annotation, ms),
Yan's avatar
Yan committed
279
        QtCore.Qt.CTRL + QtCore.Qt.Key_O)
Yan's avatar
Yan committed
280
281
    file_menu.addAction('&Quit', main_window.close,
                        QtCore.Qt.CTRL + QtCore.Qt.Key_Q)
Yan's avatar
Yan committed
282
283
284
285
286
287

    main_widget = QtWidgets.QWidget(main_window)
    main_window.setCentralWidget(main_widget)

    layout = QtWidgets.QVBoxLayout(main_widget)
    layout.addWidget(mpl_canvas)
288
    mpl_canvas.setFocus()
Yan's avatar
Yan committed
289

Yan's avatar
Yan committed
290
291
292
293

    if len(sys.argv) == 2:
        raw_file=sys.argv[1]
        ds['chrom_dat'],ds['masses'],ds['matrix'] = load_raw(raw_file)
Yan's avatar
Yan committed
294
295
        populate(chromatogram, spectrum, ds, timeSelector, annotation,
                 ms)
Yan's avatar
Yan committed
296
297
298
    else:
        spectrum_plot(spectrum, [0], [0], annotation)
        chrom_plot(chromatogram, [0], [0])
Yan's avatar
Yan committed
299
300
301
302

    main_window.show()
    sys.exit(app.exec_())