Commits
Neal Schweighart authored and Ville Suoranta committed ab4ded97aab Merge
149 149 | ## |
150 150 | ## |
151 151 | ## ============================================================================= |
152 152 | ## ============================================================================= |
153 153 | |
154 154 | |
155 155 | import glob, pylab, os, datetime |
156 156 | import numpy as np |
157 157 | from matplotlib import rc |
158 158 | import matplotlib.pyplot as plt |
159 + | import ftplib |
159 160 | |
160 161 | from casatasks.private.casa_transition import is_CASA6 |
161 162 | if is_CASA6: |
162 163 | from casatools import table, quanta, coordsys, image, measures |
163 164 | |
164 165 | tb = table() |
165 166 | qa = quanta() |
166 167 | cs = coordsys() |
167 168 | ia = image() |
168 169 | me = measures() |
301 302 | |
302 303 | ## Set up the days for which we need to go get TEC files |
303 304 | day_list = [] |
304 305 | next_day = begin_day |
305 306 | for iter in range(call_num): |
306 307 | day_list.append(next_day) |
307 308 | next_day = qa.time(str(t_min[0]+86400.*(iter+1))+'s',form='ymd')[0][:10] |
308 309 | |
309 310 | print('IGS files required for: '+str(day_list)) |
310 311 | |
312 + | # plot name declared here to avoid potential reference before assignment error |
313 + | plot_name = '' |
311 314 | ## Runs the IGS methods |
312 315 | if tec_server == 'IGS': |
313 316 | ymd_date_num = 0 |
314 317 | #array = [] |
315 318 | lo=0 |
316 319 | hi=0 |
317 320 | for ymd_date in day_list: |
318 321 | points_long,points_lat,ref_long,ref_lat,incr_long,incr_lat,incr_time,num_maps,tec_array,tec_type = get_IGS_TEC(ymd_date) |
319 322 | |
320 323 | ## Fill a new array with all the full set of TEC/DTEC values for all days in the observation set. |
486 489 | ## released ~ 2-3 days after data is collected. If this file is |
487 490 | ## unavailable, it will try to retrieve the JPL Rapid Product (JPRG), |
488 491 | ## released ~ 1 day after data is collected. While the 'uncompress' command |
489 492 | ## is not necessary, it is the most straightforward on Linux. |
490 493 | ## |
491 494 | ## ========================================================================= |
492 495 | |
493 496 | #CDDIS = 'ftp://cddis.gsfc.nasa.gov/gnss/products/ionex/' # pre-2020Nov01 version |
494 497 | CDDIS = 'ftp://gdc.cddis.eosdis.nasa.gov/gps/products/ionex/' # new, more secure ftp-ssl server (2020Nov01) |
495 498 | file_location = CDDIS+str(year)+'/'+str(dayofyear)+'/' |
496 - | curlcmd='curl -u anonymous:casa-feedback@nrao.edu --ftp-ssl-reqd ' |
499 + | #curlcmd='curl -u anonymous:casa-feedback@nrao.edu --ftp-ssl-reqd ' |
500 + | |
497 501 | |
498 502 | ## The name of the IONEX file you require. |
499 503 | igs_file='igsg'+str(dayofyear)+'0.'+str(year)[2:4]+'i' if ( gpsweek<2238) else 'IGS0OPSFIN_'+str(year)+str(dayofyear)+'0000_01D_02H_GIM.INX' |
500 504 | get_file=igs_file + ('.Z' if gpsweek<2238 else '.gz') |
501 505 | |
502 506 | print('\nFor '+ymd_date+', the required IGS file is called: '+igs_file) |
503 507 | |
504 508 | if len(glob.glob(igs_file))<1: # file does not yet exist locally |
505 - | print('Attempting retrieval of IGS Final product file: '+str(get_file)) |
506 - | get_path = file_location+get_file |
507 - | if test_IONEX_connection(get_path): |
508 - | os.system(curlcmd+get_path+' > '+workDir+get_file) |
509 - | os.system('gunzip '+get_file) |
510 - | else: |
511 - | # IGS final product file does not exist; try rapid product file |
512 - | igs_file='igrg'+str(dayofyear)+'0.'+str(year)[2:4]+'i' if ( gpsweek<2238) else 'IGS0OPSRAP_'+str(year)+str(dayofyear)+'0000_01D_02H_GIM.INX' |
513 - | get_file=igs_file + ('.Z' if gpsweek<2238 else '.gz') |
514 - | |
515 - | if len(glob.glob(igs_file))<1: # file does not yet exist |
516 - | print('Attempting retrieval of IGS Rapid product file: '+str(get_file)) |
509 + | #ftps login and navigation |
510 + | try: |
511 + | ftps = ftplib.FTP_TLS(host = 'gdc.cddis.eosdis.nasa.gov') # ftp-ssl version |
512 + | ftps.login(user='anonymous', passwd='casa-feedback@nrao.edu') |
513 + | ftps.prot_p() |
514 + | ftps.cwd('gps/products/ionex/'+str(year)+'/'+str(dayofyear)+'/') |
515 + | |
516 + | if len(glob.glob(igs_file))<1: # file does not yet exist locally |
517 + | print('Attempting retrieval of IGS Final product file: '+str(get_file)) |
517 518 | get_path = file_location+get_file |
518 519 | if test_IONEX_connection(get_path): |
519 - | os.system(curlcmd+get_path+' > '+workDir+get_file) |
520 + | with open(get_file, 'wb') as fp: |
521 + | ftps.retrbinary("RETR " + get_file, fp.write) |
520 522 | os.system('gunzip '+get_file) |
523 + | print('FILENAME: ', get_file) |
521 524 | else: |
522 - | #igs_file = igs_file.replace('igr','jpr') |
523 - | igs_file= 'jprg'+str(dayofyear)+'0.'+str(year)[2:4]+'i' if (gpsweek<2274.5) else 'JPL0OPSRAP_'+str(year)+str(dayofyear)+'0000_01D_02H_GIM.INX' |
524 - | get_file=igs_file + ('.Z' if gpsweek<2274.5 else '.gz') |
525 - | |
525 + | # IGS final product file does not exist; try rapid product file |
526 + | igs_file='igrg'+str(dayofyear)+'0.'+str(year)[2:4]+'i' if ( gpsweek<2238) else 'IGS0OPSRAP_'+str(year)+str(dayofyear)+'0000_01D_02H_GIM.INX' |
527 + | get_file=igs_file + ('.Z' if gpsweek<2238 else '.gz') |
528 + | |
526 529 | if len(glob.glob(igs_file))<1: # file does not yet exist |
527 - | print('Attempting retrieval of JPL Rapid product file: '+str(igs_file)) |
530 + | print('Attempting retrieval of IGS Rapid product file: '+str(get_file)) |
528 531 | get_path = file_location+get_file |
529 532 | if test_IONEX_connection(get_path): |
530 - | os.system(curlcmd+get_path+' > '+workDir+get_file) |
533 + | with open(get_file, 'wb') as fp: |
534 + | ftps.retrbinary("RETR " + get_file, fp.write) |
531 535 | os.system('gunzip '+get_file) |
536 + | |
532 537 | else: |
533 - | print('\nNo data products available. You may try to manually'+\ |
534 - | ' download the products at:\n'+\ |
535 - | 'ftp://gdc.cddis.eosdis.nasa.gov/gps/products/ionex\n') |
536 - | return 0,0,0,0,0,0,0,0,[0],'' |
538 + | #igs_file = igs_file.replace('igr','jpr') |
539 + | igs_file= 'jprg'+str(dayofyear)+'0.'+str(year)[2:4]+'i' if (gpsweek<2274.5) else 'JPL0OPSRAP_'+str(year)+str(dayofyear)+'0000_01D_02H_GIM.INX' |
540 + | get_file=igs_file + ('.Z' if gpsweek<2274.5 else '.gz') |
541 + | |
542 + | if len(glob.glob(igs_file))<1: # file does not yet exist |
543 + | print('Attempting retrieval of JPL Rapid product file: '+str(igs_file)) |
544 + | get_path = file_location+get_file |
545 + | if test_IONEX_connection(get_path): |
546 + | with open(get_file, 'wb') as fp: |
547 + | ftps.retrbinary("RETR " + get_file, fp.write) |
548 + | os.system('gunzip '+get_file) |
549 + | else: |
550 + | print('\nNo data products available. You may try to manually'+\ |
551 + | ' download the products at:\n'+\ |
552 + | 'ftp://gdc.cddis.eosdis.nasa.gov/gps/products/ionex\n') |
553 + | return 0,0,0,0,0,0,0,0,[0],'' |
554 + | else: |
555 + | print('JPL Rapid product file: '+igs_file+' already available in current working directory.') |
537 556 | else: |
538 - | print('JPL Rapid product file: '+igs_file+' already available in current working directory.') |
557 + | print('IGS Rapid product file: '+igs_file+' already available in current working directory.') |
539 558 | else: |
540 - | print('IGS Rapid product file: '+igs_file+' already available in current working directory.') |
559 + | print('IGS Final product file: '+igs_file+' already available in current working directory.') |
560 + | |
561 + | ftps.quit() |
562 + | except ftplib.all_errors as e: |
563 + | print('Failed to connect to ftp://gdc.cddis.eosdis.nasa.gov/gps/products/ionex/ with error: ', e) |
541 564 | else: |
542 565 | print('IGS Final product file: '+igs_file+' already available in current working directory.') |
543 566 | |
544 567 | if igs_file.startswith('igs') or igs_file.startswith('IGS0OPSFIN'): |
545 568 | tec_type = 'IGS_Final_Product' |
546 569 | elif igs_file.startswith('igr') or igs_file.startswith('IGS0OPSRAP'): |
547 570 | tec_type = 'IGS_Rapid_Product' |
548 571 | elif igs_file.startswith('jpr') or igs_file.startswith('JPL0OPSRAP'): |
549 - | tec_type = 'JPL_Rapid_Product' |
572 + | tec_type = 'JPL_Rapid_Product' |
550 573 | |
551 574 | ## ========================================================================= |
552 575 | ## |
553 576 | ## The following section reads the lines of the ionex file for 1 day |
554 577 | ## (13 maps total) into an array a[]. It also retrieves the thin-shell |
555 578 | ## ionosphere height used by IGS, the lat./long. spacing, etc. for use |
556 579 | ## later in this script. |
557 580 | ## |
558 581 | ## ========================================================================= |
559 582 | print('Transforming IONEX data format to a TEC/DTEC array for '+str(ymd_date)) |