Make your own NOAA / GOES familiar X-ray flux 3days plot in Python

Make your own NOAA / GOES familiar X-ray flux 3days plot in Python

** Note: It contains fairly specialized content, such as the data used. I will leave the explanation of physics and terms to the appropriate literature. For those who came to see from the search etc., the explanation of some modules is given from [4. Explanation of each part of the code](# 4-Explanation of each part of the code), so please refer to that. ** ** I know people in the same industry reinventing the wheel, but if there is a more innovative method, please teach me. However, the departure from IDL to Python is one major premise.

1. Background

1.1. What and why do you bother to make it yourself

I will make the following figure by myself, which I should have seen somewhere if I was doing research related to the sun. ** The figure you made is in [5. Actually draw](# 5-Actually draw), so please take a look at it alone. ** ** goes_xrays_20170906.png

Recently ~~ (I haven't seen it recently so I don't know when) ~~, GOES Xray flux page is interactive It changed to a graph. This is not a problem, but ――It is true that the value is easy to see, but honestly there is no opportunity to need real-time information on what time and minute the flux was. --Since X-ray time profiles in this format have been used for a long time, many people are relieved to see images in this format. --Past 3days plot is Page above or Solar Monitor You can get it from, but if you paste it directly on a slide etc., you will get angry because the line is thin --The operation of GOES-15 and 14 will end with the operation of GOES-16, so I don't know if the same format diagram will be updated in the future. There are many things to think about, so the purpose of this post is that if you study matplotlib, you will be able to make it yourself.

1.2. Prerequisite knowledge

In this article, in this article 2Fqiita-image-store.s3.amazonaws.com%2F0%2F222849%2F1e675340-79ff-4aa1-1891-b79d3c369eb3.png?ixlib=rb-1.2.2&auto=format&gif-q=60&q=75&s=3ea10339ee8bdcf6cb7feba540112e) I often use the word ʻaxesthat appears in. If your goal is to understand the following content and code, we recommend that you read theFigure and ʻAxes with a quick understanding.

In addition, as described in the preface, technical terms are not explained. Please note. I will explain the modules etc. without technical terms as much as possible.

2. Creating a program

2.1. Module import

This time, I imported various modules that I use for the first time in order to stick to the appearance.

import numpy as np
import datetime,calendar,requests
import pandas as pd
from io import StringIO

import matplotlib
%matplotlib inline ##For drawing on jupyter notebook

import matplotlib.pyplot as plt
import matplotlib.ticker as ticker
import matplotlib.transforms as transforms
from matplotlib.dates import DateFormatter
from matplotlib.dates import date2num
from matplotlib.offsetbox import AnchoredOffsetbox, TextArea, HPacker, VPacker

2.2. Where to put the data

Basically, all the data up to GOES 15 is open to the public. For example, in September 2017, GOES-15 X-ray 1-minute average data is in csv format [here](https://satdat.ngdc.noaa.gov/sem/goes/data/avg/2017/09/goes15 It is located at /csv/g15_xrs_1m_20170901_20170930.csv). This time, we will proceed so that you can create a graph without downloading the file locally.

At this point, we will assume that you know the GOES primary and GOES secondary numbers when you want to create the graph. I think that the IDL command contained primary information, so I will add it someday if possible.

To create a 3days plot, first reference the original data. The 1-minute average X-ray data of GOES is published in csv format every month, so be careful when referring to the data on a monthly basis.

def read_csv_skiprows(url):
    resp = requests.get(url)
    data_text = resp.text
    
    is_skip, skip_row = False,1
    for line in data_text.split('\n'):##Divide each line'data:'search for
        if 'data:' in line:
            is_skip = True
            break
        skip_row += 1

    if not is_skip:
        skip_row = 0
    return pd.read_csv(StringIO(data_text),skiprows=skip_row)

def get_goes_csv(goesnum,date,trange=3):
    ##First date on the plot:tstart
    if type(date)==list:
        ts = datetime.date(date[0],date[1],date[2])
    elif (type(date)==datetime.date)or(type(date)==datetime.datetime):
        ts = date
    ##Last date on the plot:tend
    te = ts+datetime.timedelta(days=trange)
    url = []
    if ts.month!=te.month: ##If the graph spans months
        for i,t in enumerate([ts,te]): ##Get index with enumerate
            ##Get the last day of each month
            last_day = calendar.monthrange(t.year,t.month)[1]
            url.append('https://satdat.ngdc.noaa.gov/sem/goes/data/avg/{:04}/{:02}/goes{:02}/csv/g{:02}_xrs_1m_{:04}{:02}{:02}_{:02}{:02}{:02}.csv'.format(t.year,t.month,goesnum,goesnum,t.year,t.month,1,t.year,t.month,last_day))
            if i==0:
                d = read_csv_skiprows(url[i])
            else: ##This time, it is assumed that the month will not be crossed more than once.
                df = pd.concat([d,read_csv_skiprows(url[i])])
            
    else:
        last_day = calendar.monthrange(ts.year,ts.month)[1]
        url.append('https://satdat.ngdc.noaa.gov/sem/goes/data/avg/{:04}/{:02}/goes{:02}/csv/g{:02}_xrs_1m_{:04}{:02}{:02}_{:02}{:02}{:02}.csv'\
                    .format(ts.year,ts.month,goesnum,goesnum,ts.year,ts.month,1,ts.year,ts.month,last_day))
        df = read_csv_skiprows(url[0])
    
    ##Since the date column is str type, change it to a type that can process time
    df.time_tag = pd.to_datetime(df.time_tag)
    ##Returns data between ts and te
    df_obj = df[(df.time_tag>=ts)&(df.time_tag<=te)]
    return df_obj

Although it is said to be in csv format, it can be a bit annoying, so here As you can see, the data you want starts from the bottom, and the first one contains a description of the data. Read_csv of pandas can jump to the target line by specifying the argument of skiprows or header when there is such a description, but what is troublesome is the file (or observation The number of lines to reach the data varies depending on the time).
However, since there is a description 'data:' just before the data, find the number of lines to skip with this as the target. This method was previously taught at teratail.

With the above processing, GOES X-ray data can be acquired. If you want to make your own GOES X-ray plot for the time being, add GOES data to Sunpy [^ 1] Since there is a command to get, you can also use that. [^ 1]: This part dares to refer to the old version of the document. Because, for some reason, the figure is not displayed in the latest version document and it is intuitive to understand. Because it's difficult.

However, the drawback is that you can only get data from one satellite at a time (although it may be an overstatement as it is not a problem in most cases). For example, in the flare on September 6, 2017, the data of flare that occurred at 8:57 is just It has fallen in GOES-15, and it is necessary to supplement it with the data of GOES-13. Considering this case, it may be better to go to NOAA.

3. Completed code

From the conclusion, I was able to output with the following code. It's very long. Partially explained in [4. Explanation of each part of the code](# 4-Explanation of each part of the code). If you want a graph for the time being, skip to [5. Actually draw](# 5-Actually draw).

plt.rcParams['font.family'] = 'Simplex'
plt.rcParams['font.size'] = 10
plt.rcParams['figure.dpi'] = 100

plt.rcParams['axes.linewidth'] = 1

plt.rcParams['lines.antialiased'] = False
plt.rcParams['text.antialiased'] = False
plot_antialiased = False

##All in pixel units
wdict = dict(
    grid  = 1, ##Grid line width
    spine = 1, ##Frame line width
    tick  = 1, ##Scale line width
    line  = 0.01, ##Line width of the graph
)
##Convert pixel to point
def px2pt(px):
    dpi = plt.rcParams['figure.dpi']
    return px*72/dpi

def plot_goes(t, trange ,goes ,savename, wdict, width=6.4, ylim=(10**-9,10**-2)):
    '''
    t      :Plot start date(type=datetime.date or list)
             (e.g. datetime.date(2017,9,4) or [2017,9,4]) 
    trange :Period to plot Daily(type=int) (e.g. trange=3 > 3 days plot)
    goes   :GOES satellite number list format to acquire(e.g. [15,13])
             index =0 is primary(1-8A>red, 0.5-4A>blue)
             index =1 is secondary(1-8A>yellow, 0.5-4A>purple)
    width  :Plot width default=6.4
             if it is 'extend', extend the width with trange
    ylim   :y-axis range default=(10**-9,10**-2)
Probably not necessary to change, but designed to be okay to change
    '''
    if type(t)==list:
        t = datetime.date(t[0],t[1],t[2])
    elif (type(t)==datetime.date)or(type(t)==datetime.datetime):
        pass
    else:
        print('t should be a list or a datetime object.')
    
    fs_def = plt.rcParams['font.size']
    dpi    = plt.rcParams['figure.dpi']
    y0_fig = plt.rcParams['figure.subplot.bottom']
    y1_fig = plt.rcParams['figure.subplot.top']
    
    ##Width setting
    if width=='extend':
        width = (dpi*6.4 + (dpi*6.4 - fs_def*10)/3*(trange-3))/dpi
    elif type(width)==float:
        width = width
    else:
        print('width is \'extend\' or float')

    ##Left and right margins are font size* 10 (Percentage to figure)
    margin = fs_def/(100*width)*10

    ##Width per day(Percentage to figure)
    w_day  = (1-2*margin)/trange
    ax = []

    ##frame(axes)Creation
    fig = plt.figure(figsize=(width,4.8))
    for i in range(trange):
        t_s = t + datetime.timedelta(days=i)
        t_e = t_s + datetime.timedelta(days=1)

        if i==0: ##First day
            ax.append(fig.add_axes([margin+i*w_day,y0_fig,w_day,(y1_fig-y0_fig)],\
                                   yscale='log', xlim=[t_s,t_e], ylim=ylim, zorder=-i))
            new_xticks = date2num([t_s,t_e])

        else: ##After the second day
            ax.append(fig.add_axes([margin+i*w_day,y0_fig,w_day,(y1_fig-y0_fig)],\
                                   yscale='log', xlim=[t_s,t_e], ylim=ylim, yticklabels=[], zorder=-i))
            new_xticks = date2num([t_e])

        ax[i].xaxis.set_major_locator(ticker.FixedLocator(new_xticks))
        ax[i].xaxis.set_major_formatter(DateFormatter('%b %d'))
        ax[i].xaxis.set_minor_locator(matplotlib.dates.HourLocator(interval=3))
        ax[i].yaxis.grid(lw=px2pt(wdict['grid']),color='black')

        ax[i].tick_params(axis='x', which='major', bottom=True, top=True, direction='in', length=6, width=px2pt(wdict['tick']), pad=fs_def)
        ax[i].tick_params(axis='x', which='minor', bottom=True, top=True, direction='in', length=3, width=px2pt(wdict['tick']))

        for j in ax[i].spines.keys(): ##Adjusting the line width of the frame
            ax[i].spines[j].set_linewidth(px2pt(wdict['spine']))
        ##Adjusting the scale
        if i==0: ##First day
            ax[i].tick_params(axis='y', which='both', right=True, left=True, direction='in', length=6, width=px2pt(wdict['tick']))
            ax[i].spines['right'].set_visible(False)
        elif i==trange-1: ##Last day
            ax[i].tick_params(axis='y', which='both', right=True, left=False, direction='in', length=6, width=px2pt(wdict['tick']))
            ax[i].spines['left'].set_visible(False)            
        else: ##Other
            ax[i].tick_params(axis='y', which='both', right=True, left=False, direction='in', length=6, width=px2pt(wdict['tick']))
            ax[i].spines['left'].set_visible(False)
            ax[i].spines['right'].set_visible(False)


    ## Flare class text
    f_pos = ax[-1].get_position()

    trans = transforms.blended_transform_factory(fig.transFigure,transforms.IdentityTransform())
    X10_y_display   = ax[-1].transData.transform((date2num(t_e),10**-3))[1]
    X_y_display     = ax[-1].transData.transform((date2num(t_e),10**-4))[1]
    w = X10_y_display - X_y_display

    ##Relative size of font size to the entire figure
    ##Used to allow margins when placing class text so that it does not stick to the outer frame
    cls_text_size = plt.rcParams['font.size']/(plt.rcParams['figure.figsize'][0])/72

    for i,cls in enumerate(['X', 'M', 'C', 'B', 'A']):
        pos = [f_pos.x1,X10_y_display - w/2 - i*w]
        if (pos[1]>=ax[-1].transData.transform((date2num(t_e),ax[-1].get_ylim()[0]))[1])&\
           (pos[1]<=ax[-1].transData.transform((date2num(t_e),ax[-1].get_ylim()[1]))[1]):
            ax[-1].text(pos[0]+cls_text_size/5,pos[1],cls,horizontalalignment='left',verticalalignment='center',transform=trans)

    ## GOES number and wavelength
    ybox_b = TextArea('GOES{} 0.5-4.0 A'.format(goes[0]), textprops=dict(color='#001DFF', size='large' ,rotation=90,ha='left',va='bottom')) #Blue
    ybox_r = TextArea('GOES{} 1.0-8.0 A'.format(goes[0]), textprops=dict(color='#FF0000', size='large' ,rotation=90,ha='left',va='bottom')) #Red
    ybox_primary = VPacker(children=[ybox_r, ybox_b],align='bottom', pad=0, sep=2*fs_def)
    if len(goes)==2:
        ybox_p = TextArea('GOES{} 0.5-4.0 A'.format(goes[1]), textprops=dict(color='#5400A8', size='large' ,rotation=90,ha='left',va='bottom')) #purple
        ybox_y = TextArea('GOES{} 1.0-8.0 A'.format(goes[1]), textprops=dict(color='#FFA500', size='large' ,rotation=90,ha='left',va='bottom')) #yellow
        ybox_secondary = VPacker(children=[ybox_y, ybox_p],align='bottom', pad=0, sep=2*fs_def)
        ybox = HPacker(children=[ybox_primary,ybox_secondary],align='bottom', pad=0, sep=fs_def)
        anchored_ybox = AnchoredOffsetbox(loc=5, child=ybox, pad=0., frameon=False, bbox_to_anchor=(1-margin*0.25, 0.5), 
                                          bbox_transform=fig.transFigure, borderpad=0.)
    else:
        ybox = ybox_primary
        anchored_ybox = AnchoredOffsetbox(loc=5, child=ybox, pad=0., frameon=False, bbox_to_anchor=(1-margin*0.5, 0.5), 
                                          bbox_transform=fig.transFigure, borderpad=0.)
    ax[-1].add_artist(anchored_ybox)

    ##Axis label creation(dummy axes)
    ax_dummy = fig.add_axes([margin,y0_fig,1-margin*2,y1_fig-y0_fig],frameon=False)
    ax_dummy.tick_params(labelcolor='none', top='off', bottom='off', left='off', right='off')
    ax_dummy.set_xlabel('Universal Time',size='large')
    ax_dummy.set_ylabel('Watts m$^{-2}$',size='large')
    ax_dummy.set_title('GOES Xray Flux (1-minute data)',size='x-large')

    ax_dummy.text(1, 1-cls_text_size/5, 'Begin: {} {} {} 0000 UTC'.format(t.year,t.strftime('%b'),t.day),
            horizontalalignment='right',
            verticalalignment='top',
            transform=fig.transFigure)
    ## X-ray data plot
    color = iter(['#FF0000','#001DFF','#FFA500','#5400A8']) ## red, blue, yellow, purple 
    zo = iter([4,3,2,1])
    for g in goes:
        wavelength = iter(['1.0-8.0 A','0.5-4.0 A'])
        data = get_goes_csv(g,t,trange)
        xray_long  = np.ma.masked_where((data.B_AVG<=ylim[0])|(data.B_AVG>=ylim[1]),data.B_AVG)
        xray_short = np.ma.masked_where((data.A_AVG<=ylim[0])|(data.A_AVG>=ylim[1]),data.A_AVG)
        ax[0].plot(data.time_tag,xray_long, lw=px2pt(wdict['line']),c=next(color),label='GOES{} 1.0-8.0 A'.format(g),aa=plot_antialiased,clip_on=False,zorder=next(zo))
        ax[0].plot(data.time_tag,xray_short,lw=px2pt(wdict['line']),c=next(color),label='GOES{} 0.5-4.0 A'.format(g),aa=plot_antialiased,clip_on=False,zorder=next(zo))
        
    ## save graph
    fig.savefig(savename)

4. Explanation of each part of the code

I will introduce each part. If you want to be able to use it for the time being, please skip it.

4.1. Fonts

According to here, the font used in IDL is Simplex Roman. This time, I was particular about the appearance, so I downloaded it from Free Font Site and used it. It was not installed as standard in Mac OS, so if you want to stick to it, please download it. If you don't need it

plt.rcParams['font.family'] = 'Arial'

I think that the graph will be easy to see if you set it as such.

4.2. Aliasing

plt.rcParams['lines.antialiased'] = False
plt.rcParams['text.antialiased'] = False
plot_antialiased = False

When matplotlib outputs in bitmap format such as png or jpg, it automatically complements and anti-aliases so that there is no jaggedness. If you draw a diagram normally, you don't have to worry about this, but this time I was particular about it. If you set all the above parameters to True, antialiasing will be performed, and I think that is safer. In addition, it seems that there is no way to cancel tick antialiasing.

4.3. Pixel and point conversion

def px2pt(px):
    dpi = plt.rcParams['figure.dpi']
    return px*72/dpi

This is also a wasteful commitment. In the original image, the width of every line was specified as 1 pixel, so I created it for that. Because, the unit when specifying the width etc. with matplotlib is all point. 1 point is 1/72 of 1 inch, and the relationship between inch and pixel changes depending on the dpi, so that is corrected.

4.4. Width setting

if width=='extend':
    width = (dpi*6.4 + (dpi*6.4 - fs_def*10)/3*(trange-3))/dpi
elif type(width)==float:
    width = width
else:
    print('width is \'extend\' or float')

The original image had 480 pixels in height and 640 pixels in width, so the basics are that (dpi is 100). However, if you want to create a plot longer than 3 days, you can specify width =='extend' to increase the horizontal size while keeping the width per day the same as for the 3days plot. Did. Also, if you normally specify the width with the float type, it is set to that size, but operation cannot be guaranteed. If you think about it carefully, if you create a plot shorter than 3 days, an error will be thrown, so I will describe it later ([6. Draw simply](# 6-Draw simply)).

4.5. Margin settings

##Left and right margins are font size* 10 (Percentage to figure)
margin = fs_def/(100*width)*10

##Width per day(Percentage to figure)
w_day  = (1-2*margin)/trange

Set the left and right margins. As you can see in the comment out, I set it to about 10 times the font size.

4.6. Frame setting

##frame(axes)Creation
fig = plt.figure(figsize=(width,4.8))
for i in range(trange):
    t_s = t + datetime.timedelta(days=i)
    t_e = t_s + datetime.timedelta(days=1)

    if i==0: ##First day
        ax.append(fig.add_axes([margin+i*w_day,y0_fig,w_day,(y1_fig-y0_fig)],\
                               yscale='log', xlim=[t_s,t_e], ylim=ylim, zorder=-i))
        new_xticks = date2num([t_s,t_e])

    else: ##After the second day
        ax.append(fig.add_axes([margin+i*w_day,y0_fig,w_day,(y1_fig-y0_fig)],\
                               yscale='log', xlim=[t_s,t_e], ylim=ylim, yticklabels=[], zorder=-i))
        new_xticks = date2num([t_e])

    ax[i].xaxis.set_major_locator(ticker.FixedLocator(new_xticks))
    ax[i].xaxis.set_major_formatter(DateFormatter('%b %d'))
    ax[i].xaxis.set_minor_locator(matplotlib.dates.HourLocator(interval=3))
    ax[i].yaxis.grid(lw=px2pt(grid_width),color='black')

    ax[i].tick_params(axis='x', which='major', bottom=True, top=True, direction='in', length=6, width=px2pt(tick_width), pad=fs_def)
    ax[i].tick_params(axis='x', which='minor', bottom=True, top=True, direction='in', length=3, width=px2pt(tick_width))

    for j in ax[i].spines.keys(): ##Adjusting the line width of the frame
        ax[i].spines[j].set_linewidth(px2pt(spine_width))
    ##Adjusting the scale
    if i==0: ##First day
        ax[i].tick_params(axis='y', which='both', right=True, left=True, direction='in', length=6, width=px2pt(tick_width))
        ax[i].spines['right'].set_visible(False)
    elif i==trange-1: ##Last day
        ax[i].tick_params(axis='y', which='both', right=True, left=False, direction='in', length=6, width=px2pt(tick_width))
        ax[i].spines['left'].set_visible(False)            
    else: ##Other
        ax[i].tick_params(axis='y', which='both', right=True, left=False, direction='in', length=6, width=px2pt(tick_width))
        ax[i].spines['left'].set_visible(False)
        ax[i].spines['right'].set_visible(False)

Probably the most time-consuming part of the finest details. After all, it seems that matplotlib doesn't have the memory inside the original graph, the means to draw it. inside_ticks.png

I thought about forcibly writing a horizontal bar that separated the range with ʻax.hline () , but I didn't want to do much brute force. So, when I asked about synchronization that can handle IDL, it seems that IDL has a function to add an axis anywhere. ~~ Give it to matplotlib. ~~ I couldn't help it, so I decided to deal with it by splitting ʻaxes. Specifically, by specifying the coordinates with ʻadd_axes`, the area is divided as shown in the figure below, and by changing the display of the border (it seems to be called spine in matplotlib), the frame can be created safely. Has completed. multiple_axes.png

4.6. Flare class text

## Flare class text
##Get the coordinates of the leftmost axes
f_pos = ax[-1].get_position()

trans = transforms.blended_transform_factory(fig.transFigure,transforms.IdentityTransform())
X10_y_display   = ax[-1].transData.transform((date2num(t_e),10**-3))[1]
X_y_display     = ax[-1].transData.transform((date2num(t_e),10**-4))[1]
w = X10_y_display - X_y_display

##Relative size of font size to the entire figure
##Used to allow margins when placing class text so that it does not stick to the outer frame
cls_text_size = plt.rcParams['font.size']/(plt.rcParams['figure.figsize'][0])/72

for i,cls in enumerate(['X', 'M', 'C', 'B', 'A']):
    pos = [f_pos.x1,X10_y_display - w/2 - i*w]
    if (pos[1]>=ax[-1].transData.transform((date2num(t_e),ax[-1].get_ylim()[0]))[1])&\
       (pos[1]<=ax[-1].transData.transform((date2num(t_e),ax[-1].get_ylim()[1]))[1]):
        ax[-1].text(pos[0]+cls_text_size/5,pos[1],cls,horizontalalignment='left',verticalalignment='center',transform=trans)

This is the least readable part of this code. I'm sorry.
The flare class is determined by the value on the vertical axis, so you need to put a letter in ʻaxes.text at the correct value. There are several ways to set the coordinates when placing ʻaxes.text,

  1. Coordinates corresponding to the pixels of the entire figure
  2. Relative coordinates with (0,0) in the lower left and (1,1) in the upper right.
  3. In the axes (frame) in the figure, the lower left is (0,0) and the upper right is (1,1).
  4. Coordinates corresponding to the values that the figure has

However, since the horizontal axis of this figure is the time, I gave up (4) (is it possible using date2num?). Therefore, I designed it so that characters can be placed in the correct position outside the border while using blended_transform_factory which can handle a mixture of (1) to (3). In particular,

trans = transforms.blended_transform_factory(fig.transFigure,transforms.IdentityTransform())

By doing so, the x coordinate refers to (2; figure.transFigure) and the y coordinate refers to (4;transforms.IdentityTransform ()) so that text can be placed. There seems to be a mess of improvement here ...

4.7. Creating a Legend

## GOES number and wavelength

##Set by rotating the text of each color 90 degrees
ybox_b = TextArea('GOES{} 0.5-4.0 A'.format(goes[0]), textprops=dict(color='#001DFF', size='large' ,rotation=90,ha='left',va='bottom')) #Blue
ybox_r = TextArea('GOES{} 1.0-8.0 A'.format(goes[0]), textprops=dict(color='#FF0000', size='large' ,rotation=90,ha='left',va='bottom')) #Red

##Connect vertically with VPacker
ybox_primary = VPacker(children=[ybox_r, ybox_b],align='bottom', pad=0, sep=2*fs_def)

##If you want to plot both primary and secondary
##Create secondary in the same way and connect horizontally with HPacker
if len(goes)==2:
    ybox_p = TextArea('GOES{} 0.5-4.0 A'.format(goes[1]), textprops=dict(color='#5400A8', size='large' ,rotation=90,ha='left',va='bottom')) #purple
    ybox_y = TextArea('GOES{} 1.0-8.0 A'.format(goes[1]), textprops=dict(color='#FFA500', size='large' ,rotation=90,ha='left',va='bottom')) #yellow
    ybox_secondary = VPacker(children=[ybox_y, ybox_p],align='bottom', pad=0, sep=2*fs_def)
    ybox = HPacker(children=[ybox_primary,ybox_secondary],align='bottom', pad=0, sep=fs_def)
    ##Installation in Anchored Offsetbox The installation location is determined for each loc number.
    ## loc=5 is'center right'In other words, the right center of the figure
    anchored_ybox = AnchoredOffsetbox(loc=5, child=ybox, pad=0., frameon=False, bbox_to_anchor=(1-margin*0.25, 0.5), bbox_transform=fig.transFigure, borderpad=0.)
else: ##If only primary, just primary
    ybox = ybox_primary
    anchored_ybox = AnchoredOffsetbox(loc=5, child=ybox, pad=0., frameon=False, bbox_to_anchor=(1-margin*0.5, 0.5), bbox_transform=fig.transFigure, borderpad=0.)
ax[-1].add_artist(anchored_ybox)

Readability is also low here.

In the original image, the legend is rotated 90 degrees to the right outside of the frame. In matplotlib, you can put the legend out of the frame, but it doesn't seem to do anything to rotate it. In addition, it seems that you can't change the color in ʻaxes.text, so you need to put different colored texts separately. This is very difficult to do just by using coordinate specifications. Therefore, I used VPacker and HPackerofmatplotlib.offsetbox`, which can be placed in relative positions, to reproduce the legend rotated outside the figure. It's a pretty brute force technique.

First, generate each text with TextArea, and concatenate them using VPacker and HPacker. The completed pseudo-legend was installed at the right end of the figure using ʻAnchored Offsetbox`.

The place to install is decided for each number of loc, and loc = 5 means to install at the right end with center right.

4.8. Creating axis labels

##Axis label creation(dummy axes)
ax_dummy = fig.add_axes([margin,y0_fig,1-margin*2,y1_fig-y0_fig],frameon=False)
ax_dummy.tick_params(labelcolor='none', top='off', bottom='off', left='off', right='off')
ax_dummy.set_xlabel('Universal Time',size='large')
ax_dummy.set_ylabel('Watts m$^{-2}$',size='large')
ax_dummy.set_title('GOES Xray Flux (1-minute data)',size='x-large')

ax_dummy.text(1, 1-cls_text_size/5, 'Begin: {} {} {} 0000 UTC'.format(t.year,t.strftime('%b'),t.day),
        horizontalalignment='right',
        verticalalignment='top',
        transform=fig.transFigure)

I used the method of displaying ʻaxes separately to display the scale inside the figure ~~ (rough technique) ~~, but thanks to that, the label setting position on the horizontal axis could not be determined. .. So I made a big dummy ʻaxes with nothing displayed and set the axis label for it. dummy_axes.png

Also, I gave up the detailed reproduction of the characters displayed on the top, bottom, left, and right of the original image (because the characters at the bottom are not necessary). However, if there is no start time of the data in the upper right, the information of the year of the data will not be known, so I decided to display this part in a hurry.

4.9. Plot of X-ray data

## X-ray data plot
color = iter(['#FF0000','#001DFF','#FFA500','#5400A8']) ## red, blue, yellow, purple 
zo = iter([4,3,2,1])
for g in goes:
    wavelength = iter(['1.0-8.0 A','0.5-4.0 A'])
    data = get_goes_csv(g,t,trange)
    xray_long  = np.ma.masked_where((data.B_AVG<=ylim[0])|(data.B_AVG>=ylim[1]),data.B_AVG)
    xray_short = np.ma.masked_where((data.A_AVG<=ylim[0])|(data.A_AVG>=ylim[1]),data.A_AVG)
    ax[0].plot(data.time_tag,xray_long, lw=px2pt(line_width),c=next(color),label='GOES{} 1.0-8.0 A'.format(g),aa=plot_antialiased,clip_on=False,zorder=next(zo))
    ax[0].plot(data.time_tag,xray_short,lw=px2pt(line_width),c=next(color),label='GOES{} 0.5-4.0 A'.format(g),aa=plot_antialiased,clip_on=False,zorder=next(zo))

I finally came to the plot of X-ray data.

Of the columns in the GOES 1-minute average data csv file, you need 'time_tag',' A_AVG', and 'B_AVG' to plot the X-ray flux. 'time_tag' stores the time per minute as it is, and'A_AVG' and 'B_AVG' store the average X-ray data corresponding to 0.05 --0.4 nm and 0.1-0.8 nm, respectively. It has been. These are taken out based on the DataFrame grammar of pandas.

zorder is the drawing order. The larger this value is, the more it is displayed in the foreground. In the original image, the value of GOES primary was displayed in the foreground, so I set the number to be so.

Another thing I'm doing is specifying clip_on = False when plotting. Since we have prepared ʻaxes for the date this time, it is necessary to call ʻaxes.plot for that number, but if you specify clip_on = False, it will be plotted outside the prepared ʻaxes. Therefore, the data for all days is plotted with ʻaxes on the first day. Furthermore, by performing this clip_on = False, the plot can be placed on the border.

In addition, the value is masked using np.ma.masked_where. The data used this time stores -99999 when the data is lost, and if this is plotted on the log scale as it is, a figure like a well-shaped potential will be created.

Furthermore, if you draw in the range of ylim = (10 ** -9, 10 ** -2) set by default, there is basically no problem, but if you change this, clip_on = False Due to the effect of setting, there is a possibility that it will be plotted in the vertical direction. for that reason,

xray_long  = np.ma.masked_where((data.B_AVG<=ylim[0])|(data.B_AVG>=ylim[1]),data.B_AVG)

In this way, all data that comes out of the set y-axis range, including -99999, is masked. This will leave the masked plot blank without connecting.

Normally, storing np.nan also makes the plot blank, so I don't think there are many examples of using numpy.ma, but when drawing a completely different graph, np. Sometimes it was more convenient to put on a mask instead of storing nan. Therefore, I am using numpy.ma here for the purpose of introduction.

5. Actually draw

** Note: Simplex fonts may not be available on each PC by default. ** ** I will draw using the above function. Basically, the required parameters are declared first,

plt.rcParams['font.family'] = 'Simplex'
plt.rcParams['font.size'] = 10
plt.rcParams['figure.dpi'] = 100

plt.rcParams['axes.linewidth'] = 1

plt.rcParams['lines.antialiased'] = False
plt.rcParams['text.antialiased'] = False
plot_antialiased = False

##All in pixel units
wdict = dict(
    grid  = 1, ##Grid line width
    spine = 1, ##Frame line width
    tick  = 1, ##Scale line width
    line  = 0.01, ##Line width of the graph
)

plot_goes([2017,9,4], 3 ,[15,13] ,'goes.png', wdict, width=6.4, ylim=(10**-9,10**-2))

You can plot it like this. goes.png It looks like this when compared side by side. compare.png

Isn't it quite reproducible? ?? Since the downloaded Simplex font is quite thin, I think that part of the characters are not displayed due to aliasing, and that it was possible to reproduce almost everything except the character information on the top, bottom, left, and right. ~~ It seems that the size of the frame will have to be tuned manually. ~~ This completes copying.

You should be able to draw without problems even if you change each input a little.

plot_goes([2017,9,4], 5 ,[15] ,'goes_extend.png', wdict, width='extend', ylim=(10**-9,10**-2))

goes_extend.png

6. Simple drawing

By [5. Actually drawing](# 5-Actually drawing), the reproduction of the original family is almost completed. From here on, let's make it a little more practical (simple) diagram.

The top two elements where the code for goes_plot () above is cumbersome are:

--Attempting to display the scale in the figure --I tried to rotate the legend 90 degrees

I think. I don't think you need to be so particular about posters and presentation slides, so I'll try to eliminate them and simplify the code. In addition, the above function goes_plot () will fetch data from the server each time it is executed, which can be time consuming and overaccessible. Therefore, we will reuse the downloaded data and make some modifications so that the details of the figure can be adjusted. (Of course, it's best to download the data file locally and load it.)

6.1. Code

plt.rcParams['font.family'] = 'Arial'
plt.rcParams['font.size'] = 11
plt.rcParams['figure.dpi'] = 100

plt.rcParams['axes.linewidth'] = 1

plt.rcParams['lines.antialiased'] = True
plt.rcParams['text.antialiased'] = True
plot_antialiased = True

##All in pixel units
wdict = dict(
    grid  = 1, ##Grid line width
    spine = 1, ##Frame line width
    tick  = 1, ##Scale line width
    line  = 2, ##Line width of the graph
)
##Convert pixel to point
def px2pt(px):
    dpi = plt.rcParams['figure.dpi']
    return px*72/dpi

def plot_goes_simple(t, trange ,goes ,savename, wdict, width=6.4, ylim=(10**-9,10**-2)):
    '''
    t      :Plot start date(type=datetime.date or list)
             (e.g. datetime.date(2017,9,4) or [2017,9,4]) 
    trange :Period to plot Daily(type=int) (e.g. trange=3 > 3 days plot)
    goes   :GOES satellite number list format to acquire(e.g. [15,13])
             index =0 is primary(1-8A>red, 0.5-4A>blue)
             index =1 is secondary(1-8A>yellow, 0.5-4A>purple)
    width  :Plot width default=6.4
             if it is 'extend', extend the width with trange
    ylim   :y-axis range default=(10**-9,10**-2)
Probably not necessary to change, but designed to be okay to change
    '''
    if type(t)==list:
        t = datetime.date(t[0],t[1],t[2])
    elif (type(t)==datetime.date)or(type(t)==datetime.datetime):
        pass
    else:
        print('t should be a list or a datetime object.')
    
    fs_def = plt.rcParams['font.size']
    dpi    = plt.rcParams['figure.dpi']
    y0_fig = plt.rcParams['figure.subplot.bottom']
    y1_fig = plt.rcParams['figure.subplot.top']
    
    ##Width setting
    if width=='extend':
        width = (dpi*6.4 + (dpi*6.4 - fs_def*10)/3*(trange-3))/dpi
    elif type(width)==float:
        width = width
    else:
        print('width is \'extend\' or float')

    ##Left and right margins are font size* 10 (Percentage to figure)
    margin = fs_def/(100*width)*10

    ##Width per day(Percentage to figure)
    w_day  = (1-2*margin)/trange
    ax = []
    
    ##frame(axes)Creation
    t_s = t
    t_e = t_s + datetime.timedelta(days=trange)
    
    fig = plt.figure(figsize=(width,4.8))
    ##If you put a normal legend in the frame, add forcibly_Without specifying with axes
    ##Default fig.add_I think subplot is fine
    #ax = fig.add_axes([margin,y0_fig,1-2*margin,y1_fig-y0_fig],yscale='log', xlim=[t_s,t_e], ylim=ylim)
    ax = fig.add_subplot(111,yscale='log', xlim=[t_s,t_e], ylim=ylim)
    new_xticks = date2num([t_s+datetime.timedelta(days=i) for i in range(trange+1)])

    ax.xaxis.set_major_locator(ticker.FixedLocator(new_xticks))
    ax.xaxis.set_major_formatter(DateFormatter('%b %d'))
    ax.xaxis.set_minor_locator(matplotlib.dates.HourLocator(interval=3))
    ax.yaxis.grid(lw=px2pt(wdict['grid']),color='black')

    ax.tick_params(axis='x', which='major', bottom=True, top=True, direction='in', length=6, width=px2pt(wdict['tick']), pad=fs_def)
    ax.tick_params(axis='x', which='minor', bottom=True, top=True, direction='in', length=3, width=px2pt(wdict['tick']))
    ax.tick_params(axis='y', which='both', right=True, left=True, direction='in', length=6, width=px2pt(wdict['tick']))
    
    for j in ax.spines.keys(): ##Adjusting the line width of the frame
        ax.spines[j].set_linewidth(px2pt(wdict['spine']))

    ## Flare class text
    f_pos = ax.get_position()

    trans = transforms.blended_transform_factory(fig.transFigure,transforms.IdentityTransform())
    X10_y_display   = ax.transData.transform((date2num(t_e),10**-3))[1]
    X_y_display     = ax.transData.transform((date2num(t_e),10**-4))[1]
    w = X10_y_display - X_y_display

    ##Relative size of font size to the entire figure
    ##Used to allow margins when placing class text so that it does not stick to the outer frame
    cls_text_size = plt.rcParams['font.size']/(plt.rcParams['figure.figsize'][0])/72
    
    for i,cls in enumerate(['X', 'M', 'C', 'B', 'A']):
        pos = [f_pos.x1,X10_y_display - w/2 - i*w]
        if (pos[1]>=ax.transData.transform((date2num(t_e),ax.get_ylim()[0]))[1])&\
           (pos[1]<=ax.transData.transform((date2num(t_e),ax.get_ylim()[1]))[1]):
            ax.text(pos[0]+cls_text_size/2,pos[1],cls,horizontalalignment='center',verticalalignment='center',transform=trans)

    ##If you uncomment here, you can put the legend outside the frame
    '''
    ## GOES number and wavelength
    ybox_b = TextArea('GOES{} 0.5-4.0 A'.format(goes[0]), textprops=dict(color='#001DFF', size='large' ,rotation=90,ha='left',va='bottom')) #Blue
    ybox_r = TextArea('GOES{} 1.0-8.0 A'.format(goes[0]), textprops=dict(color='#FF0000', size='large' ,rotation=90,ha='left',va='bottom')) #Red
    ybox_primary = VPacker(children=[ybox_r, ybox_b],align='bottom', pad=0, sep=2*fs_def)
    if len(goes)==2:
        ybox_p = TextArea('GOES{} 0.5-4.0 A'.format(goes[1]), textprops=dict(color='#5400A8', size='large' ,rotation=90,ha='left',va='bottom')) #purple
        ybox_y = TextArea('GOES{} 1.0-8.0 A'.format(goes[1]), textprops=dict(color='#FFA500', size='large' ,rotation=90,ha='left',va='bottom')) #yellow
        ybox_secondary = VPacker(children=[ybox_y, ybox_p],align='bottom', pad=0, sep=2*fs_def)
        ybox = HPacker(children=[ybox_primary,ybox_secondary],align='bottom', pad=0, sep=fs_def)
        anchored_ybox = AnchoredOffsetbox(loc=5, child=ybox, pad=0., frameon=False, bbox_to_anchor=(1-margin*0.25, 0.5), 
                                          bbox_transform=fig.transFigure, borderpad=0.)
    else:
        ybox = ybox_primary
        anchored_ybox = AnchoredOffsetbox(loc=5, child=ybox, pad=0., frameon=False, bbox_to_anchor=(1-margin*0.5, 0.5), 
                                          bbox_transform=fig.transFigure, borderpad=0.)
    ax.add_artist(anchored_ybox)
    '''
    
    ##Axis label creation
    ax.set_xlabel('Universal Time',size='large')
    ax.set_ylabel('Watts m$^{-2}$',size='large')
    ax.set_title('GOES Xray Flux (1-minute data)',size='x-large')

    ax.text(1, 1-cls_text_size/5, 'Begin: {} {} {} 0000 UTC'.format(t.year,t.strftime('%b'),t.day),
            horizontalalignment='right',
            verticalalignment='top',
            transform=fig.transFigure)

    ## X-ray data plot
    color = iter(['#FF0000','#001DFF','#FFA500','#5400A8']) ## red, blue, yellow, purple 
    zo = iter([4,3,2,1])

    for i,g in enumerate(goes):
        data = data_temp[i]
        xray_long  = np.ma.masked_where((data.B_AVG<=ylim[0])|(data.B_AVG>=ylim[1]),data.B_AVG)
        xray_short = np.ma.masked_where((data.A_AVG<=ylim[0])|(data.A_AVG>=ylim[1]),data.A_AVG)
        ax.plot(data.time_tag,xray_long, lw=px2pt(wdict['line']),c=next(color),label='GOES{} 1.0-8.0 A'.format(g),aa=plot_antialiased,clip_on=False,zorder=next(zo))
        ax.plot(data.time_tag,xray_short,lw=px2pt(wdict['line']),c=next(color),label='GOES{} 0.5-4.0 A'.format(g),aa=plot_antialiased,clip_on=False,zorder=next(zo))
        
    ## simple legend
    leg = ax.legend(handlelength=0, handletextpad=0, ncol=2)
    for line, text in zip(leg.get_lines(), leg.get_texts()):
        text.set_color(line.get_color())
    
    ## save graph
    fig.savefig(savename)

The notation around the variable ʻax has changed slightly because it is no longer necessary to put multiple ʻaxes. Please note the difference from plot_goes () when arranging.

## simple legend
leg = ax.legend(handlelength=0, handletextpad=0, ncol=2)
for line, text in zip(leg.get_lines(), leg.get_texts()):
    text.set_color(line.get_color())

So, I've made the Legend markers invisible and changed the color of the Legend text to the color that corresponds to each plot.

6.2. Drawing

I've tried a few, but the following parameters seem to be just right.

plt.rcParams['font.family'] = 'Arial'
plt.rcParams['font.size'] = 11
plt.rcParams['figure.dpi'] = 100

plt.rcParams['axes.linewidth'] = 1

plt.rcParams['lines.antialiased'] = True
plt.rcParams['text.antialiased'] = True
plot_antialiased = True

##All in pixel units
wdict = dict(
    grid  = 1, ##Grid line width
    spine = 1, ##Frame line width
    tick  = 1, ##Scale line width
    line  = 2, ##Line width of the graph
)

Store the X-ray data in data_temp and plot it as follows.

t = [2017,9,4]
trange = 3
goes = [15,13]

data_temp = []
for g in goes:
    data_temp.append(get_goes_csv(g,t,trange))

plot_goes_simple(t, trange, goes, 'goes_simple.png', wdict, width=6.4, ylim=(10**-9,10**-2))

goes_simple.png

If it looks like this, I think it's relatively easy to see even if you put it on a poster.

7. Conclusion

If I was particular about it, I couldn't stop. ʻAxes` I learned a lot about things around me. If you have the opportunity, please use it for papers, presentation slides, posters, etc. If it helps, even when I meet you somewhere, if you say "it was useful", I'm very happy that you feel like "I did it".

If you have any questions or improvements, please feel free to contact us. GOES-16 has real-time data released in json format, but you only have the data for the last week. Will past data be released soon? I'll be in trouble if I don't do it. ~~ Once the format is decided, I would like to make it so that it can be supported.

References

matplotlib in general matplotlib.figure.Figure — Matplotlib 3.1.2 documentation Basic knowledge of matplotlib that I wanted to know early, or the story of an artist who can adjust the appearance --Qiita A very summary of matplotlib-Qiita matplotlib reverse lookup memo (simple graph) --Qiita matplotlib reverse lookup memo (frame edition) --Qiita

around tick, spine Set the axis scale of the graph of time series data with matplotlib.dates | Analysis Train Matplotlib scale and scale label, scale line Create a graph with borders removed with matplotlib --Qiita

text placement and coordinate specification (transform) Text properties and layout — Matplotlib 3.1.2 documentation Transformations Tutorial — Matplotlib 3.1.2 documentation

Setting axis labels for multiple subplots (setting dummy axes) python - pyplot axes labels for subplots - Stack Overflow

Adjust legend marker and legend text color python - How do you just show the text label in plot legend? (e.g. remove a label's line in the legend) - Stack Overflow Option to set the text color in legend to be same as the line · Issue #10720 · matplotlib/matplotlib · GitHub

Creating a legend rotated 90 degrees (matplotlib.offsetbox) python - Matplotlib: y-axis label with multiple colors - Stack Overflow matplotlib.offsetbox — Matplotlib 3.1.2 documentation

numpy.ma The numpy.ma module — NumPy v1.17 Manual

Recommended Posts

Make your own NOAA / GOES familiar X-ray flux 3days plot in Python
[Python] logging in your own module
Create your own Linux commands in Python
[LLDB] Create your own command in Python
Easily use your own functions in Python
Get your own IP address in Python
Make your own module quickly with setuptools (python)
Make a joyplot-like plot of R in python
Import your own modules in Grasshopper's Python development
Create your own Big Data in Python for validation
Create your own Random Dot Stereogram (RDS) in Python.
Try to improve your own intro quiz in Python
Use the CASA Toolkit in your own Python environment
[Road to intermediate Python] Define in in your own class
Try sorting your own objects with priority queue in Python
Let's make a number guessing game in your own language!
Don't make test.py in Python!
Make a bookmarklet in Python
Make Opencv available in Python
Make python segfault in 2 lines
Plot geographic information in Python
X-ray diffraction simulation in Python
Create your own graph structure class and its drawing in python
Specify your own class in class argument and return type annotation in Python