download_service.py 10.4 KB
Newer Older
Lozac'h Loic's avatar
Lozac'h Loic committed
1
#!/usr/bin/python3
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
2
import os, sys, argparse, shlex, glob, time
Lozac'h Loic's avatar
Lozac'h Loic committed
3
import cdsapi
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
4
from subprocess import Popen, PIPE, STDOUT
Lozac'h Loic's avatar
Lozac'h Loic committed
5
6
from osgeo import gdal, osr, ogr

Lozac'h Loic's avatar
Lozac'h Loic committed
7

Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
8
9
10
11
12
def stream_process(process):
    go = process.poll() is None
    for line in process.stdout:
        print(line.strip())
    return go
Lozac'h Loic's avatar
Lozac'h Loic committed
13
14
15

def process_command(cmd):
    print("Starting : "+" ".join(cmd))
Lozac'h Loic's avatar
Lozac'h Loic committed
16
17
18
    process = Popen(" ".join(cmd), stdout=PIPE, stderr=STDOUT, shell=True)
    while stream_process(process):
        time.sleep(0.1)
Lozac'h Loic's avatar
Lozac'h Loic committed
19
20
21

def download_cds(m_year, m_month, m_day, m_time, m_latmax, m_lonmin, m_latmin, m_lonmax, outfile):
    c = cdsapi.Client()
Lozac'h Loic's avatar
Lozac'h Loic committed
22
    requ = (
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
23
    'reanalysis-era5-land',
Lozac'h Loic's avatar
Lozac'h Loic committed
24
    {
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
25
        'format': 'grib',
Lozac'h Loic's avatar
Lozac'h Loic committed
26
27
28
29
30
31
32
33
34
35
36
        'variable': '2m_temperature',
        'year': m_year,#AAAA
        'month': m_month,#MM
        'day': m_day,#DD
        'time': m_time,#HH:00
        'area': [
            m_latmax, m_lonmin, m_latmin,
            m_lonmax,
        ],
    },
    outfile)
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
37
38
39
40
41
    try:
        c.retrieve(requ[0],requ[1],requ[2])
    except Exception as err:
        print(err)
        return requ
Lozac'h Loic's avatar
Lozac'h Loic committed
42
43
    
    return requ
Lozac'h Loic's avatar
Lozac'h Loic committed
44
    
Lozac'h Loic's avatar
Lozac'h Loic committed
45
46
47
48
49
50
51
52
def download_cds_inst(m_year, m_month, m_day, m_time, m_latmax, m_lonmin, m_latmin, m_lonmax, outfile):
    c = cdsapi.Client()
    requ = (
    'reanalysis-era5-pressure-levels',
    {
        'product_type': 'reanalysis',
        'format': 'grib',
        'variable': 'temperature',
Lozac'h Loic's avatar
Lozac'h Loic committed
53
        'pressure_level': '1000',
Lozac'h Loic's avatar
Lozac'h Loic committed
54
55
56
57
58
59
60
61
62
63
        'year': m_year,#AAAA
        'month': m_month,#MM
        'day': m_day,#DD
        'time': m_time,#HH:00
        'area': [
            m_latmax, m_lonmin, m_latmin,
            m_lonmax,
        ],
    },
    outfile)
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
64
65
66
67
    try:
        c.retrieve(requ[0],requ[1],requ[2])
    except Exception as err:
        print(err)
Lozac'h Loic's avatar
Lozac'h Loic committed
68
        requ[0]=err
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
69
        return requ
Lozac'h Loic's avatar
Lozac'h Loic committed
70
71
72
    
    return requ

Lozac'h Loic's avatar
Lozac'h Loic committed
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
def search_files(directory='.', resolution='S1', extension='SAFE', fictype='d'):
    images=[]
    extension = extension.lower()
    resolution = resolution.lower()
    for dirpath, dirnames, files in os.walk(directory):
        if fictype == 'f':
            for name in files:
    #            print(os.path.join(dirpath, name) + " test")
                if extension and name.lower().endswith(extension) and name.lower().find(resolution) >= 0 :
                    
    #                print(os.path.join(dirpath, name) + " OK")
                    abspath = os.path.abspath(os.path.join(dirpath, name))
                    images.append(abspath)
        elif fictype == 'd':
            for dirname in dirnames:
    #            print(os.path.join(dirpath, name) + " test")
                if extension and dirname.lower().endswith(extension) and dirname.lower().find(resolution) >= 0 :
                    
    #                print(os.path.join(dirpath, name) + " OK")
                    abspath = os.path.abspath(os.path.join(dirpath, dirname))
                    images.append(abspath)
        else:
            print("search_files type error")
            exit()
            
    return images

def round_s1time(s1time):
    m_hour = int(s1time[:2])
    m_minu = int(s1time[2:4])
    m_seco = int(s1time[4:6])
    
    if m_seco >= 30 :
        m_minu += 1;
    if m_minu >= 30 :
        m_hour += 1;
    if m_hour == 24:
        m_hour -= 1;
    st_hour = str(m_hour)
    if m_hour < 10 :
        st_hour = "0"+st_hour
    return st_hour+":00"

Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
116
def round_extend(points, supinf):
Lozac'h Loic's avatar
Lozac'h Loic committed
117
    
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
    for pt in points:
        if supinf == "sup":
            pt.SetPoint(0,int(pt.GetX()*10)/10+0.1,int(pt.GetY()*10)/10+0.1)
        elif supinf == "inf":
            pt.SetPoint(0,int(pt.GetX()*10)/10-0.1,int(pt.GetY()*10)/10-0.1)
        
#     val=val*10
#     pe = int(val)
#     deci = val/pe
#     if supinf == "sup":
#         if deci < 0.25 :
#             deci = 0.25
#         elif 0.25 <= deci < 0.5:
#             deci = 0.5
#         elif 0.5 <= deci < 0.75:
#             deci = 0.75
#         else:
#             deci = 0.
#             pe += 1
#     elif supinf == "inf":
#         if deci < 0.25 :
#             deci = 0.
#         elif 0.25 <= deci < 0.5:
#             deci = 0.25
#         elif 0.5 <= deci < 0.75:
#             deci = 0.5
#         else:
#             deci = 0.75
#     
#     return (pe+deci)/10
Lozac'h Loic's avatar
Lozac'h Loic committed
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
    
def get_ref_extend(imgref):
    raster = gdal.Open(imgref)
    proj = osr.SpatialReference(wkt=raster.GetProjection())
    #print("ref: "+raster.GetProjection())
    upx, xres, xskew, upy, yskew, yres = raster.GetGeoTransform()
    cols = raster.RasterXSize
    rows = raster.RasterYSize
     
    ulx = upx + 0*xres + 0*xskew
    uly = upy + 0*yskew + 0*yres
     
    llx = upx + 0*xres + rows*xskew
    lly = upy + 0*yskew + rows*yres
     
    lrx = upx + cols*xres + rows*xskew
    lry = upy + cols*yskew + rows*yres
     
    urx = upx + cols*xres + 0*xskew
    ury = upy + cols*yskew + 0*yres
    
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
169
    pointUL, pointLR = ogr.Geometry(ogr.wkbPoint), ogr.Geometry(ogr.wkbPoint)
Lozac'h Loic's avatar
Lozac'h Loic committed
170
    pointLL, pointUR = ogr.Geometry(ogr.wkbPoint), ogr.Geometry(ogr.wkbPoint)
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
171
172
    pointUL.AddPoint(ulx, uly)
    pointLR.AddPoint(lrx, lry)
Lozac'h Loic's avatar
Lozac'h Loic committed
173
174
175
176
177
178
179
180
    pointLL.AddPoint(llx, lly)
    pointUR.AddPoint(urx, ury)
    
    spatialRef = osr.SpatialReference()
    spatialRef.ImportFromEPSG(4326) 

    if not proj.IsSame(spatialRef) :
        coordTransform = osr.CoordinateTransformation(proj,spatialRef)
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
181
182
        pointUL.Transform(coordTransform)
        pointLR.Transform(coordTransform)
Lozac'h Loic's avatar
Lozac'h Loic committed
183
184
        pointLL.Transform(coordTransform)
        pointUR.Transform(coordTransform)
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
185
186
187
188
    
    print(str(pointUR))
    print(str(pointUL))
    print(str(pointLL))
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
189
    print(str(pointLR))
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
190
    
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
191
192
193
194
    latmax = str(max([pointUL.GetY(),pointLL.GetY(),pointUR.GetY(),pointLR.GetY()])+0.2)
    latmin = str(min([pointUL.GetY(),pointLL.GetY(),pointUR.GetY(),pointLR.GetY()])-0.2)
    lonmax = str(max([pointUL.GetX(),pointLL.GetX(),pointUR.GetX(),pointLR.GetX()])+0.2)
    lonmin = str(min([pointUL.GetX(),pointLL.GetX(),pointUR.GetX(),pointLR.GetX()])-0.2)
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
195
    
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
196
197
#     pointLL.SetPoint(0,round_extend(pointLL.GetX(), "inf"),round_extend(pointLL.GetY(), "inf"))
#     pointUR.SetPoint(0,round_extend(pointUR.GetX(), "sup"),round_extend(pointUR.GetY(), "sup"))
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
198
    print([latmax, lonmin, latmin, lonmax])
Lozac'h Loic's avatar
Lozac'h Loic committed
199
    
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
200
    return [latmax, lonmin, latmin, lonmax]
Lozac'h Loic's avatar
Lozac'h Loic committed
201
202
203
204
205
206
207
208
209
210
211
212

if __name__ == "__main__":
    
    # Make parser object
    parser = argparse.ArgumentParser(description=
        """
        Download services
        """)

    subparsers = parser.add_subparsers(help='Choose server', dest="pipeline")

    # Short pipeline
Lozac'h Loic's avatar
Lozac'h Loic committed
213
214
215
216
217
218
219
220
221
222
    list_parser = subparsers.add_parser('eodag', help="Download S1 through Eodag API")
    list_parser.add_argument('-downdir', action='store', required=True, help="Download directory")
    list_parser.add_argument('-tile_str', action='store', required=True, help="S2 Tile name, ex. T30UVU")
    list_parser.add_argument('-startdate', action='store', required=True, help="Starting date S1, ex. 2020-08-31")
    list_parser.add_argument('-enddate', action='store', required=True, help="Ending date S1, ex. 2020-09-31")
    
    list_parser = subparsers.add_parser('eodagrestart', help="Restart download S1 through Eodag API using Eodag SearchResult geojson files")
    list_parser.add_argument('-downdir', action='store', required=True, help="Download directory")
    list_parser.add_argument('-remaining_geojson', action='store', required=True, help="Remaining geojson file created from eodag pipeline")
    list_parser.add_argument('-alternate_geojson', action='store', required=False, help="[Optional] Alternate geojson file created from eodag pipeline")
Lozac'h Loic's avatar
Lozac'h Loic committed
223
224
225
226
227
    
    list_parser = subparsers.add_parser('cds', help="Download ERA5 temperature image file according to S1 dates")
    list_parser.add_argument('-s1dir', action='store', required=True, help="Calibrated Sentinel1 directory")
    list_parser.add_argument('-ref', action='store', required=True, help="AOI image reference")
    list_parser.add_argument('-zone', action='store', required=True, help='Geographic zone reference for output files naming')
Lozac'h Loic's avatar
Lozac'h Loic committed
228
    list_parser.add_argument('-gap', choices=['3m', '6d'], required=False, default='land', help="[Optional] Query ERA5-Land reanalysis collection with 3 month gap (default '3m') OR ERA5 basic reanalysis collection with 6 days gap ('6d') ")
Lozac'h Loic's avatar
Lozac'h Loic committed
229
230
231
232
    list_parser.add_argument('-outdir', action='store', required=True, help='Output directory')
        
    args=parser.parse_args()
    
Lozac'h Loic's avatar
Lozac'h Loic committed
233
234
235
236
237
238
239
    if args.pipeline == 'eodag' :
        cmd=["python", "manage.py", "runscript", "freezeDetectWrapper", "--traceback", "--script-args",
              "12", "downdir@"+args.downdir, "tile_str@"+args.tile_str, "startdate@"+args.startdate, "enddate@"+args.enddate
            ]
        process_command(cmd)
    elif args.pipeline == 'eodagrestart' :
        cmd=["python", "manage.py", "runscript", "freezeDetectWrapper", "--traceback", "--script-args",
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
240
              "13", "remaining_geojson@"+args.remaining_geojson, "downdir@"+args.downdir, "alternate_geojson@"+str(args.alternate_geojson)
Lozac'h Loic's avatar
Lozac'h Loic committed
241
242
243
            ]
        process_command(cmd)
        
Lozac'h Loic's avatar
Lozac'h Loic committed
244
245
246
247
248
    elif args.pipeline == 'cds' :
        outabs = os.path.abspath(args.outdir)
        if not os.path.exists(outabs):
            os.makedirs(outabs)
        
Lozac'h Loic's avatar
Lozac'h Loic committed
249
#         s1files = search_files(args.s1dir, 'S1', 'tif', fictype='f')
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
250
        s1files = glob.glob(os.path.join(args.s1dir, 'S1*.TIF'))
Lozac'h Loic's avatar
debug    
Lozac'h Loic committed
251
        extend = get_ref_extend(args.ref)
Lozac'h Loic's avatar
Lozac'h Loic committed
252
        errfile = os.path.join(args.outdir, "download_error.log")
Lozac'h Loic's avatar
Lozac'h Loic committed
253
        with open(errfile, 'w') as err:
Lozac'h Loic's avatar
Lozac'h Loic committed
254
255
256
257
258
            for s1 in s1files:
                stdatetime = os.path.basename(s1)[:-4].split("_")[4]
                stdate = stdatetime.split("T")[0]
                sttime = stdatetime.split("T")[1]
                rdtime = round_s1time(sttime) 
Lozac'h Loic's avatar
Lozac'h Loic committed
259
260
261
                
                if args.gap == '3m':
                    outfile = os.path.join(args.outdir,"ERA5-2m_"+args.zone+"_"+stdatetime+".grib")
Lozac'h Loic's avatar
Lozac'h Loic committed
262
                    print("Downloading CDS : "+outfile)
Lozac'h Loic's avatar
Lozac'h Loic committed
263
264
                    request = download_cds(stdate[:4], stdate[4:6], stdate[6:8], rdtime, extend[0], extend[1], extend[2], extend[3], outfile)
                else:
Lozac'h Loic's avatar
Lozac'h Loic committed
265
                    outfile = os.path.join(args.outdir,"ERA5-1000hpa_"+args.zone+"_"+stdatetime+".grib")
Lozac'h Loic's avatar
Lozac'h Loic committed
266
                    print("Downloading CDS : "+outfile)
Lozac'h Loic's avatar
Lozac'h Loic committed
267
                    request = download_cds_inst(stdate[:4], stdate[4:6], stdate[6:8], rdtime, extend[0], extend[1], extend[2], extend[3], outfile)
Lozac'h Loic's avatar
Lozac'h Loic committed
268
269
270
                
                if not os.path.exists(outfile):
                    err.write(request[2]+";"+request[1]+"\n")
Lozac'h Loic's avatar
Lozac'h Loic committed
271
                    err.write(request[0]+"\n")
Lozac'h Loic's avatar
Lozac'h Loic committed
272
273
274
275
276
277
278
279