Project

General

Profile

Download (9.76 KB) Statistics
| Branch: | Revision:
1
# -*- coding: utf-8 -*-
2
"""
3
Created on Mon May 13 10:02:02 2013
4

    
5
@author: parmentier
6
"""
7

    
8
#!/usr/bin/python
9
#
10
# Early version of assorted Python code to download MODIS 11A1 (daily)
11
# data associated with specified dates and tiles, then interface with
12
# GRASS to calculate temporally averaged LST in a way that accounts for
13
# QC flags.
14
#
15
# TODO:
16
#  - functionalize to encapsulate high level procedural steps
17
#  - create a 'main' function to allow script to be run as a command
18
#    that can accept arguments?
19
#  - put downloaded data somewhere other than working directory?
20
#  - record all downloads to a log file?
21
#
22
# Created by Jim Regetz
23
# NCEAS on 16-May-2012
24
# Modified by Benoit Parmentier
25
#NCEAS on 18-October-2012
26
#NCEAS on 19-January-2013
27
#NCEAS on 14-May-2013
28
# 
29

    
30
import os, glob
31
import datetime, calendar
32
import ftplib
33
#import grass.script as gs
34
import argparse
35
import errno
36
#from datetime import date
37
#from datetime import datetime
38

    
39

    
40
#------------------
41
# helper functions
42
#------------------
43

    
44
def yj_to_ymd(year, doy):
45
    """Return date as e.g. '2000.03.05' based on specified year and
46
    numeric day-of-year (doy) """
47
    date = datetime.datetime.strptime('%d%03d' % (year, doy), '%Y%j')
48
    return date.strftime('%Y.%m.%d')
49

    
50
def get_doy_range(year, month):
51
    """Determine starting and ending numeric day-of-year (doy)
52
    asscociated with the specified month and year.
53

    
54
    Arguments:
55
    year -- four-digit integer year
56
    month -- integer month (1-12)
57

    
58
    Returns tuple of start and end doy for that month/year.
59
    """
60
    last_day_of_month = calendar.monthrange(year, month)[1]
61
    start_doy = int(datetime.datetime.strptime('%d.%02d.%02d' % (year,
62
        month, 1), '%Y.%m.%d').strftime('%j'))
63
    end_doy = int(datetime.datetime.strptime('%d.%02d.%02d' % (year,
64
        month, last_day_of_month), '%Y.%m.%d').strftime('%j'))
65
    return (int(start_doy), int(end_doy))
66

    
67
#added by Benoit, to create a list of raster images for particular month (2001-2010)
68

    
69
def date_sequence(year_start,month_start,day_start,year_end,month_end,day_end):
70
    """Create a date sequence for a start date and end date defined by year, month, date
71
    Arguments:
72
    year_start -- starting year
73
 
74
    Returns list of absolute paths to the downloaded HDF files.Note the format yyyydoy
75
    """
76
    from datetime import date
77
    from dateutil.rrule import rrule, DAILY
78

    
79
    a = date(year_start, month_start,day_start)
80
    b = date(year_end, month_end, day_end)
81
    date_seq=list()
82
    for dt in rrule(DAILY, dtstart=a, until=b):
83
        print dt.strftime("%Y%j")    
84
        date_seq.append(dt.strftime("%Y%j"))
85
        
86
    return(date_seq)
87

    
88
# quick function to return list of dirs in wd
89
def list_contents(ftp):
90
    """Parse ftp directory listing into list of names of the files
91
    and/or directories contained therein. May or may not be robust in
92
    general, but seems to work fine for LP DAAP ftp server."""
93
    listing = []
94
    ftp.dir(listing.append)
95
    contents = [item.split()[-1] for item in listing[1:]]
96
    return contents
97

    
98
def download_mod11a1(destdir, tile, start_doy, end_doy, year, ver=5):
99
    """Download into destination directory the MODIS 11A1 HDF files
100
    matching a given tile, year, and day range. If for some unexpected
101
    reason there are multiple matches for a given day, only the first is
102
    used. If no matches, the day is skipped with a warning message.
103

    
104
    Arguments:
105
    destdir -- path to local destination directory
106
    tile -- tile identifier (e.g., 'h08v04')
107
    start_doy -- integer start of day range (0-366)
108
    end_doy -- integer end of day range (0-366)
109
    year -- integer year (>=2000)
110
    ver -- MODIS version (4 or 5)
111

    
112
    Returns list of absolute paths to the downloaded HDF files.
113
    """
114
    # connect to data pool and change to the right directory
115
    ftp = ftplib.FTP('e4ftl01.cr.usgs.gov')
116
    ftp.login('anonymous', '')
117
    ftp.cwd('MOLT/MOD11A1.%03d' % ver)
118
    # make list of daily directories to traverse
119
    available_days = list_contents(ftp)
120
    desired_days = [yj_to_ymd(year, x) for x in range(start_doy, end_doy+1)]
121
    days_to_get = filter(lambda day: day in desired_days, available_days)
122
    if len(days_to_get) < len(desired_days):
123
        missing_days = [day for day in desired_days if day not in days_to_get]
124
        print 'skipping %d day(s) not found on server:' % len(missing_days)
125
        print '\n'.join(missing_days)
126
    # get each tile in turn
127
    hdfs = list()
128
    for day in days_to_get:
129
        ftp.cwd(day)
130
        files_to_get = [file for file in list_contents(ftp)
131
            if tile in file and file[-3:]=='hdf']
132
        if len(files_to_get)==0:
133
            # no file found -- print message and move on
134
            print 'no hdf found on server for tile', tile, 'on', day
135
            ftp.cwd('..') #added by Benoit
136
            continue
137
        elif 1<len(files_to_get):
138
            # multiple files found! -- just use the first...
139
            print 'multiple hdfs found on server for tile', tile, 'on', day
140
        file = files_to_get[0]
141
        local_file = os.path.join(destdir, file)
142
        ftp.retrbinary('RETR %s' % file, open(local_file, 'wb').write)
143
        hdfs.append(os.path.abspath(local_file))
144
        ftp.cwd('..')
145
    # politely disconnect
146
    ftp.quit()
147
    # return list of downloaded paths
148
    return hdfs
149

    
150
def get_hdf_paths(hdfdir, tile, start_doy, end_doy, year):
151
    """Look in supplied directory to find the MODIS 11A1 HDF files
152
    matching a given tile, year, and day range. If for some unexpected
153
    reason there are multiple matches for a given day, only the first is
154
    used. If no matches, the day is skipped with a warning message.
155

    
156
    Arguments:
157
    hdfdir -- path to directory containing the desired HDFs
158
    tile -- tile identifier (e.g., 'h08v04')
159
    start_doy -- integer start of day range (0-366)
160
    end_doy -- integer end of day range (0-366)
161
    year -- integer year (>=2000)
162

    
163
    Returns list of absolute paths to the located HDF files.
164
    """
165
    hdfs = list()
166
    for doy in range(start_doy, end_doy+1):
167
        fileglob = 'MOD11A1.A%d%03d.%s*hdf' % (year, doy, tile)
168
        pathglob = os.path.join(hdfdir, fileglob)
169
        files = glob.glob(pathglob)
170
        if len(files)==0:
171
            # no file found -- print message and move on
172
            print 'cannot access %s: no such file' % pathglob
173
            continue
174
        elif 1<len(files):
175
            # multiple files found! -- just use the first...
176
            print 'multiple hdfs found for tile', tile, 'on', day
177
        hdfs.append(os.path.abspath(files[0]))
178
    return hdfs
179

    
180
def create_dir_and_check_existence(path):
181
    try:
182
        os.makedirs(path)
183
    except OSError as exception:
184
        if exception.errno != errno.EEXIST:
185
            raise
186
            
187
def main():
188
    #from sys import argv                       # example client code
189
    #myargs = getopts(argv)
190
    
191
    #import argparse
192
    parser = argparse.ArgumentParser()
193
    parser.add_argument("tiles", type=str, help="list of Modis tiles")
194
    parser.add_argument("start_year", type=int, help="start year")
195
    parser.add_argument("end_year", type=int, help="end year")
196
    parser.add_argument("start_month", type=int, help="start month")
197
    parser.add_argument("end_month", type=int, help="end month")
198
    parser.add_argument("hdfdir", type=str, help="destination/source directory for hdf file")
199
    parser.add_argument("night", type=int, help="night")
200
    parser.add_argument("download", type=int, help="out_suffix")
201
    parser.add_argument("out_suffix", type=str, help="out_suffix")
202

    
203
    myargs = parser.parse_args()
204

    
205
    tiles = myargs.tiles #These tiles correspond to Venezuela.
206
    start_year = myargs.start_year
207
    end_year = myargs.end_year 
208
    end_month = myargs.end_month
209
    start_month= myargs.start_month
210
    hdfdir =  myargs.hdfdir
211
    night= myargs.night    # if 1 then produce night climatology
212
    out_suffix= myargs.out_suffix #"_03192013"
213
    download=myargs.download# if 1 then download files
214
    
215
    ################## First step: download data ###################
216
    # Added tile loop 
217
    year_list=range(start_year,end_year+1) #list of year to loop through
218
    #print 'Number of arguments:', len(myargs), 'arguments.'
219
    #print 'Argument List:', str(myargs)
220
    
221
    tiles =tiles.split(",") #Create a list from string
222
    print 'this is the result',str(year_list)
223
    print(myargs)
224
    print(tiles)
225
 
226
 
227
    #tiles = ['h11v08','h11v07','h12v07','h12v08','h10v07','h10v08'] #These tiles correspond to Venezuela.
228
    #start_year = 2001
229
    #end_year = 2010
230
    #end_month=12
231
    #start_month=1
232
    #hdfdir =  '/home/parmentier/Data/benoit_test' #destination file where hdf files are stored locally after download.
233
    #night=1    # if 1 then produce night climatology
234
    #out_suffix="_03192013"
235
    #download=1  # if 1 then download files
236
    
237
    ################## First step: download data ###################
238
    # Added tile loop 
239
    year_list=range(start_year,end_year+1) #list of year to loop through
240
    #for testing...
241
    #year_list=[2001,2002]
242
    #hdfdir = '/home/parmentier/Data/benoit_test2' 
243
    #tiles = ['h10v07','h10v08','h12v07']
244
    if download==1:
245
        for tile in tiles:
246
            outDir = os.path.join(hdfdir,tile)
247
            create_dir_and_check_existence(outDir)
248
            
249
            for year in year_list:
250
                start_doy = 1
251
                end_doy = 365
252
                if calendar.isleap(year):
253
                    end_doy=366
254
                    
255
                hdfs = download_mod11a1(outDir, tile, start_doy, end_doy, year)
256

    
257
if __name__ == '__main__':
258
    main()
259
#    from sys import argv                       # example client code
260
#    myargs = getopts(argv)
261
#    if '-a' in myargs:
262
#        print(myargs['-a'])
263
#    if '-b' in myargs:
264
#        print(myargs['-b'])
265
#    print(myargs)
266
    
267
#i = 1
268
#while i < arg_length:
269
#    arg = argv[i]
270
#    if arg =='-a':
271
#        a     = argv[i]
272
#    elif arg =='-b':
273
#        b     = argv[i]
274

    
275
#x = a + b
276

    
277
#print x
(14-14/22)