diff --git a/api/GetAverageData.py b/api/GetAverageData.py
index bc5491c4390f7fc6d8ff8d3645b664c8563dcfa5..fb6a32682c0735b7da150c2ea4639e1a37de2320 100644
--- a/api/GetAverageData.py
+++ b/api/GetAverageData.py
@@ -1,6 +1,5 @@
 import configparser
-
-import psycopg2
+import numpy as np
 from psycopg2 import sql
 import numpy as np
 
@@ -12,6 +11,25 @@ param_postgres = cfg["POSTGRES"]
 param_interpol = cfg["INTERPOLATION"]
 
 
+def get_average_of_multiple_years(cursor, years):
+    avg_strings = " "
+    where_sql = " WHERE lat IS NOT NULL "
+    and_strings = ""
+    n = int(years[1]) - int(years[0])
+    for year in range(int(years[0]), int(years[1])+1):
+        avg_string = ' AVG ("{}") + '.format(str(year))
+        and_string = """ AND "{}" != 'NaN' """.format(str(year))
+
+        avg_strings += avg_string
+        and_strings += and_string
+    avg_strings = avg_strings[:-2]
+
+    query = """SELECT station_id, ROUND(({}) / {}, 1), transparent FROM stations WHERE file IS NULL GROUP BY station_id, transparent ORDER BY station_id ASC;""".format(avg_strings, n)
+    print(query)
+    cursor.execute(query)
+    return cursor.fetchall()
+
+
 # Getting all available year columns from database
 def get_year_columns(cursor):
     columns = []
diff --git a/api/api.py b/api/api.py
index 6ee89fce13d41d974a5c8f287726fd2042779908..a415431cb01865b37cf6e7a3a78c859f82611927 100644
--- a/api/api.py
+++ b/api/api.py
@@ -1,10 +1,11 @@
 import psycopg2
-from flask import Flask, jsonify, request
-from flask_cors import cross_origin
+from flask import Flask, jsonify, request, send_from_directory
+# from flask_cors import cross_origin
 from psycopg2 import sql
 import configparser
-from GetAverageData import get_interpolation_data_for_point
-import SQLPandasTools as s2pTool
+from GetAverageData import get_interpolation_data_for_point, get_year_columns, get_average_of_multiple_years
+from dataacquisition import write_raster
+# import SQLPandasTools as s2pTool
 
 cfg = configparser.ConfigParser()
 cfg.read('../config.ini')
@@ -94,22 +95,38 @@ def getStandardQuery():
     
     return query
 
+
+@app.route('/raster', methods=['GET'])
+def get_raster():
+    if 'years' in request.args:
+        years = request.args['years'].split('-')
+        if int(years[1]) < int(years[0]):
+            years = [years[1], years[0]]
+        with psycopg2.connect(database=param_postgres["dbName"], user=param_postgres["user"], password=param_postgres["password"], host=param_postgres["host"], port=param_postgres["port"]) as connection:
+            with connection.cursor() as cursor:
+                # all_years = get_year_columns(cursor)
+                average_data = get_average_of_multiple_years(cursor, years)
+                write_raster.write_raster(average_data)
+        return send_from_directory('D:/Uni/Master/01_SS2021/Automatisierte_Geodatenprozessierung/temperaturverteilung/dataacquisition/output', filename='myraster.tif', as_attachment=True)
+        # return 'Läuft, Brudi'
+
 @app.route('/annualMean', methods=['GET'])
 @cross_origin()
 def annualMean():
-    
+
     query = getStandardQuery()
-    
+
     with psycopg2.connect(database=param_postgres["dbName"], user=param_postgres["user"], password=param_postgres["password"], host=param_postgres["host"], port=param_postgres["port"]) as connection:
         with connection.cursor() as cursor:
-            
+
             cursor.execute(query)
-            
+
             results = cursor.fetchall()[0][0]
-            
+
             return jsonify(s2pTool.determineAnnualMean(results))
-    
+
     return "{}"
 
+
 if __name__ == '__main__':
     app.run(host='127.0.0.1', port=42000)
diff --git a/dataacquisition/sandbox_raster.py b/dataacquisition/sandbox_raster.py
index 53fc68ed6f12f10b2d292172cf8ec052f6ee8d80..2bfc98883ff445ba3b1b35d853e0146ed55ac536 100644
--- a/dataacquisition/sandbox_raster.py
+++ b/dataacquisition/sandbox_raster.py
@@ -16,28 +16,39 @@ with psycopg2.connect(database=db_name, user=user, password=pw, host=param_postg
     with connection.cursor() as cursor:
         cursor.execute('SELECT station_id, "2018" FROM stations WHERE file IS NULL ORDER BY station_id ASC;')
         pixel_array = []
+        pixel_array_alpha = []
         results = cursor.fetchall()
         for j in range(0, 36):
             row_array = []
+            row_array_alpha = []
             for i, station_id in enumerate(results):
                 if i % 36 == 0:
-                    row_array.append(results[i + j][1])                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                 
+                    value = results[i + j][1]
+                    row_array.append(value)
+                    row_array_alpha.append(value)
+                    np_row_array = np.array(row_array)
+                    np_row_array_alpha = np.array(row_array_alpha)
+                    np_row_array_alpha = np.where(np_row_array_alpha < 8, 0, 255)
                     # print(row_array)
-            pixel_array.append(row_array)
+            pixel_array.append(np_row_array)
+            pixel_array_alpha.append(np_row_array_alpha)
         np_pixel_array = np.array(pixel_array)
-        print(np_pixel_array)
-        
+        np_pixel_array_alpha = np.array(pixel_array_alpha)
+
         xmin, ymin, xmax, ymax = [5.01, 47.15, 14.81, 55.33]
         nrows, ncols = np.shape(np_pixel_array)
         xres = (xmax - xmin) / float(ncols)
         yres = (ymax - ymin) / float(nrows)
         geotransform = (xmin, xres, 0, ymax, 0, -yres)
 
-        output_raster = gdal.GetDriverByName('GTiff').Create('./output/myraster.tif', ncols, nrows, 1, gdal.GDT_Float32)  # Open the file
+        output_raster = gdal.GetDriverByName('GTiff').Create('D:/Uni/Master/01_SS2021/Automatisierte_Geodatenprozessierung/temperaturverteilung/dataacquisition/output/myraster.tif', ncols, nrows, 4, gdal.GDT_Float32)  # Open the file
         output_raster.SetGeoTransform(geotransform)
         srs = osr.SpatialReference()
         srs.ImportFromEPSG(4326)
         output_raster.SetProjection(srs.ExportToWkt())
         output_raster.GetRasterBand(1).WriteArray(np_pixel_array)
+        # output_raster.GetRasterBand(2).WriteArray(np_pixel_array)
+        # output_raster.GetRasterBand(3).WriteArray(np_pixel_array)
+        output_raster.GetRasterBand(4).WriteArray(np_pixel_array_alpha)
 
-        output_raster.FlushCache()
\ No newline at end of file
+        output_raster.FlushCache()
diff --git a/dataacquisition/sandbox_raster.py.orig b/dataacquisition/sandbox_raster.py.orig
new file mode 100644
index 0000000000000000000000000000000000000000..af76b563313e98bfcdc8941d86a255e7df044032
--- /dev/null
+++ b/dataacquisition/sandbox_raster.py.orig
@@ -0,0 +1,58 @@
+import configparser
+import numpy as np
+import psycopg2
+import pandas as pd
+from osgeo import gdal, osr
+
+cfg = configparser.ConfigParser()
+cfg.read('../config.ini')
+assert "POSTGRES" in cfg, "missing POSTGRES in config.ini"
+param_postgres = cfg["POSTGRES"]
+param_interpol = cfg["INTERPOLATION"]
+
+db_name, user, pw, host, port = param_postgres['dbName'], param_postgres["user"], param_postgres["password"], param_postgres["host"], param_postgres["port"]
+
+with psycopg2.connect(database=db_name, user=user, password=pw, host=param_postgres["host"], port=port) as connection:
+    with connection.cursor() as cursor:
+        cursor.execute('SELECT station_id, "2018" FROM stations WHERE file IS NULL ORDER BY station_id ASC;')
+        pixel_array = []
+        pixel_array_alpha = []
+        results = cursor.fetchall()
+        for j in range(0, 36):
+            row_array = []
+            row_array_alpha = []
+            for i, station_id in enumerate(results):
+                if i % 36 == 0:
+<<<<<<< HEAD
+                    row_array.append(results[i + j][1])                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                 
+=======
+                    value = results[i + j][1]
+                    row_array.append(value)
+                    row_array_alpha.append(value)
+                    np_row_array = np.array(row_array)
+                    np_row_array_alpha = np.array(row_array_alpha)
+                    np_row_array_alpha = np.where(np_row_array_alpha < 8, 0, 255)
+>>>>>>> 94001e7882d184c3101f8663bc43b41ad9cb4e98
+                    # print(row_array)
+            pixel_array.append(np_row_array)
+            pixel_array_alpha.append(np_row_array_alpha)
+        np_pixel_array = np.array(pixel_array)
+        np_pixel_array_alpha = np.array(pixel_array_alpha)
+
+        xmin, ymin, xmax, ymax = [5.01, 47.15, 14.81, 55.33]
+        nrows, ncols = np.shape(np_pixel_array)
+        xres = (xmax - xmin) / float(ncols)
+        yres = (ymax - ymin) / float(nrows)
+        geotransform = (xmin, xres, 0, ymax, 0, -yres)
+
+        output_raster = gdal.GetDriverByName('GTiff').Create('D:/Uni/Master/01_SS2021/Automatisierte_Geodatenprozessierung/temperaturverteilung/dataacquisition/output/myraster.tif', ncols, nrows, 4, gdal.GDT_Float32)  # Open the file
+        output_raster.SetGeoTransform(geotransform)
+        srs = osr.SpatialReference()
+        srs.ImportFromEPSG(4326)
+        output_raster.SetProjection(srs.ExportToWkt())
+        output_raster.GetRasterBand(1).WriteArray(np_pixel_array)
+        # output_raster.GetRasterBand(2).WriteArray(np_pixel_array)
+        # output_raster.GetRasterBand(3).WriteArray(np_pixel_array)
+        output_raster.GetRasterBand(4).WriteArray(np_pixel_array_alpha)
+
+        output_raster.FlushCache()
diff --git a/dataacquisition/write_raster.py b/dataacquisition/write_raster.py
new file mode 100644
index 0000000000000000000000000000000000000000..e3c4d79b466583b982b9b21fa3f361ac5411e00b
--- /dev/null
+++ b/dataacquisition/write_raster.py
@@ -0,0 +1,46 @@
+import numpy as np
+from osgeo import gdal, osr
+
+
+def write_raster(data):
+    pixel_array = []
+    pixel_array_alpha = []
+    for j in range(0, 36):
+        row_array = []
+        row_array_alpha = []
+        for i, station_id in enumerate(data):
+            if i % 36 == 0:
+                value = data[i + j][1]
+                value = 0 if not value else value
+                value = 0 if str(value) == 'NaN' else value
+                alpha_value = data[i + j][2]
+                alpha_value = 0 if alpha_value else 255
+                # print(value)
+                row_array.append(value)
+                row_array_alpha.append(alpha_value)
+                np_row_array = np.array(row_array)
+                np_row_array_alpha = np.array(row_array_alpha)
+                # np_row_array_alpha = np.where(np_row_array_alpha < 8, 0, 255)
+                # print(row_array)
+        pixel_array.append(np_row_array)
+        pixel_array_alpha.append(np_row_array_alpha)
+    np_pixel_array = np.array(pixel_array)
+    np_pixel_array_alpha = np.array(pixel_array_alpha)
+
+    xmin, ymin, xmax, ymax = [5.01, 47.15, 14.81, 55.33]
+    nrows, ncols = np.shape(np_pixel_array)
+    xres = (xmax - xmin) / float(ncols)
+    yres = (ymax - ymin) / float(nrows)
+    geotransform = (xmin, xres, 0, ymax, 0, -yres)
+
+    output_raster = gdal.GetDriverByName('GTiff').Create('D:/Uni/Master/01_SS2021/Automatisierte_Geodatenprozessierung/temperaturverteilung/dataacquisition/output/myraster.tif', ncols, nrows, 4, gdal.GDT_Float32)  # Open the file
+    output_raster.SetGeoTransform(geotransform)
+    srs = osr.SpatialReference()
+    srs.ImportFromEPSG(4326)
+    output_raster.SetProjection(srs.ExportToWkt())
+    output_raster.GetRasterBand(1).WriteArray(np_pixel_array)
+    # output_raster.GetRasterBand(2).WriteArray(np_pixel_array)
+    # output_raster.GetRasterBand(3).WriteArray(np_pixel_array)
+    output_raster.GetRasterBand(4).WriteArray(np_pixel_array_alpha)
+
+    output_raster.FlushCache()