Search code examples
pandasdatetimedata-processing

How do I select data for a period of time say a month for plotting?


I want to be able to visualise and analyse data for a specific period of time. The data is in different files so i read each file and appended them into an array. The data I have is in Julien date format so i used a function to convert it to datetime so that i am able to plot the data. However the function is now giving an error:

  1. File "D:\DATA\TEC DATA\2018\amco\RES\CMN_append_df_stack_merge.py", line 162, in jd_to_date jd = jd + 0.5

    TypeError: can only concatenate str (not "float") to str

  2. It was not giving this error before however the plot had lines crisscrossing everywhere which i suspect is because of the datetime not being formatted correctly.

Here is the code i used

*- coding: utf-8 -*-
"""
Created on Mon Mar  9 17:51:41 2020

@author: user
"""

import glob as glob

import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
from datetime import datetime
"""
Functions for converting dates to/from JD and MJD. Assumes dates are historical
dates, including the transition from the Julian calendar to the Gregorian
calendar in 1582. No support for proleptic Gregorian/Julian calendars.
:Author: Matt Davis
:Website: http://github.com/jiffyclub
"""

import math
import datetime as dt

# Note: The Python datetime module assumes an infinitely valid Gregorian calendar.
#       The Gregorian calendar took effect after 10-15-1582 and the dates 10-05 through
#       10-14-1582 never occurred. Python datetime objects will produce incorrect
#       time deltas if one date is from before 10-15-1582.

def mjd_to_jd(mjd):
    """
    Convert Modified Julian Day to Julian Day.

    Parameters
    ----------
    mjd : float
        Modified Julian Day

    Returns
    -------
    jd : float
        Julian Day


    """
    return mjd + 2400000.5


def jd_to_mjd(jd):
    """
    Convert Julian Day to Modified Julian Day

    Parameters
    ----------
    jd : float
        Julian Day

    Returns
    -------
    mjd : float
        Modified Julian Day

    """
    return jd - 2400000.5


def date_to_jd(year,month,day):
    """
    Convert a date to Julian Day.

    Algorithm from 'Practical Astronomy with your Calculator or Spreadsheet', 
        4th ed., Duffet-Smith and Zwart, 2011.

    Parameters
    ----------
    year : int
        Year as integer. Years preceding 1 A.D. should be 0 or negative.
        The year before 1 A.D. is 0, 10 B.C. is year -9.

    month : int
        Month as integer, Jan = 1, Feb. = 2, etc.

    day : float
        Day, may contain fractional part.

    Returns
    -------
    jd : float
        Julian Day

    Examples
    --------
    Convert 6 a.m., February 17, 1985 to Julian Day

    >>> date_to_jd(1985,2,17.25)
    2446113.75

    """
    if month == 1 or month == 2:
        yearp = year - 1
        monthp = month + 12
    else:
        yearp = year
        monthp = month

    # this checks where we are in relation to October 15, 1582, the beginning
    # of the Gregorian calendar.
    if ((year < 1582) or
        (year == 1582 and month < 10) or
        (year == 1582 and month == 10 and day < 15)):
        # before start of Gregorian calendar
        B = 0
    else:
        # after start of Gregorian calendar
        A = math.trunc(yearp / 100.)
        B = 2 - A + math.trunc(A / 4.)

    if yearp < 0:
        C = math.trunc((365.25 * yearp) - 0.75)
    else:
        C = math.trunc(365.25 * yearp)

    D = math.trunc(30.6001 * (monthp + 1))

    jd = B + C + D + day + 1720994.5

    return jd


def jd_to_date(jd):
    """
    Convert Julian Day to date.

    Algorithm from 'Practical Astronomy with your Calculator or Spreadsheet', 
        4th ed., Duffet-Smith and Zwart, 2011.

    Parameters
    ----------
    jd : float
        Julian Day

    Returns
    -------
    year : int
        Year as integer. Years preceding 1 A.D. should be 0 or negative.
        The year before 1 A.D. is 0, 10 B.C. is year -9.

    month : int
        Month as integer, Jan = 1, Feb. = 2, etc.

    day : float
        Day, may contain fractional part.

    Examples
    --------
    Convert Julian Day 2446113.75 to year, month, and day.

    >>> jd_to_date(2446113.75)
    (1985, 2, 17.25)

    """
    jd = jd + 0.5

    F, I = math.modf(jd)
    I = int(I)

    A = math.trunc((I - 1867216.25)/36524.25)

    if I > 2299160:
        B = I + 1 + A - math.trunc(A / 4.)
    else:
        B = I

    C = B + 1524

    D = math.trunc((C - 122.1) / 365.25)

    E = math.trunc(365.25 * D)

    G = math.trunc((C - E) / 30.6001)

    day = C - E + F - math.trunc(30.6001 * G)

    if G < 13.5:
        month = G - 1
    else:
        month = G - 13

    if month > 2.5:
        year = D - 4716
    else:
        year = D - 4715

    return year, month, day


def hmsm_to_days(hour=0,min=0,sec=0,micro=0):
    """
    Convert hours, minutes, seconds, and microseconds to fractional days.

    Parameters
    ----------
    hour : int, optional
        Hour number. Defaults to 0.

    min : int, optional
        Minute number. Defaults to 0.

    sec : int, optional
        Second number. Defaults to 0.

    micro : int, optional
        Microsecond number. Defaults to 0.

    Returns
    -------
    days : float
        Fractional days.

    Examples
    --------
    >>> hmsm_to_days(hour=6)
    0.25

    """
    days = sec + (micro / 1.e6)

    days = min + (days / 60.)

    days = hour + (days / 60.)

    return days / 24.


def days_to_hmsm(days):
    """
    Convert fractional days to hours, minutes, seconds, and microseconds.
    Precision beyond microseconds is rounded to the nearest microsecond.

    Parameters
    ----------
    days : float
        A fractional number of days. Must be less than 1.

    Returns
    -------
    hour : int
        Hour number.

    min : int
        Minute number.

    sec : int
        Second number.

    micro : int
        Microsecond number.

    Raises
    ------
    ValueError
        If `days` is >= 1.

    Examples
    --------
    >>> days_to_hmsm(0.1)
    (2, 24, 0, 0)

    """
    hours = days * 24.
    hours, hour = math.modf(hours)

    mins = hours * 60.
    mins, min = math.modf(mins)

    secs = mins * 60.
    secs, sec = math.modf(secs)

    micro = round(secs * 1.e6)

    return int(hour), int(min), int(sec), int(micro)


def datetime_to_jd(date):
    """
    Convert a `datetime.datetime` object to Julian Day.

    Parameters
    ----------
    date : `datetime.datetime` instance

    Returns
    -------
    jd : float
        Julian day.

    Examples
    --------
    >>> d = datetime.datetime(1985,2,17,6)  
    >>> d
    datetime.datetime(1985, 2, 17, 6, 0)
    >>> jdutil.datetime_to_jd(d)
    2446113.75

    """
    days = date.day + hmsm_to_days(date.hour,date.minute,date.second,date.microsecond)

    return date_to_jd(date.year,date.month,days)


def jd_to_datetime(jd):
    """
    Convert a Julian Day to an `jdutil.datetime` object.

    Parameters
    ----------
    jd : float
        Julian day.

    Returns
    -------
    dt : `jdutil.datetime` object
        `jdutil.datetime` equivalent of Julian day.

   --------------------------------------------------------------------------------------------------------------------------------------------------------
    --------
    >>> jd_to_datetime(2446113.75)
    datetime(1985, 2, 17, 6, 0)

    """
    year, month, day = jd_to_date(jd)

    frac_days,day = math.modf(day)
    day = int(day)

    hour,min,sec,micro = days_to_hmsm(frac_days)

    return datetime(year,month,day,hour,min,sec,micro)
'---------------------------------------------------------------------------------------------------------------------------------------------------------'
#"""
#JULIEN DAY CONVERTOR 

def timedelta_to_days(td):
    """
    Convert a `datetime.timedelta` object to a total number of days.

    Parameters
    ----------
    td : `datetime.timedelta` instance

    Returns
    -------
    days : float
        Total number of days in the `datetime.timedelta` object.

    Examples
    --------
    >>> td = datetime.timedelta(4.5)
    >>> td
    datetime.timedelta(4, 43200)
    >>> timedelta_to_days(td)
    4.5

    """
    seconds_in_day = 24. * 3600.

    days = td.days + (td.seconds + (td.microseconds * 10.e6)) / seconds_in_day

    return days


"""
'Start of main Code---------------------------------------------------------------------------------------------------'
"""
#########

files = glob.glob('*.Cmn')
files = files[:3]
np_array_values = []
for file in files:
         df = pd.read_csv(file, delimiter="\t", sep ='\t', skiprows = 5, names = ["Jdate" ,'Time' ,'PRN' ,'Az','Ele','Lat', 'Lon' ,'Stec', 'Vtec', 'S4'])
         #df.set_index('Jdatet')



         #Appending Each data frame to the Giant Array np-array_values then Stacking them

         np_array_values.append (df)
         merge_values = np.vstack(np_array_values) 
         #CONVERTING ARRAY BACK TO DATA FRAME
         TEC_data =pd.DataFrame(merge_values)

         Vtec = TEC_data.loc[:,7]
         jdate = TEC_data.loc[:,0]
         Time = TEC_data.loc[:,1]
         month_name = file[:-7]
         STATION_NAME = file[:4]

df.replace('-99.000', np.nan)

pos= np.where(jdate)[0]
pos1=np.where(Vtec[pos]<-20.0)[0]

Vtec[pos1]='nan'

fulldate = []

for i in jdate:
    a = jd_to_datetime(i)
    fulldate.append(a)


         #plt.show()
plt.plot(fulldate,  Vtec)
 #plt.xlim(0, 24)
 #plt.xticks(np.arange(0, 26, 2))
plt.ylabel('TECU')

plt.xlabel('Date')

#plt.grid(axis='both')

plt.title("Station : " + STATION_NAME.upper())t**

Here is a sample of the data. The columns are: nknown_station, "E:\GPS DATA\rbmc2\2018\amco\amco3630.18o"

-4.87199 294.66602 75.87480

Jdatet Time PRN Az Ele Lat Lon Stec Vtec S4 2458481.500000 -24.000000 1 198.34 23.37 -10.70 292.70 18.62 13.70 -99.000 2458481.500347 0.008333 1 198.21 23.53 -10.67 292.73 18.65 13.74 -99.000 2458481.500694 0.016667 1 198.07 23.69 -10.64 292.75 18.76 13.84 -99.000 2458481.501042 0.025000 1 197.94 23.85 -10.61 292.78 18.68 13.83 -99.000 2458481.501389 0.033333 1 197.81 24.01 -10.58 292.80 18.60 13.83 -99.000 2458481.501736 0.041667 1 197.68 24.17 -10.55 292.83 18.53 13.83 -99.000 2458481.502083 0.050000 1 197.54 24.33 -10.52 292.85 18.53 13.86 -99.000 2458481.502431 0.058333 1 197.41 24.49 -10.49 292.88 18.51 13.88 -99.000 2458481.502778 0.066667 1 197.28 24.65 -10.46 292.90 18.66 14.00 -99.000 2458481.503125 0.075000 1 197.15 24.81 -10.43 292.92 18.78 14.09 0.238 2458481.503472 0.083333 1 197.01 24.98 -10.40 292.95 18.55 14.01 -99.000 2458481.503819 0.091667 1 196.88 25.14 -10.37 292.97 18.39 13.96 -99.000 2458481.504167 0.100000 1 196.75 25.30 -10.34 292.99 18.33 13.97 -99.000 2458481.504514 0.108333 1 196.62 25.47 -10.31 293.02 18.20 13.94 -99.000 2458481.504861 0.116667 1 196.49 25.63 -10.28 293.04 17.61 13.67 -99.000 2458481.505208 0.125000 1 196.36 25.80 -10.25 293.06 16.74 13.25 -99.000 2458481.505556 0.133333 1 196.23 25.96 -10.22 293.09 16.06 12.92 -99.000 2458481.505903 0.141667 1 196.10 26.13 -10.19 293.11 15.46 12.64 -99.000 2458481.506250 0.150000 1 195.97 26.30 -10.16 293.13 14.77 12.31 -99.000 2458481.506597 0.158333 1 195.84 26.46 -10.13 293.15 14.42 12.15 0.127


Solution

  • I added the list containing the datetime object as an index to the Dataframe. Apparently to recognise the dates you need to use plt.plot_date instead of using the usual plt.plot. Datetime index will take care of the rest.