Calculate Daily Averages / Sums from hourly data
Jump to navigation
Jump to search
Copy and paste the text in the box into the visual basic editor that comes with Excel. To get access to personal.xls in the visual basic editor check out this article: Save Macros in Excel
This macro and collection of subroutines creates daily averages from hourly data. It was written originally for data suitable for use with MicroMet but wider application are totally possible. What's handy about this collection of macros is that Air temperature is averaged, precipitation is totalized, wind speed is converted to vector format and back to radial format and relative humidity is averaged correctly, too.
Any questions check with Bob.
Option Explicit
Private Sub daily_avg_sum_etc()
' worksheet variables
Dim inputrow
Dim outputrow
Dim inputWKS
Dim outputWKS
Dim tempCOL
Dim rhCOl
Dim wsCOL
Dim wdCOL
Dim precipCOL
' date and time variables
Dim inputdate
Dim olddate
Dim outputYEAR
Dim outputDAY
Dim outputMONTH
Dim outputTIME
' position variables
Dim easting
Dim northing
Dim stationID
Dim ELevation
' environmental input variables
Dim inputWS As Double
Dim inputWD As Double
Dim inputAT As Double
Dim inputRH As Double
Dim inputPRECIP As Double
' environmental sums for daily means
Dim sumAT As Double
Dim sumVP As Double
Dim sumxWS As Double
Dim sumyWS As Double
Dim sumPRECIP As Double
Dim countAT As Double
Dim countVP As Double
Dim countWS As Double
Dim countPRECIP As Double
Dim zoo As Double
' environmental output variables
Dim outputWS As Double
Dim outputWD As Double
Dim outputAT As Double
Dim outputVP As Double
Dim outputRH As Double
Dim outputPRECIP As Double
' things to do:
' 1) split date string
' 2) add an arbitrary hour to worksheet
' 3) add station ID, easting, & northing to worksheet
' 4) convert station data to metric and output to worksheet
' 5) convert dewpoint & temperature to RH
'''''''''''''''''''''''''''''''''''''''''''''
things to change before each station process
'''''''''''''''''''''''''''''''''''''''''''''
easting = 516394
northing = 7256371
northing = northing - 7000000 ' for micromet format
ELevation = 95
stationID = 102
outputrow = 4
inputrow = 1
''''''''''''''''''''''''
set worksheet variables
''''''''''''''''''''''''
inputWKS = 2
outputWKS = 1
tempCOL = 7
rhCOl = 10
wsCOL = 14
wdCOL = 15
precipCOL = 17
'''''''''''''''''''''''''''''''
Start running through the data
'''''''''''''''''''''''''''''''
olddate = Format(Worksheets(inputWKS).Cells(inputrow, 1).Value, "mm/dd/yyyy")
Do While Worksheets(inputWKS).Cells(inputrow, 1).Value <> Empty
inputdate = Format(Worksheets(inputWKS).Cells(inputrow, 1).Value, "mm/dd/yyyy")
If olddate <> inputdate Then
'''''''''''''''''''''''''''''''''''''''''''''''''''''''''
do final formatting and calculations prior to outputting
'''''''''''''''''''''''''''''''''''''''''''''''''''''''''
outputYEAR = Right(CStr(olddate), 4)
outputMONTH = Left(CStr(olddate), 2)
outputDAY = Mid(CStr(olddate), 4, 2)
outputTIME = 1
If countAT <> 0 Then outputAT = sensor_avg(sumAT, countAT) Else outputAT = -9999
If countVP <> 0 Then
outputVP = sensor_avg(sumVP, countVP)
outputRH = VaporpressuretoRH(outputVP, outputAT) * 100
If outputRH > 100 Then outputRH = 100
Else:
outputRH = -9999
End If
' note to self: adjust WS to 10m before outputtting
' also, use two functions, one to convert ux, uy into WS and one to convert into WD.
' if there is only WS then output only WS also.
If countWS <> 0 And sumxWS <> 0 And sumyWS <> 0 Then
' there's stuff to do.
sumxWS = sumxWS / countWS
sumyWS = sumyWS / countWS
outputWS = uxy_to_WS(sumxWS, sumyWS)
If outputMONTH < 6 Or outputMONTH > 8 Then
zoo = 0.01
Else
zoo = 0.03
End If
outputWS = elevateWS(outputWS, 3, 10, zoo)
outputWD = uxy_to_WD(sumxWS, sumyWS)
Else
outputWS = -9999
outputWD = -9999
End If
If countPRECIP <> 0 Then outputPRECIP = sumPRECIP Else outputPRECIP = -9999
''''''''''''''''''''''''''
dump stuff to spreadsheet
''''''''''''''''''''''''''
Worksheets(outputWKS).Cells(outputrow, 1).Value = outputYEAR
Worksheets(outputWKS).Cells(outputrow, 2).Value = outputMONTH
Worksheets(outputWKS).Cells(outputrow, 3).Value = outputDAY
Worksheets(outputWKS).Cells(outputrow, 4).Value = outputTIME
Worksheets(outputWKS).Cells(outputrow, 5).Value = stationID
Worksheets(outputWKS).Cells(outputrow, 6).Value = easting
Worksheets(outputWKS).Cells(outputrow, 7).Value = northing
Worksheets(outputWKS).Cells(outputrow, 8).Value = ELevation
Worksheets(outputWKS).Cells(outputrow, 9).Value = outputAT
Worksheets(outputWKS).Cells(outputrow, 10).Value = outputRH
Worksheets(outputWKS).Cells(outputrow, 11).Value = outputWS
Worksheets(outputWKS).Cells(outputrow, 12).Value = outputWD
Worksheets(outputWKS).Cells(outputrow, 13).Value = outputPRECIP
'''''''''''''''''''''''''''''''
reset values and increment row
'''''''''''''''''''''''''''''''
sumAT = 0
sumVP = 0
sumxWS = 0
sumyWS = 0
outputPRECIP = 0
countAT = 0
countVP = 0
countWS = 0
countPRECIP = 0
outputrow = outputrow + 1
olddate = inputdate
End If
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
okay with the outputting out of the way we're ready to start reading in the data.
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
If tempCOL <> 6999 Then inputAT = Worksheets(inputWKS).Cells(inputrow, tempCOL).Value Else inputAT = -9999
If rhCOl <> 6999 Then inputRH = Worksheets(inputWKS).Cells(inputrow, rhCOl).Value Else inputRH = -9999
If wsCOL <> 6999 Then inputWS = Worksheets(inputWKS).Cells(inputrow, wsCOL).Value Else inputWS = -9999
If wdCOL <> 6999 Then inputWD = Worksheets(inputWKS).Cells(inputrow, wdCOL).Value Else inputWD = -9999
If precipCOL <> 6999 Then inputPRECIP = Worksheets(inputWKS).Cells(inputrow, precipCOL).Value Else inputPRECIP = -9999
If Abs(inputAT) <> 6999 And inputAT <> -9999 Then
sumAT = sumAT + inputAT
countAT = countAT + 1
End If
If Abs(inputRH) <> 6999 And inputRH <> -9999 And Abs(inputAT) <> 6999 And inputAT <> -9999 Then
sumVP = sumVP + RHtoVaporpressure(inputRH, inputAT)
countVP = countVP + 1
End If
If Abs(inputWS) <> 6999 And inputWS <> -9999 And Abs(inputWD) <> 6999 And inputWD <> -9999 Then
sumxWS = sumxWS + WS_to_ux(inputWS, inputWD)
sumyWS = sumyWS + WS_to_uy(inputWS, inputWD)
countWS = countWS + 1
End If
If Abs(inputPRECIP) <> 6999 And inputPRECIP <> -9999 Then
sumPRECIP = sumPRECIP + inputPRECIP
countPRECIP = countPRECIP + 1
End If
' next time step
inputrow = inputrow + 1
Loop
End Sub
Public Function RHtoVaporpressure(inputRH As Double, inputAT As Double) ' Equations come from Dingman's Physical Hydrology second edition. ' pp 586 - 587 ' eq D-7: ' e* = 0.611 * exp( (17.3 * T) / ( T + 237.3 ) ) ' e* [kPa], T [Deg C] ' eq D-10: ' Wa = ea / e*a ' Wa = Relative Humidity [%] ' ea = Vapor Pressure [kPa] ' e*a = Saturation Vapor Pressure [kPa] Dim ea As Double Dim e_star_a As Double e_star_a = 0.611 * Exp((17.3 * inputAT) / (inputAT + 237.3)) ea = e_star_a * inputRH / 100 RHtoVaporpressure = ea End Function
Private Function VaporpressuretoRH(VaporPressure As Double, AirTemp As Double) ' Equations come from Dingman's Physical Hydrology second edition. ' pp 586 - 587 ' eq D-7: ' e* = 0.611 * exp( (17.3 * T) / ( T + 237.3 ) ) ' e* [kPa], T [Deg C] ' eq D-10: ' Wa = ea / e*a ' Wa = Relative Humidity [%] ' ea = Vapor Pressure [kPa] ' e*a = Saturation Vapor Pressure [kPa] Dim ea As Double Dim e_star_a As Double e_star_a = 0.611 * Exp((17.3 * AirTemp) / (AirTemp + 237.3)) VaporpressuretoRH = VaporPressure / e_star_a End Function
Public Function sensor_avg(total As Double, count As Double)
' this function just returns an average
If count <> 0 Then
sensor_avg = total / count
Else:
total = -9999
End If
End Function
Public Function WS_to_ux(inputWS, inputWD) ' assumes WD is in degrees. Dim ux ux = inputWS * Cos(inputWD * 3.14592 / 180) WS_to_ux = ux End Function
Public Function WS_to_uy(inputWS, inputWD) ' assumes WD is in degrees. Dim uy uy = inputWS * Sin(inputWD * 3.141592 / 180) WS_to_uy = uy End Function
Public Function uxy_to_WS(inputux, inputuy) ' go from vector components back to radial. uxy_to_WS = (inputux ^ 2 + inputuy ^ 2) ^ 0.5 End Function
Public Function uxy_to_WD(inputux, inputuy) Dim WD ' go from vector components back to radial WD = Atn(inputuy / inputux) * 180 / 3.141592 ' now dump it in the right quandrant.. If inputux > 0 And inputuy > 0 Then WD = WD + 0 If inputux < 0 And inputuy > 0 Then WD = WD + 180 If inputux < 0 And inputuy < 0 Then WD = WD + 180 If inputux > 0 And inputuy < 0 Then WD = WD + 360 uxy_to_WD = WD End Function
Public Function elevateWS(inputWS, wsheight, reportheight, zoo) ' U(h) = U(H) * ( ln(h/z0) / ln(H/z0) ) ' h = target WS height in meters ' H = anemometer height in meters ' z0 = 0.01 September - May ' z0 = 0.03 June - August elevateWS = inputWS * (Log(reportheight / zoo) / Log(wsheight / zoo)) End Function
