Commit 608e2c83 authored by Jonathan  Minz's avatar Jonathan Minz
Browse files

ORIGINAL: Python script for plotting Doppler Lidar data share by Dr. Oliver...

ORIGINAL: Python script for plotting Doppler Lidar data share by Dr. Oliver Branch. Created by Dr. Florian Spaeth. 
parent 1c22db8b
Loading
Loading
Loading
Loading

Phython_hpl_plots.py

0 → 100644
+476 −0
Original line number Diff line number Diff line
#!/usr/bin/env python3
# -*- coding: utf-8 -*-

import matplotlib
matplotlib.use('Agg') # use non gui backend

import numpy as np
#from netCDF4 import Dataset
import os
from datetime import datetime, timedelta
#import matplotlib.pyplot as plt
from matplotlib import pylab as plt
import sys

import datetime
from datetime import datetime
from datetime import timedelta 
import matplotlib.dates as mdates
import matplotlib.pyplot as plt
from numpy import meshgrid
import matplotlib.colors as mc
from matplotlib.colors import ListedColormap, LinearSegmentedColormap 
from mpl_toolkits.axes_grid1.axes_divider import make_axes_locatable


#function imports .hpl file and stores information in dictionary data_temp
def hpl2dict(file_path):
    #import hpl files into intercal storage
    with open(file_path, 'r') as text_file:
        lines=text_file.readlines()
    
    #write lines into Dictionary
    data_temp=dict()
    
    header_n=17 #length of header
    data_temp['filename']=lines[0].split()[-1]
    data_temp['system_id']=int(lines[1].split()[-1])
    data_temp['number_of_gates']=int(lines[2].split()[-1])
    data_temp['range_gate_length_m']=float(lines[3].split()[-1])
    data_temp['gate_length_pts']=int(lines[4].split()[-1])
    data_temp['pulses_per_ray']=int(lines[5].split()[-1])
    data_temp['number_of_waypoints_in_file']=int(lines[6].split()[-1])
    rays_n=(len(lines)-header_n)/(data_temp['number_of_gates']+1)
    if not rays_n.is_integer():
        raise Exception('Number of lines does not match expected format')
    data_temp['no_of_rays_in_file']=int(rays_n)
    data_temp['scan_type']=' '.join(lines[7].split()[2:])
    data_temp['focus_range']=lines[8].split()[-1]
    data_temp['start_time']=' '.join(lines[9].split()[-2:])
    data_temp['resolution']=('%s %s' % (lines[10].split()[-1],'m s-1'))
    data_temp['range_gates']=np.arange(0,data_temp['number_of_gates'])
    data_temp['center_of_gates']=(data_temp['range_gates']+0.5)*data_temp['range_gate_length_m']
    
    #dimensions of data set
    gates_n=data_temp['number_of_gates']
    rays_n=data_temp['no_of_rays_in_file']
    
    #keys for measurement variables are predefined as symetric numpy arrays filled with NaN values
    data_temp['radial_velocity']=np.full([gates_n,rays_n],np.nan) #m/s
    data_temp['intensity']=np.full([gates_n,rays_n],np.nan) #SNR+1
    data_temp['beta']=np.full([gates_n,rays_n],np.nan) #m-1 sr-1
    data_temp['elevation']=np.full(rays_n,np.nan) #degrees
    data_temp['azimuth']=np.full(rays_n,np.nan) #degrees
    data_temp['decimal_time']=np.full(rays_n,np.nan) #hours
    data_temp['pitch']=np.full(rays_n,np.nan) #degrees
    data_temp['roll']=np.full(rays_n,np.nan) #degrees
    for ri in range(0,rays_n): #loop through rays
        lines_temp=lines[header_n+(ri*gates_n)+ri+1:header_n+(ri*gates_n)+gates_n+ri+1]
        header_temp=np.asarray(lines[header_n+(ri*gates_n)+ri].split(),dtype=float)
        data_temp['decimal_time'][ri]=header_temp[0]
        data_temp['azimuth'][ri]=header_temp[1]
        data_temp['elevation'][ri]=header_temp[2]
        data_temp['pitch'][ri]=header_temp[3]
        data_temp['roll'][ri]=header_temp[4]
        for gi in range(0,gates_n): #loop through range gates
            line_temp=np.asarray(lines_temp[gi].split(),dtype=float)
            data_temp['radial_velocity'][gi,ri]=line_temp[1]
            data_temp['intensity'][gi,ri]=line_temp[2]
            data_temp['beta'][gi,ri]=line_temp[3]    
        
    return data_temp


'''
#=======================
# select data files
#=======================

# print count of specified arguments
print("Sie haben %d Argument(e) angegeben." % len(sys.argv))

# repeat for each argument
for argument in sys.argv:
  # print current argument
  print(argument)
'''

def file2plot(datetime_now, hour_in, loop_start):
    #hour_in = int(12)
    #hour_in = int(sys.argv[1])
#    datetime_now = datetime(year=2020, month=7, day=1, hour=hour_in, minute=3, second=23)
    #datetime_now = datetime(year=2020, month=7, day=1, hour=12, minute=3, second=23)
    print('Now:   ', datetime_now)

    #/ipmserver/ag_fernerkundung/Projekte/LAFO/DL/DL_UAE
    #/data/data_shared/DopplerLidar/DL_UAE/Proc/
    #/data/data_shared/DopplerLidar/DL_UAE/Quicklooks
    #/data/data_personal/quicklooks/DL/metek_server

    path1=os.path.join('/data/data_shared/DopplerLidar/DL_156/Proc', (datetime_now - timedelta(hours=3)).strftime("%Y"), (datetime_now - timedelta(hours=3)).strftime("%Y%m"), (datetime_now - timedelta(hours=3)).strftime("%Y%m%d"))
    file1 = 'Stare_156_' + (datetime_now - timedelta(hours=3)).strftime("%Y%m%d_%H") + '.hpl'
    file_path1 = os.path.join(path1, file1)
    path2=os.path.join('/data/data_shared/DopplerLidar/DL_156/Proc', (datetime_now - timedelta(hours=2)).strftime("%Y"), (datetime_now - timedelta(hours=2)).strftime("%Y%m"), (datetime_now - timedelta(hours=2)).strftime("%Y%m%d"))
    file2 = 'Stare_156_' + (datetime_now - timedelta(hours=2)).strftime("%Y%m%d_%H") + '.hpl'
    file_path2 = os.path.join(path2, file2)
    path3=os.path.join('/data/data_shared/DopplerLidar/DL_156/Proc', (datetime_now - timedelta(hours=1)).strftime("%Y"), (datetime_now - timedelta(hours=1)).strftime("%Y%m"), (datetime_now - timedelta(hours=1)).strftime("%Y%m%d"))
    file3 = 'Stare_156_' + (datetime_now - timedelta(hours=1)).strftime("%Y%m%d_%H") + '.hpl'
    file_path3 = os.path.join(path3, file3)
    path_out1=os.path.join('/data/data_personal/quicklooks/DL', datetime_now.strftime("%Y"), datetime_now.strftime("%Y%m"), datetime_now.strftime("%Y%m%d"))
    file_out1='LAFO_'+datetime_now.strftime("%Y%m%d_%H")+'_DLw.png'
    file_path_out1 = os.path.join(path_out1, file_out1)
    #path_out2=os.path.join('/data/data_shared/DopplerLidar/DL_UAE/Quicklooks', datetime_now.strftime("%Y"), datetime_now.strftime("%m%d"))
    path_out2=os.path.join('/data/data_personal/quicklooks/DL/metek_server', datetime_now.strftime("%Y"), datetime_now.strftime("%m%d"))

    file_out2=datetime_now.strftime("%H")+'01.LAFO_DLw.png'
    file_path_out2 = os.path.join(path_out2, file_out2)
    file_out3=datetime_now.strftime("%H")+'01.LAFO_DLbeta.png'
    file_path_out3 = os.path.join(path_out2, file_out3)
    file_out4=datetime_now.strftime("%H")+'01.LAFO_DLsnr.png'
    file_path_out4 = os.path.join(path_out2, file_out4)
    file_out5=datetime_now.strftime("%H")+'01.LAFO_DLw_4km.png'
    file_path_out5 = os.path.join(path_out2, file_out5)
    file_out6=datetime_now.strftime("%H")+'01.LAFO_DLbeta_4km.png'
    file_path_out6 = os.path.join(path_out2, file_out6)
    print('file1: ', file_path1)
    print('file2: ', file_path2)
    print('file3: ', file_path3)
    print('file_out1: ', file_path_out1)
    print('file_out2: ', file_path_out2)
    print('file_out3: ', file_path_out3)
    print('file_out4: ', file_path_out4)
    print('file_out5: ', file_path_out5)
    print('file_out6: ', file_path_out6)
    
    # quit python script
    #sys.exit(0)

    file_exists = 1
    if not os.path.exists(file_path1):
        #raise Exception('%s cannot be found' %os.path.basename(file_path1))
        file_exists = 0
        print('%s cannot be found' %os.path.basename(file_path1))
    
    if not os.path.exists(file_path2):
        #raise Exception('%s cannot be found' %os.path.basename(file_path2))
        file_exists = 0
        print('%s cannot be found' %os.path.basename(file_path2))

    if not os.path.exists(file_path3):
        #raise Exception('%s cannot be found' %os.path.basename(file_path3))
        file_exists = 0
        print('%s cannot be found' %os.path.basename(file_path3))

    if file_exists == 1:
        if not os.path.isdir(path_out1):
            #raise Exception('%s cannot be found' %os.path.basename(file_path2))
            print('Folder does not exist.')
            os.makedirs(path_out1)
            print('Quicklook folder created: %s' %os.path.basename(path_out1))

        if not os.path.isdir(path_out2):
            #raise Exception('%s cannot be found' %os.path.basename(file_path2))
            print('Folder does not exist.')
            os.makedirs(path_out2)
            print('Quicklook folder created: %s' %os.path.basename(path_out2))

        # read data from files
        data_temp1=hpl2dict(file_path1)
        print(data_temp1['filename'] + ' loaded.')
        
        data_temp2=hpl2dict(file_path2)
        print(data_temp2['filename'] + ' loaded.')
        
        data_temp3=hpl2dict(file_path3)
        print(data_temp3['filename'] + ' loaded.')
        
        
        # merge data_temp1, data_temp2, and date_temp3
        gates_n1=data_temp1['number_of_gates']
        rays_n1=data_temp1['no_of_rays_in_file']
        
        gates_n2=data_temp2['number_of_gates']
        rays_n2=data_temp2['no_of_rays_in_file']
        
        gates_n3=data_temp3['number_of_gates']
        rays_n3=data_temp3['no_of_rays_in_file']
        
        data_12=dict()
        
        data_12['filename']='Stare_115_' + (datetime_now - timedelta(hours=3)).strftime("%Y%m%d_%H") + '-' + (datetime_now - timedelta(hours=1)).strftime("%H") + '.hpl'
        data_12['system_id']=data_temp1['system_id']
        data_12['number_of_gates']=data_temp1['number_of_gates']
        data_12['range_gate_length_m']=data_temp1['range_gate_length_m']
        data_12['gate_length_pts']=data_temp1['gate_length_pts']
        data_12['pulses_per_ray']=data_temp1['pulses_per_ray']
        data_12['number_of_waypoints_in_file']=data_temp1['number_of_waypoints_in_file']
        data_12['no_of_rays_in_file']=int(rays_n1+rays_n2+rays_n3)
        data_12['scan_type']=data_temp1['scan_type']
        data_12['focus_range']=data_temp1['focus_range']
        data_12['start_time']=data_temp1['start_time']
        data_12['resolution']=data_temp1['resolution']
        data_12['range_gates']=np.arange(0,data_temp1['number_of_gates'])
        data_12['center_of_gates']=(data_temp1['range_gates']+0.5)*data_temp1['range_gate_length_m']
        
        #dimensions of data set
        gates_n=data_12['number_of_gates']
        rays_n=data_12['no_of_rays_in_file']
        
        #keys for measurement variables are predefined as symetric numpy arrays filled with NaN values
        data_12['radial_velocity']=np.full([gates_n,rays_n],np.nan) #m/s
        data_12['intensity']=np.full([gates_n,rays_n],np.nan) #SNR+1
        data_12['beta']=np.full([gates_n,rays_n],np.nan) #m-1 sr-1
        data_12['elevation']=np.full(rays_n,np.nan) #degrees
        data_12['azimuth']=np.full(rays_n,np.nan) #degrees
        data_12['decimal_time']=np.full(rays_n,np.nan) #hours
        data_12['date']=np.full(rays_n,np.nan) #date string
        data_12['pitch']=np.full(rays_n,np.nan) #degrees
        data_12['roll']=np.full(rays_n,np.nan) #degrees
        
        data_12['elevation'][0:rays_n1]=data_temp1['elevation'] #degrees
        data_12['elevation'][rays_n1:rays_n1+rays_n2]=data_temp2['elevation'] #degrees
        data_12['elevation'][rays_n1+rays_n2:]=data_temp3['elevation'] #degrees
        data_12['azimuth'][0:rays_n1]=data_temp1['azimuth'] #degrees
        data_12['azimuth'][rays_n1:rays_n1+rays_n2]=data_temp2['azimuth'] #degrees
        data_12['azimuth'][rays_n1+rays_n2:]=data_temp3['azimuth'] #degrees
        data_12['decimal_time'][0:rays_n1]=data_temp1['decimal_time'] #hours
        data_12['decimal_time'][rays_n1:rays_n1+rays_n2]=data_temp2['decimal_time'] #hours
        data_12['decimal_time'][rays_n1+rays_n2:]=data_temp3['decimal_time'] #hours
        data_12['pitch'][0:rays_n1]=data_temp1['pitch'] #degrees
        data_12['pitch'][rays_n1:rays_n1+rays_n2]=data_temp2['pitch'] #degrees
        data_12['pitch'][rays_n1+rays_n2:]=data_temp3['pitch'] #degrees
        data_12['roll'][0:rays_n1]=data_temp1['roll'] #degrees
        data_12['roll'][rays_n1:rays_n1+rays_n2]=data_temp2['roll'] #degrees
        data_12['roll'][rays_n1+rays_n2:]=data_temp3['roll'] #degrees
        
        data_12['radial_velocity'][:,0:rays_n1]=data_temp1['radial_velocity'] #m/s
        data_12['radial_velocity'][:,rays_n1:rays_n1+rays_n2]=data_temp2['radial_velocity'] #m/s
        data_12['radial_velocity'][:,rays_n1+rays_n2:]=data_temp3['radial_velocity'] #m/s
        data_12['intensity'][:,0:rays_n1]=data_temp1['intensity'] #m/s
        data_12['intensity'][:,rays_n1:rays_n1+rays_n2]=data_temp2['intensity'] #m/s
        data_12['intensity'][:,rays_n1+rays_n2:]=data_temp3['intensity'] #m/s
        data_12['beta'][:,0:rays_n1]=data_temp1['beta'] #m/s
        data_12['beta'][:,rays_n1:rays_n1+rays_n2]=data_temp2['beta'] #m/s
        data_12['beta'][:,rays_n1+rays_n2:]=data_temp3['beta'] #m/s
        
        #=======================
        # PLOTTING
        #=======================
        cmap = LinearSegmentedColormap.from_list("mycmap", ["lime","blue", "white", "red", "yellow"])
        
        
#        plt.rc('font', **{'family': 'serif', 'serif': ['Computer Modern']})
#        plt.rc('text', usetex=True)
#        plt.rc('text')
#        plt.rc('axes', labelsize=14)
#        plt.rc('xtick', labelsize=14)
#        plt.rc('ytick', labelsize=14)
        
        # masks data
        data_12['radial_velocity'] = np.ma.masked_less(data_12['radial_velocity'], -3)
        data_12['radial_velocity'] = np.ma.masked_greater(data_12['radial_velocity'], 3)
        data_12['beta'] = np.ma.masked_less(data_12['beta'], 1e-7)
        data_12['beta'] = np.ma.masked_greater(data_12['beta'], 1e-4)
        data_12['intensity'] = np.ma.masked_less(data_12['intensity'], 0.8)
        data_12['intensity'] = np.ma.masked_greater(data_12['intensity'], 1.6)
        
        temp1 = [data_temp1['start_time']]*len(data_temp1['decimal_time'])
        temp2 = [data_temp2['start_time']]*len(data_temp2['decimal_time'])
        temp3 = [data_temp3['start_time']]*len(data_temp3['decimal_time'])
        temp= temp1+temp2+temp3

        t_start = data_12['start_time'][:4]+'-' + data_12['start_time'][4:6] +\
                    '-' + data_12['start_time'][6:8]
        time = data_12['decimal_time']
        height = data_12['center_of_gates']
        
        Hours = time
        Minutes = 60 * (Hours % 1)
        Seconds = 60 * (Minutes % 1)
        
        timestamps = []
        for n in range(len(time)):
            day = temp[n][:4]+'-' + temp[n][4:6] +'-' + temp[n][6:8]+"%d:%02d:%02d" % (Hours[n], Minutes[n], Seconds[n])
            day1 = temp[n][:4]+'-' + temp[n][4:6] +'-' + temp[n][6:8]+"%d:%02d:%02d" % (Hours[0], Minutes[0], Seconds[0])
#            day = t_start+"%d:%02d:%02d" % (Hours[n], Minutes[n], Seconds[n])
#            day1 = t_start+"%d:%02d:%02d" % (Hours[0], Minutes[0], Seconds[0])
#            timestamps.append(datetime.datetime.strptime(day, '%Y-%m-%d%H:%M:%S'))
            timestamps.append(datetime.strptime(day, '%Y-%m-%d%H:%M:%S'))
        
        # plotting
        
        fig = plt.figure(figsize=[10, 6.2], dpi=100)
#        fig = plt.figure(figsize=[10, 6.2], dpi=75)
        
        xx, yy = meshgrid(timestamps, height)
        data_doppler = data_12['radial_velocity']
        data_beta = data_12['beta']
        data_snr = data_12['intensity']
        
        pc = plt.pcolormesh(xx, yy/1000, data_doppler, cmap=cmap)
        plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%H:%M'))
#        plt.yticks(ticks=[0, 2, 4, 6, 8, 10], labels=['0', '2', '4', '6', '8', '10'])	# 'ticks' and 'labels' in new python version not allowed any more
        plt.yticks([0, 2, 4, 6, 8, 10], ['0', '2', '4', '6', '8', '10'])
        plt.tick_params('both', length=8, width=0.5, which='major')
        plt.ylabel('Range (km)')
        plt.xlabel('Time (UTC)')
        ax2_divider = make_axes_locatable(plt.gca())
        cax2 = ax2_divider.append_axes("top", size=0.12, pad="20%")
        clb = fig.colorbar(pc, cax=cax2, orientation="horizontal")
        clb.ax.set_title('Doppler Wind Velocity, m/s', size=24)
#        plt.gcf().text(0.1, 0.78, 'Date:'+'\t'+'\t'+t_start, fontsize=14)
#        plt.gcf().text(0.3, 0.78, 'Time:'+'\t'+'\t'+day1[10:15] + '\t' + '\t' + 'to' +
#                       '\t' + '\t'+day[0:10]+' '+day[9:15], fontsize=14)
        plt.gcf().text(0.1, 0.78, t_start+' '+day1[10:15], fontsize=14)
        plt.gcf().text(0.3, 0.78, 'to'+' '+day[0:10]+' '+day[10:15], fontsize=14)
        plt.gcf().text(0.75, 0.78, 'Az = '+str(int(data_12['azimuth'][0])),
                       fontsize=14)
        plt.gcf().text(0.9, 0.78, 'El = '+str(int(data_12['elevation'][0])),
                       fontsize=14)
        plt.tight_layout()
        #       plt.show()
        #       plt.savefig('Matplotlib_save_plot-2.png')
        # plt.savefig(file_path_out1)
        # print('Quicklook created: ' + file_path_out1)
        plt.savefig(file_path_out2)
        print('Quicklook created: ' + file_path_out2)
        print('Doppler Quicklook created: ', datetime.utcnow())
        
#       plot beta data
        fig = plt.figure(figsize=[10, 6.2], dpi=100)
#        fig = plt.figure(figsize=[10, 6.2], dpi=75)
        pc = plt.pcolormesh(xx, yy/1000, data_beta, cmap=plt.get_cmap('rainbow'))
        plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%H:%M'))
#        plt.yticks(ticks=[0, 2, 4, 6, 8, 10], labels=['0', '2', '4', '6', '8', '10'])	# 'ticks' and 'labels' in new python version not allowed any more
        plt.yticks([0, 2, 4, 6, 8, 10], ['0', '2', '4', '6', '8', '10'])
        plt.tick_params('both', length=8, width=0.5, which='major')
        plt.ylabel('Range (km)')
        plt.xlabel('Time (UTC)')
        ax2_divider = make_axes_locatable(plt.gca())
        cax2 = ax2_divider.append_axes("top", size=0.12, pad="20%")
        clb = fig.colorbar(pc, cax=cax2, orientation="horizontal")
        clb.ax.set_title('Backscatter Signal', size=24)
        plt.gcf().text(0.1, 0.78, t_start+' '+day1[10:15], fontsize=14)
        plt.gcf().text(0.3, 0.78, 'to'+' '+day[0:10]+' '+day[10:15], fontsize=14)
        plt.gcf().text(0.75, 0.78, 'Az = '+str(int(data_12['azimuth'][0])),
                       fontsize=14)
        plt.gcf().text(0.9, 0.78, 'El = '+str(int(data_12['elevation'][0])),
                       fontsize=14)
        plt.tight_layout()
        plt.savefig(file_path_out3)
        print('Quicklook created: ' + file_path_out3)
        print('beta    Quicklook created: ', datetime.utcnow())

#       plot snr data
        fig = plt.figure(figsize=[10, 6.2], dpi=100)
#        fig = plt.figure(figsize=[10, 6.2], dpi=75)
        pc = plt.pcolormesh(xx, yy/1000, data_snr, cmap=plt.get_cmap('rainbow'))
        plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%H:%M'))
#        plt.yticks(ticks=[0, 2, 4, 6, 8, 10], labels=['0', '2', '4', '6', '8', '10'])	# 'ticks' and 'labels' in new python version not allowed any more
        plt.yticks([0, 2, 4, 6, 8, 10], ['0', '2', '4', '6', '8', '10'])
        plt.tick_params('both', length=8, width=0.5, which='major')
        plt.ylabel('Range (km)')
        plt.xlabel('Time (UTC)')
        ax2_divider = make_axes_locatable(plt.gca())
        cax2 = ax2_divider.append_axes("top", size=0.12, pad="20%")
        clb = fig.colorbar(pc, cax=cax2, orientation="horizontal")
        clb.ax.set_title('SNR+1', size=24)
        plt.gcf().text(0.1, 0.78, t_start+' '+day1[10:15], fontsize=14)
        plt.gcf().text(0.3, 0.78, 'to'+' '+day[0:10]+' '+day[10:15], fontsize=14)
        plt.gcf().text(0.75, 0.78, 'Az = '+str(int(data_12['azimuth'][0])),
                       fontsize=14)
        plt.gcf().text(0.9, 0.78, 'El = '+str(int(data_12['elevation'][0])),
                       fontsize=14)
        plt.tight_layout()
        plt.savefig(file_path_out4)
        print('Quicklook created: ' + file_path_out4)
        print('snr     Quicklook created: ', datetime.utcnow())

#       plot Doppler data up to 4 km
        fig = plt.figure(figsize=[10, 6.2], dpi=100)
#        fig = plt.figure(figsize=[10, 6.2], dpi=75)
        pc = plt.pcolormesh(xx[0:150,:], yy[0:150,:]/1000, data_doppler[0:150,:], cmap=cmap)
        plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%H:%M'))
#        plt.yticks(ticks=[0, 1, 2, 3, 4], labels=['0', '1', '2', '3', '4'])	# 'ticks' and 'labels' in new python version not allowed any more
        plt.yticks([0, 1, 2, 3, 4], ['0', '1', '2', '3', '4'])
        plt.tick_params('both', length=8, width=0.5, which='major')
        plt.ylabel('Range (km)')
        plt.xlabel('Time (UTC)')
        ax2_divider = make_axes_locatable(plt.gca())
        cax2 = ax2_divider.append_axes("top", size=0.12, pad="20%")
        clb = fig.colorbar(pc, cax=cax2, orientation="horizontal")
        clb.ax.set_title('Doppler Wind Velocity, m/s', size=24)
        plt.gcf().text(0.1, 0.78, t_start+' '+day1[10:15], fontsize=14)
        plt.gcf().text(0.3, 0.78, 'to'+' '+day[0:10]+' '+day[10:15], fontsize=14)
        plt.gcf().text(0.75, 0.78, 'Az = '+str(int(data_12['azimuth'][0])),
                       fontsize=14)
        plt.gcf().text(0.9, 0.78, 'El = '+str(int(data_12['elevation'][0])),
                       fontsize=14)
        plt.tight_layout()
        plt.savefig(file_path_out5)
        print('Quicklook created: ' + file_path_out5)
        print('Doppler Quicklook created: ', datetime.utcnow())

#       plot beta data up to 4 km
        fig = plt.figure(figsize=[10, 6.2], dpi=100)
#        fig = plt.figure(figsize=[10, 6.2], dpi=75)
        pc = plt.pcolormesh(xx[0:150,:], yy[0:150,:]/1000, data_beta[0:150,:], cmap=plt.get_cmap('rainbow'))
        plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%H:%M'))
#        plt.yticks(ticks=[0, 1, 2, 3, 4], labels=['0', '1', '2', '3', '4'])	# 'ticks' and 'labels' in new python version not allowed any more
        plt.yticks([0, 1, 2, 3, 4], ['0', '1', '2', '3', '4'])
        plt.tick_params('both', length=8, width=0.5, which='major')
        plt.ylabel('Range (km)')
        plt.xlabel('Time (UTC)')
        ax2_divider = make_axes_locatable(plt.gca())
        cax2 = ax2_divider.append_axes("top", size=0.12, pad="20%")
        clb = fig.colorbar(pc, cax=cax2, orientation="horizontal")
        clb.ax.set_title('Backscatter Signal', size=24)
        plt.gcf().text(0.1, 0.78, t_start+' '+day1[10:15], fontsize=14)
        plt.gcf().text(0.3, 0.78, 'to'+' '+day[0:10]+' '+day[10:15], fontsize=14)
        plt.gcf().text(0.75, 0.78, 'Az = '+str(int(data_12['azimuth'][0])),
                       fontsize=14)
        plt.gcf().text(0.9, 0.78, 'El = '+str(int(data_12['elevation'][0])),
                       fontsize=14)
        plt.tight_layout()
        plt.savefig(file_path_out6)
        print('Quicklook created: ' + file_path_out6)
        print('beta    Quicklook created: ', datetime.utcnow())

    return

# crontab -e
# 3 * * * * python3 /ipmserver/ag_fernerkundung/Projekte/LAFO/DL/DL_UAE/hpl2plot_time3h.py


# print count of specified arguments
print("Sie haben %d Argument(e) angegeben." % len(sys.argv))

# repeat for each argument
for argument in sys.argv:
  # print current argument
  print(argument)


# current UTC time
datetime_now = datetime.utcnow()

# modify datetime_now
#datetime_now = datetime(year=2020, month=11, day=16, hour=1, minute=3, second=23)
#datetime_now=datetime_now - timedelta(hours=2)

print('Current UTC date/time: ', datetime_now)

hh=datetime_now.strftime("%H")
loop_start=hh
file2plot(datetime_now, hh, loop_start)

'''
datetime_now = datetime(year=2020, month=12, day=5, hour=0, minute=3, second=23)
loop_start = 0
loop_end = 24
for hh in range(loop_start, loop_end):
    print(hh, datetime_now + timedelta(hours=hh))
#    datetime_now=datetime_now - timedelta(hours=hh)
#    file2plot(datetime_now, hh, loop_start)
    file2plot(datetime_now + timedelta(hours=hh), hh, loop_start)
'''