# NS41 Satellite Data # Michael D Cayton # December 22, 2018 # Last update 12-23-18 ### sum of col V81 through 99 ### 21 not 240 remove # number of rows, sum # midified Bessal function of the second kind. #### Load the libraries ### library(tidyverse) library(RCurl) library(lubridate) # Set the url for the ftp site working directory. # This is where the origianal data is stored. url <- 'ftp://ftp.ngdc.noaa.gov/STP/space-weather/satellite-data/satellite-systems/gps/data/ns41/' # Create a vector of file names from the directory # This will allow us to extract the dates and will be used to create a master file. filename <- getURL(url, dirlistonly = TRUE) filename # Clean up the file list destnames <- strsplit(filename, "\r*\n")[[1]] destnames # Create data frame for the data. ns41 = NULL ns41 <- tbl_df(ns41) # For loop to access each file, combine into one dataframe, and add a column for the filename date l <- length(destnames) data <- NULL i <- NULL for (i in 1:l) { data <- read.table(file = paste(url, destnames[i], sep="")) data$Date <- as.Date(str_sub(destnames[i], start=6, end=11), "%y%m%d") ns41 <- rbind(ns41, data) print(paste(i, "of", l, "complete")) } # Save table. write_csv(ns41, path="/Users/mdcay/Documents/Happy Dad/ns41.csv") #Import raw data as needed. ns41 <- read_csv(file="/Users/mdcay/Documents/Happy Dad/ns41.csv") ns41 <- tbl_df(ns41) # Select for V21 = 240 ns41 <- ns41 %>% filter(V21 == 240) # Create table, grouped by unique week, and sum of V81 through V99 ns41_sums <- ns41 %>% group_by(Date) %>% select(Date,V25, V82:V100) %>% summarise(n=n(), V24=sum(V25), V81=sum(V82), V82=sum(V83), V83=sum(V84), V84=sum(V85), V85=sum(V86), V86=sum(V87), V87=sum(V88), V88=sum(V89), V89=sum(V90), V90=sum(V91), V91=sum(V92), V92=sum(V93), V93=sum(V94), V94=sum(V95), V95=sum(V96), V96=sum(V97), V97=sum(V98), V98=sum(V99), V99=sum(V100) ) write_csv(ns41_sums, path='ns41_sums.csv') ns41_sums <- read_csv(file="/Users/mcayton/Documents/Happy Dad/ns41_sums.csv") # Create a list of the energy levels. These can be found in V100 through V130. energies <- c(.6, .7, .8, .9, 1, 1.25, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, 5, 6, 7, 8, 9, 10) # Use Modified Bessel function of the second kind. Each V81-V99 row is an interval and go with the energy list. x <- ns41_sums %>% select(-Date, -n) nu <- as.numeric(x[1,]) bessel <- besselK(nu, energies) bessel ggplot() + geom_line(aes(nu, energies)) ### sum of col V81 through 99 ### 21 not 240 remove # number of rows, sum # midified Bessal function of the second kind. ### sum of col V81 through 99 ### 21 not 240 remove # number of rows, sum ### JSON with table info ### # test <- fromJSON('{ # "decimal_day": { "DESCRIPTION":"Decimal Day of year", # "DIMENSION": [1], # "UNITS": "days", # "START_COLUMN": 0 # }, # # "Geographic_Latitude": { "DESCRIPTION":"Geographic_Latitude", # "DIMENSION": [1], # "UNITS": "degrees", # "START_COLUMN": 1 # }, # "Geographic_Longitude": { "DESCRIPTION":"Geographic_Longitude", # "DIMENSION": [1], # "UNITS": "degrees", # "START_COLUMN": 2 # }, # "Rad_Re": { "DESCRIPTION":"Radius (earth radii)", # "DIMENSION": [1], # "UNITS": "R_E", # "START_COLUMN": 3 # }, # "rate_electron_measured": { "DESCRIPTION":"rate_electron_measured", # "DIMENSION": [8], # "UNITS": "hertz", # "START_COLUMN": 4 # }, # "rate_proton_measured": { "DESCRIPTION":"rate_proton_measured", # "DIMENSION": [8], # "UNITS": "hertz", # "START_COLUMN": 12 # }, # "collection_interval": { "DESCRIPTION":"collection_interval", # "DIMENSION": [1], # "UNITS": "seconds", # "START_COLUMN": 20 # }, # "year": { "DESCRIPTION":"year eg 2015", # "DIMENSION": [1], # "UNITS": "years", # "START_COLUMN": 21 # }, # "decimal_year": { "DESCRIPTION":"decimal_year contains fractional part of year", # "DIMENSION": [1], # "UNITS": "years", # "START_COLUMN": 22 # }, # "svn_number": { "DESCRIPTION":"SVN number", # "DIMENSION": [1], # "UNITS": "none", # "START_COLUMN": 23 # }, # "dropped_data": { "DESCRIPTION":"dropped_data=1 means ignore data", # "DIMENSION": [1], # "UNITS": "none", # "START_COLUMN": 24 # }, # "b_coord_radius": { "DESCRIPTION":"radius from earths dipole axis", # "DIMENSION": [1], # "UNITS": "R_E", # "START_COLUMN": 25 # }, # "b_coord_height": { "DESCRIPTION":"height above earths equatorial plane", # "DIMENSION": [1], # "UNITS": "R_E", # "START_COLUMN": 26 # }, # "magnetic_longitude": { "DESCRIPTION":"magnetic longitude", # "DIMENSION": [1], # "UNITS": "degrees", # "START_COLUMN": 27 # }, # "L_shell": { "DESCRIPTION":"L shell, dipole file, T89", # "DIMENSION": [1], # "UNITS": "R_E", # "START_COLUMN": 28 # }, # "bfield_ratio": { "DESCRIPTION":"Bsatellite/Bequator", # "DIMENSION": [1], # "UNITS": "none", # "START_COLUMN": 29 # }, # "local_time": { "DESCRIPTION":"Magnetic local time", # "DIMENSION": [1], # "UNITS": "h", # "START_COLUMN": 30 # }, # "b_sattelite": { "DESCRIPTION":" field at satellite", # "DIMENSION": [1], # "UNITS": "gauss", # "START_COLUMN": 31 # }, # "b_equator": { "DESCRIPTION":"B field at equator", # "DIMENSION": [1], # "UNITS": "gauss", # "START_COLUMN": 32 # }, # "diffp": { "DESCRIPTION":"not sure what this is", # "DIMENSION": [1], # "UNITS": "gauss", # "START_COLUMN": 33 # }, # "sigmap": { "DESCRIPTION":"not sure what this is", # "DIMENSION": [1], # "UNITS": "gauss", # "START_COLUMN": 34 # }, # "electron_background": { "DESCRIPTION":"electron_background estimated", # "DIMENSION": [8], # "UNITS": "hertz", # "START_COLUMN": 35 # }, # "proton_background": { "DESCRIPTION":"proton_background estimated", # "DIMENSION": [8], # "UNITS": "hertz", # "START_COLUMN": 43 # }, # "proton_activity": { "DESCRIPTION":"proton_activity=1 if proton activity is occurring", # "DIMENSION": [1], # "UNITS": "none", # "START_COLUMN": 51 # }, # "electron_temperature_fit": { "DESCRIPTION":"electron temperature from Maxwellian fit", # "DIMENSION": [1], # "UNITS": "MeV", # "START_COLUMN": 52 # }, # "electron_density_fit": { "DESCRIPTION":"electron number density fit", # "DIMENSION": [1], # "UNITS": "cm^-3", # "START_COLUMN": 53 # }, # "model_counts_electron_fit": { "DESCRIPTION":"E1-E8 rates from Maxwellian fit", # "DIMENSION": [8], # "UNITS": "hertz", # "START_COLUMN": 54 # }, # "dtc_counts_electron": { "DESCRIPTION":"E1-E8 rates dead time corrected", # "DIMENSION": [8], # "UNITS": "hertz", # "START_COLUMN": 62 # }, # "integral_flux_instrument": { "DESCRIPTION":"integral of electron flux fit above integral_flux_energy[i]", # "DIMENSION": [30], # "UNITS": "cm^-2sec^-1sr^-1", # "START_COLUMN": 70 # }, # "integral_flux_energy": { "DESCRIPTION":"energies for integral_flux_instrument", # "DIMENSION": [30], # "UNITS": "MeV", # "START_COLUMN": 100 # }, # "electron_diff_flux_energy": { "DESCRIPTION":"energies for the fluxes in electron_diff_flux_energy", # "DIMENSION": [15], # "UNITS": "MeV", # "START_COLUMN": 130 # }, # "electron_diff_flux": { "DESCRIPTION":"electron flux at energies electron_diff_flux[i]", # "DIMENSION": [15], # "UNITS": "cm^-2sec^-1sr^-1MeV^-1", # "START_COLUMN": 145 # } # }')