library(psychrolib)
SetUnitSystem("SI")


# Step 1: Create a dataframe with hourly temperature data for a full year

# Let's assume the dataframe is named "temp_data" and has two columns: "datetime" (with hourly timestamps) and "temperature" (with hourly temperature readings).
df <- read.csv("data.csv",col.names=c("DateTime", "DB", "RH","DP","P","WIND","WIND-DIR","SUN","RAIN","START"))

df=df[complete.cases(df), ]


date_format <- "%d/%m/%Y %H:%M"
df$DateTime=as.POSIXct(df$DateTime, format = date_format)
df$DB=(df$DB/10)
df$RH=(df$RH/100)
df$P=(df$P*100)
df$WB=GetTWetBulbFromRelHum(df$DB, df$RH, df$P)
# create a sample data frame

head(df,n=5)
nrow(df)
df=df[format(df$DateTime, "%m") %in% c("06", "07", "08", "09"),]
nrow(df)
total_hours=nrow(df)
# use the aggregate function to calculate the average of y values for each x value
print('hello')
length(unique(df$DB))
df2=aggregate(WB ~ DB, data = df, FUN = mean)
df2$count <- sapply(df2$DB, function(x) sum(df$DB > x))
#df2$count <- sapply(df2$DB, function(x) sum(df$DB > x))

#df2$hours=sum(df$DB > df2$DB)
df=df2

df$pct_exceeded=(df$count/total_hours)*100

head(df)
tail(df)
# Step 1: Create a dataframe with hourly temperature and relative humidity data for a full year
# Let's assume the dataframe is named "df" and has three columns: "datetime" (with hourly timestamps), "DB" (with hourly dry bulb temperature readings), and "RH" (with hourly relative humidity readings).
# Step 2: Calculate the total number of hours in a year
# Step 2: Calculate the total number of hours in a year

#print(total_hours)
#hrs=nrow(unique(df[,c("DB","WB")]))
# Step 3: Create a new dataframe named "temp_pct" with four columns: "DB", "RH", "pct_exceeded", and "RH_pct_exceeded"
#temp_pct <- data.frame(DB = df$DB, WB = df$WB, pct_exceeded = numeric(length(df$DB)))

# Step 4: For each temperature value in the "DB" column of the "df" dataframe, calculate the percentage of time that temperature is exceeded and the relative humidity for the first occurrence of that temperature, and store them in the corresponding rows of the "temp_pct" dataframe
#for (i in seq_along(df$DB)) {
#  temp <- df$DB[i]
#  hours_exceeded <- df$count
#  pct_exceeded <- hours_exceeded / total_hours
#  rh <- df$WB[df$DB == temp][1]
#  temp_pct$pct_exceeded[i] <- pct_exceeded
#  temp_pct$RH_pct_exceeded[i] <- rh
#}

# Step 5: Calculate the average relative humidity for each unique value of "pct_exceeded" and store it in the "RH_pct_exceeded" column of the "temp_pct" dataframe
#for (pct in unique(temp_pct$pct_exceeded)) {
#  avg_rh <- mean(temp_pct$RH_pct_exceeded[temp_pct$pct_exceeded == pct])
#  temp_pct$RH_pct_exceeded[temp_pct$pct_exceeded == pct] <- avg_rh
#}

# Step 6: Order the "temp_pct" dataframe by the "pct_exceeded" column, with lowest percentage first
#temp_pct <- temp_pct[order(temp_pct$pct_exceeded), ]

# Step 7: Print the first 20 rows of the "temp_pct" dataframe
#head(temp_pct,5)
#tail(temp_pct,5)
#temp_pct=head(temp_pct, 20)
df <- subset(df, pct_exceeded <= 1)

#print(temp_pct)

library(ggplot2)
library(ggpsychro)


p = ggpsychro(tdb_lim = c(-15, 40)) +
    geom_grid_relhum(alpha = 1.0, label.alpha = 1.0, label.size = 6, label.fontface = 2) +
    scale_relhum(minor_breaks = NULL) +
    geom_grid_wetbulb(size = 1.0, color = "black", alpha = 1.0, label_loc = NA) +
    scale_wetbulb(breaks = seq(25, 30, by = 5), minor_breaks = NULL) +
    geom_grid_vappres(label.size = 5) +
    scale_vappres(breaks = seq(6000, 7000, by = 500), limits = c(6000, 7000)) +
    geom_grid_specvol() +
    scale_specvol(labels = NULL) +
#    geom_line(aes(x = 20:30,  = 64000, y=1,color="orange"), stat = "wetbulb") +
    geom_point(aes(df$DB, wetbulb=df$WB),stat = "wetbulb")+
    geom_grid_enthalpy()
ggsave("chart1.pdf", plot = p, width = 6, height = 4)


# Step 8: Create a list of unique combinations of "RH_pct_exceeded" and "pct_exceeded"
#combo_list <- unique(temp_pct[, c("RH_pct_exceeded", "pct_exceeded")])