Skip to content

GitLab

D
DataEvaluationDWD
Commit 45ad43de authored by Gerrit Erichsen's avatar Gerrit Erichsen
Browse files
Options

BugFixes: The stations did not get cleared and removed correctly, so it was... 

BugFixes: The stations did not get cleared and removed correctly, so it was reset to always start at empty each loop
parent fcb7c0ae
Hide whitespace changes
Inline Side-by-side
Showing with 93 additions and 63 deletions
Code/evaluation.py
View file @ 45ad43de
......@@ -113,15 +113,15 @@ def findStationId(stations, stationId):
return -1
if __name__ == '__main__':
dataType = 'wind' #either solar, temperature or wind
dataType = 'solar' #either solar, temperature or wind
dataTypeShort = ''
dataTypeShort2 = ''
# while for map creation the loop for each year in range was quite
# successful, it had its problems with the statistics part creating
# accurate plots so for the statistic plots,
# the limits were set to one year only
startYear = 2015
endYear = 2015
startYear = 2011
endYear = 2017
path = 'D:/WetterdatenCDC'
path2 = 'D:/WetterdatenPamore'
pathFigs = 'J:/wka/user/erichsen/DataDocumentation/images/'
......@@ -159,12 +159,11 @@ if __name__ == '__main__':
fileNameH5 = dwd.getFullFileNameH5(dataType, year, path2)
fileName2H5 = dwd.getFullFileNameH5(dataType, year, path2, 2)
if yearAsInt == startYear:
lat,lon = dwd.readGrid(fileNameH5)
print("Done reading grid.")
lat,lon = dwd.readGrid(fileNameH5)
print("Done reading grid.")
stations = [] # first value of a WeatherStation (dependent on type)
stations2 = [] # second value of WeatherStation (dependent on type)
stations = [] # first value of a WeatherStation (dependent on type)
stations2 = [] # second value of WeatherStation (dependent on type)
#metaData = metaData[0:10]
......@@ -172,14 +171,14 @@ if __name__ == '__main__':
currentIdentifiers = []
for station in metaData:
currentIdentifiers.append(station[0])
deleteCounter = 0
for i in range(0, len(stations)):
if not stations[i - deleteCounter].getIdentifier() \
lengthBefore = len(stations)
for i in range(1, lengthBefore):
if not stations[lengthBefore - i].getIdentifier() \
in currentIdentifiers:
stations.remove(stations[i - deleteCounter])
if i < len(stations2):
stations2.remove(stations2[i - deleteCounter])
deleteCounter += 1
stations.remove(stations[lengthBefore - i])
if (lengthBefore - i) < len(stations2):
stations2.remove(stations2[lengthBefore - i])
i = 1
total = len(metaData)
......@@ -187,20 +186,27 @@ if __name__ == '__main__':
for station in metaData:
fileName = dwd.getFullFileName(dataType, station[0],\
station[4], station[5], path)
ind = findStationId(stations, station[0])
if fileName == '':
stations.remove(stations[ind])
if ind >= len(stations2):
stations2.remove(stations2[ind])
print("Done", i, "/", total)
i+=1
continue
ind = findStationId(stations, station[0])
if (ind < 0):
stations.append(WeatherStation(station[1], station[0], \
station[2], station[3], \
dataType))
else :
stations[ind].clear()
if len(stations2) > 0:
stations2[ind].clear()
stations[ind].readInCdcFile(fileName, year, cdcColumn)
if (len(stations[ind].getMeasuredData()) == 0):
stations.remove(stations[ind])
if ind >= len(stations2):
stations2.remove(stations2[ind])
print("Done", i, "/", total)
i+=1
continue
......@@ -211,7 +217,8 @@ if __name__ == '__main__':
stations[ind].getLongitude())
stations[ind].setLatLonIds(bLat, bLon)
if dataType == 'solar' or dataType == 'temperature':
stations[ind].readInModelFile(fileNameH5, dataTypeShort)
stations[ind].readInModelFile(fileNameH5, dataTypeShort, \
yearAsInt)
else:
height = dwd.findHeight(fileName, \
stations[ind].getIdentifier(), \
......@@ -221,10 +228,19 @@ if __name__ == '__main__':
fileNameH5u = fileNameH5 + '_' + upper + '.h5'
fileNameH5l = fileNameH5 + '_' + lower + '.h5'
stations[ind].readInModelFile(fileNameH5l, dataTypeShort, \
fileNameH5u, height)
yearAsInt, fileNameH5u, \
height)
else:
fileNameH5l = fileNameH5 + '_' + lower + '.h5'
stations[ind].readInModelFile(fileNameH5l, dataTypeShort)
stations[ind].readInModelFile(fileNameH5l, dataTypeShort, \
yearAsInt)
if len(stations[ind].modelData) == 0:
stations.remove(stations[ind])
if ind >= len(stations2):
stations2.remove(stations2[ind])
print("Done", i, "/", total)
i+=1
continue
if dataType == 'solar' or dataType == 'wind':
secondType = ''
if (dataType == 'solar'):
......@@ -238,14 +254,16 @@ if __name__ == '__main__':
station[2], station[3], \
secondType))
stations2[ind].readInCdcFile(fileName, year, cdcColumn + 1)
if (len(stations2[ind].m_measuredData) == 0):
stations2.pop()
if (len(stations2[ind].observationData) == 0):
stations.remove(stations[ind])
tations2.remove(stations2[ind])
print("Done", i, "/", total)
i+=1
continue
if dataType == 'solar':
stations2[ind].readInModelFile(fileName2H5, \
dataTypeShort2)
dataTypeShort2,\
yearAsInt)
else:
height = dwd.findHeight(fileName, \
stations[ind].getIdentifier(), \
......@@ -255,13 +273,21 @@ if __name__ == '__main__':
fileNameH5u = fileName2H5 + '_' + upper + '.h5'
fileNameH5l = fileName2H5 + '_' + lower + '.h5'
stations2[ind].readInModelFile(fileNameH5l, \
dataTypeShort2, \
dataTypeShort2,\
yearAsInt, \
fileNameH5u, \
height)
else:
fileNameH5l = fileName2H5 + '_' + lower + '.h5'
stations2[ind].readInModelFile(fileNameH5l, \
dataTypeShort2)
dataTypeShort2, \
yearAsInt)
if (len(stations2[ind].modelData) == 0):
stations.remove(stations[ind])
stations2.remove(stations2[ind])
print("Done", i, "/", total)
i+=1
continue
print("Done read", i, "/", total)
i+=1
......@@ -342,7 +368,7 @@ if __name__ == '__main__':
elif dataType == 'solar':
plt_stats1.set_ylabel('W / m$^2$')
else:
plt1.set_ylabel('m / s')
plt_stats1.set_ylabel('m / s')
plt_stats1.annotate(annotation, (0.6, 0.95 * np.max(rmse)), \
fontsize=8)
plt_stat2.boxplot(mae, labels=xLabelBox1)
......@@ -353,7 +379,7 @@ if __name__ == '__main__':
figStat.subplots_adjust(wspace=0.3, hspace=0.3)
figStat.savefig(figFileName, dpi=600, format='png')
plt.close()
fig1.clf()
figStat.clf()
if len(stations2) > 0:
figStat2= plt.figure(figsize=(widthInInch, heightInInch))
......@@ -378,7 +404,7 @@ if __name__ == '__main__':
figStat2.subplots_adjust(wspace=0.3, hspace=0.3)
figStat2.savefig(figFileName, dpi=600, format='png')
plt.close()
fig2.clf()
figStat2.clf()
#finding out, who did worst/best -----------------------------------
worstRMSE, worstMAE, worstMBE = findWorstStation(stations)
......@@ -413,7 +439,7 @@ if __name__ == '__main__':
figScatterW.subplots_adjust(wspace=0.3, hspace=0.3)
figScatterW.savefig(figFileName, dpi=600, format='png')
plt.close()
fig3.clf()
figScatterW.clf()
# plotting best ----------------------------------------------------
figFileName = pathFigs + year + '_' + dataType + '_best.png'
heightInInch2 = heightInInch
......@@ -428,9 +454,9 @@ if __name__ == '__main__':
pltSctB5 = figScatterB.add_subplot(235)
pltSctB6 = figScatterB.add_subplot(236)
else:
pltSctB1 = fig4.add_subplot(131)
pltSctB2 = fig4.add_subplot(132)
pltSctB3 = fig4.add_subplot(133)
pltSctB1 = figScatterB.add_subplot(131)
pltSctB2 = figScatterB.add_subplot(132)
pltSctB3 = figScatterB.add_subplot(133)
stations[bestRMSE].createPlots(pltSctB1)
stations[bestMAE].createPlots(pltSctB2)
stations[bestMBE].createPlots(pltSctB3)
......@@ -443,7 +469,7 @@ if __name__ == '__main__':
figScatterB.subplots_adjust(wspace=0.3, hspace=0.3)
figScatterB.savefig(figFileName, dpi=600, format='png')
plt.close()
fig4.clf()
figScatterB.clf()
#- year's statistics ---------------------------------------------------
if createMap:
......
Code/station.py
View file @ 45ad43de
......@@ -95,10 +95,10 @@ class WeatherStation:
# starting to caluclate things
def getDifferences(self):
""""(void) -> ([float])
function calculates the differences between observed and model data.
afterwards returns the resulting array. this is needed for the statistical
analysis."""
"""(void) -> ([float])
function calculates the differences between observed and model data.
afterwards returns the resulting array. this is needed for the statistical
analysis."""
differences = []
totalLength = len(self.observationData)
if (len(self.modelData) < totalLength):
......@@ -137,9 +137,9 @@ class WeatherStation:
# translating heights, as it is necessary with the wind data
def getHeightIndices(height):
""""(float) -> (str, str)
function returns indices of COSMO-DE's main planes, that lie
around the specified height"""
"""(float) -> (str, str)
function returns indices of COSMO-DE's main planes, that lie
around the specified height"""
lower = ''
upper = ''
if height <= 10.:
......@@ -177,8 +177,8 @@ class WeatherStation:
return lower, upper
def getHeightFromIndex(index):
""""(str) -> (float)
function returns height of COSMO-DE's main plane with index index"""
"""(str) -> (float)
function returns height of COSMO-DE's main plane with index index"""
if index == '50':
return 10.
elif index == '49':
......@@ -248,14 +248,12 @@ class WeatherStation:
if (int(targetYear) % 4 == 0):
supposedLength = 8784
if len(self.observationData) < supposedLength:
print(self.m_name, 'only has', len(self.observationData))
print(self.name, 'only has', len(self.observationData))
self.observationData = []
def readInModelFile(self, h5FileName, dataType, h5FileName2 = '', \
height = 10.):
def readInModelFile(self, h5FileName, dataType, targetYear, \
h5FileName2 = '', height = 10.):
self.modelData = []
if self.name == 'Brocken':
print('Brocken, alter', h5FileName, h5FileName2, height)
if h5FileName2 == '':
with h5py.File(h5FileName, 'r') as h5File:
self.modelData =h5File[dataType][:,self.latitudeId,self.longitudeId]
......@@ -272,14 +270,20 @@ class WeatherStation:
for lower, upper in zip(data1, data2):
self.modelData.append(np.interp(height, [height1, height2], \
[lower, upper]))
supposedLength = 8760
if (int(targetYear) % 4 == 0):
supposedLength = 8784
if len(self.modelData) < supposedLength:
print(self.name, 'only has', len(self.modelData))
self.modelData = []
# calculating all the statistics needed ------------------------------------
def calculateStats(self):
""""(void) -> (void)
kicks of calculation process. call-stack is:
-getDifferences()
-calculateMAE() | calculateMBE() | calculateRMSE()
afterwards operations on solar data and some min/max analysis"""
"""(void) -> (void)
kicks of calculation process. call-stack is:
-getDifferences()
-calculateMAE() | calculateMBE() | calculateRMSE()
afterwards operations on solar data and some min/max analysis"""
differences = self.getDifferences();
self.mae = WeatherStation.calculateMAE(differences);
self.mbe = WeatherStation.calculateMBE(differences);
......@@ -299,28 +303,28 @@ class WeatherStation:
self.obsAvg = valueArray.mean()
def calculateMAE(differences):
""""([float]) -> (np.array)
function returns mean of absolutes of float array"""
"""([float]) -> (np.array)
function returns mean of absolutes of float array"""
return abs(np.array(differences)).mean()
def calculateMBE(differences):
""""([float]) -> (np.array)
function returns mean of float array"""
"""([float]) -> (np.array)
function returns mean of float array"""
return np.array(differences).mean()
def calculateRMSE(differences):
""""([float]) -> (np.array)
function returns mean of square-roots of squares of float array"""
"""([float]) -> (np.array)
function returns mean of square-roots of squares of float array"""
return math.sqrt(np.square(np.array(differences)).mean())
# making plots -------------------------------------------------------------
def createPlots(self, fig, printAxis = False):
""""(matplotlib.axes, bool) -> (void)
creates a scatter plot of obersvation data and model data (x and y,
repsectively). the scatter plot is added to the matplotlib axes given
at function call.
@remark: this re-uses some code of the differences() function. may be
this can be done less verbose and more well structed"""
"""(matplotlib.axes, bool) -> (void)
creates a scatter plot of obersvation data and model data (x and y,
repsectively). the scatter plot is added to the matplotlib axes given
at function call.
@remark: this re-uses some code of the differences() function. may be
this can be done less verbose and more well structed"""
differences = self.getDifferences();
valuesObs = []
valuesObsBounds = []
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment