# import os from os import chdir, system from pathlib import Path from datetime import datetime, timedelta # from os.path import exists from re import sub import numpy as np # from glob import glob # from tools import get_sat_freq, read_gnss_v2, dt2float, lagrange, inv_quat_rotate # from tools import calculate_tec, calc_bore_sight from traceback import extract_tb from step0 import prepare_file from step1 import gather_data from step4 import calc_tec from sys import argv from process_dcb import process_dcb import multiprocessing as mp if __name__ == '__main__': if len(argv) == 3: working_leo = argv[2] nProcess = int(argv[1]) elif len(argv) == 2: working_leo = '' nProcess = int(argv[1]) else: working_leo = '' nProcess = 0 working_dir = Path(__file__).absolute().parent if '__file__' in dir() else Path.cwd() data_dirs = [# Path('/data2/xinjiaz/data_GPSRO/SCDR_unpack_STAR/cosmic2/v1.26/level1a/leoOrb/'), Path('/data2/mchen36/POD_DATA/cosmic2/LEOPOD/leoOrb/'), Path('/data2/xinjiaz/data_GPSRO/SCDR_unpack_STAR/cosmic2/v1.26/level1a/leoAtt/'), Path('/data2/xinjiaz/data_GPSRO/SCDR_unpack_STAR/cosmic2/v1.26/level1a/podCrx/')] if 'ucar' in working_leo: data_dirs = [Path('/GNSSRO_data2/xinjiaz/data_GPSRO/UCAR/cosmic2/level1b/leoOrb/'), Path('/GNSSRO_data2/xinjiaz/data_GPSRO/UCAR/cosmic2/level1a/leoAtt/'), Path('/GNSSRO_data2/xinjiaz/data_GPSRO/UCAR/cosmic2/level1a/podCrx/')] elif 'tacc' in working_leo: data_dirs = [Path('/GNSSRO_data2/xinjiaz/data_GPSRO/TACC/level1b/leoOrb/'), Path('/GNSSRO_data2/xinjiaz/data_GPSRO/TACC/level1a/leoAtt/'), Path('/GNSSRO_data2/xinjiaz/data_GPSRO/TACC/level1a/podCrx/')] root_out = working_dir / f'L0_extract_{working_leo}' root_out.mkdir(exist_ok=True) ztypes = ['leoOrb',' leoAtt', 'podCrx'] chdir(working_dir) (working_dir / 'log').mkdir(exist_ok=True) logfile = working_dir / 'log' / f'progress_rt_{working_leo}_{datetime.now():%Y%m}.txt' def log_message(msg): with open(logfile, 'a') as fid: fid.write(f"{datetime.now():%H:%M:%S} {msg}\n") log_message('- checking progress...') progfile = working_dir / f'progress_rt_{working_leo}.npz' if progfile.is_file(): with np.load(progfile, allow_pickle=True) as fid: tarfiles = fid['tarfiles'] filetypes = fid['filetypes'] fileleos = fid['fileleos'] filedates = fid['filedates'] progress = fid['progress'] timetag = fid['timetag'] lastaccess = fid['lastaccess'] else: tarfiles = np.array([],dtype=object) filetypes = np.array([],dtype='int8') fileleos = np.array([],dtype='int8') filedates = np.array([], dtype=datetime) progress = np.array([],dtype='int8') timetag = np.array([],dtype=datetime) lastaccess = np.array([datetime(2000, 1, 1)]) def update_progress(pnum,state): progress[pnum] = state timetag[pnum] = datetime.now() lastaccess[0] = datetime.now() def save_progress(): np.savez(f'progress_rt_{working_leo}.npz', tarfiles = tarfiles, filetypes = filetypes, fileleos = fileleos, filedates = filedates, progress = progress, timetag = timetag, lastaccess = lastaccess, ) # idx = np.logical_and(progress == 0, datetime.now() - timetag > timedelta(minutes=60)) # update_progress(idx, -1) # if np.sum(idx) > 0: # log_message(f"-- {np.sum(idx)} timeout processes were terminated.") # idx = progress == -2 # update_progress(idx, -1) # if np.sum(idx) > 0: # log_message(f"-- {np.sum(idx)} error processes were terminated.") if datetime.now() - lastaccess[0] < timedelta(minutes=60): log_message("-- Still progressing.") else: idx = np.logical_and(progress == 0, datetime.now() - timetag > timedelta(minutes=60)) update_progress(idx, -1) if np.sum(idx) > 0: log_message(f"-- {np.sum(idx)} timeout processes were terminated.") idx = progress == -2 update_progress(idx, -1) if np.sum(idx) > 0: log_message(f"-- {np.sum(idx)} error processes were terminated.") save_progress() manual_date = len(working_leo)>0 # set True to manually set dates to process nloop = 0 while True: new_proc = False nloop += 1 if manual_date: # manually assign dates to process process_dates = [datetime(2026, 1, 6) + timedelta(days=x) for x in range(2)] + [ datetime(2026, 3, 25) + timedelta(days=x) for x in range(2)] else: dt_today = datetime.now().replace(hour=0, minute=0, second=0, microsecond=0) process_dates = [dt_today + timedelta(days=x) for x in range(-2,2)] # 3 recent days if 'ucar' in working_leo: # idx = int(sub(r'[^\d]', r'', working_leo)) root_files = [np.hstack([list((data_dirs[0] / f'{x:%Y/%j}').glob('*_sp3')) for x in process_dates]), np.hstack([list((data_dirs[1] / f'{x:%Y/%j}').glob(f'leoAtt_*.*.*.*_txt')) for x in process_dates]), np.hstack([list((data_dirs[2] / f'{x:%Y/%j}').glob(f'podCrx_*.*.*.*.*_crx.gz')) for x in process_dates])] root_types = [np.full((len(root_files[x]),), x) for x in range(len(root_files))] root_leos = [[int(sub(r'.*_\d{4}\.\d{3}\.\d\d(\d)\.\d\d_.*', r'\1', x.name)) for x in root_files[0]], [int(sub(r'.*_\d{4}\.\d{3}\.\d\d(\d)\.\d\d_.*', r'\1', x.name)) for x in root_files[1]], [int(sub(r'.*_\d{4}\.\d{3}\.\d\d(\d)\.\d\d.\d\d_.*', r'\1', x.name)) for x in root_files[2]]] root_dates = [[datetime.strptime(x.parent.as_posix()[-8:], '%Y/%j') for x in root_files[0]], [datetime.strptime(x.parent.as_posix()[-8:], '%Y/%j') for x in root_files[1]], [datetime.strptime(x.parent.as_posix()[-8:], '%Y/%j') for x in root_files[2]]] elif 'tacc' in working_leo: root_files = [np.hstack([list((data_dirs[0] / f'{x:%Y/%j}').glob('leoOrb_*_sp3')) for x in process_dates]), np.hstack([list((data_dirs[1] / f'{x:%Y/%j}').glob('leoAtt_*_txt')) for x in process_dates]), np.hstack([list((data_dirs[2] / f'{x:%Y/%j}').glob('podCrx_*_crx.gz')) for x in process_dates])] root_types = [np.full((len(root_files[x]),), x) for x in range(len(root_files))] root_leos = [[int(sub(r'.*_\d{4}\.\d{3}\.\d\d(\d)\.\d\d_.*', r'\1', x.name)) for x in root_files[0]], [int(sub(r'.*_\d{4}\.\d{3}\.\d\d(\d)\.\d\d_.*', r'\1', x.name)) for x in root_files[1]], [int(sub(r'.*_\d{4}\.\d{3}\.\d\d(\d)\.\d\d.\d\d_.*', r'\1', x.name)) for x in root_files[2]]] root_dates = [[datetime.strptime(x.parent.as_posix()[-8:], '%Y/%j') for x in root_files[0]], [datetime.strptime(x.parent.as_posix()[-8:], '%Y/%j') for x in root_files[1]], [datetime.strptime(x.parent.as_posix()[-8:], '%Y/%j') for x in root_files[2]]] else: root_files = [np.hstack([list((data_dirs[0] / f'{x:%Y.%j}').glob('*1.SP3')) for x in process_dates]), np.hstack([list((data_dirs[1] / f'{x:%Y/%j}').glob('*_txt')) for x in process_dates]), np.hstack([list((data_dirs[2] / f'{x:%Y/%j}').glob('*_rnx')) for x in process_dates])] root_types = [np.full((len(root_files[x]),),x) for x in range(len(root_files))] root_leos = [[int(sub(r'.*_C2(\d)\d\.SP3',r'\1',x.name)) for x in root_files[0]], [int(sub(r'.*_\d{4}\.\d{3}\.\d\d(\d)\.\d\d_txt',r'\1',x.name)) for x in root_files[1]], [int(sub(r'.*\.\d\d\.\d\d\.\d\d(\d)\.\d\d\..*',r'\1',x.name)) for x in root_files[2]]] root_dates = [[datetime.strptime(x.parent.as_posix()[-8:],'%Y.%j') for x in root_files[0]], [datetime.strptime(x.parent.as_posix()[-8:],'%Y/%j') for x in root_files[1]], [datetime.strptime(x.parent.as_posix()[-8:],'%Y/%j') for x in root_files[2]]] root_files = np.hstack(root_files) root_types = np.hstack(root_types).astype('int8') root_leos = np.hstack(root_leos).astype('int8') root_dates = np.hstack(root_dates).astype(datetime) if 'tacc' in working_leo: idx = np.logical_and(root_types == 2, root_leos == 1) root_files = root_files[idx] root_types = root_types[idx] root_leos = root_leos[idx] root_dates = root_dates[idx] idx = np.array([tarfile in root_files for tarfile in tarfiles], dtype=bool) tarfiles = tarfiles[idx] filetypes = filetypes[idx] fileleos = fileleos[idx] filedates = filedates[idx] progress = progress[idx] timetag = timetag[idx] if np.any(idx==False): log_message(f'- {np.sum(idx==False):d} files are removed from the queue.') idx = np.array([root_file not in tarfiles for root_file in root_files], dtype=bool) tarfiles = np.concatenate((tarfiles, root_files[idx])) filetypes = np.concatenate((filetypes, root_types[idx])) fileleos = np.concatenate((fileleos, root_leos[idx])) filedates = np.concatenate((filedates, root_dates[idx])) progress = np.concatenate((progress, np.full((np.sum(idx),),-1,dtype='int8'))) timetag = np.concatenate((timetag, np.full((np.sum(idx),),datetime.now(),dtype=datetime))) if np.any(idx): log_message(f'- {np.sum(idx == True):d} files are added to the queue.') idx = np.argsort(tarfiles) tarfiles = tarfiles[idx] filetypes = filetypes[idx] fileleos = fileleos[idx] filedates = filedates[idx] progress = progress[idx] timetag = timetag[idx] save_progress() try: idx = np.where(progress==-1)[0] uniq_leos = np.unique(fileleos[idx]) for dt in process_dates: for leo in uniq_leos: idx_now = idx[np.logical_and(fileleos[idx] == leo, filedates[idx] == dt)] if len(idx_now)==0: continue log_message(f"-- Processing {dt:%Y-%m-%d} C2E{leo} files ({len(idx_now)}/{np.sum(progress==-1)}).") update_progress(idx_now, 0) save_progress() if True: # 'repr' not in working_leo: system('rm -r ' + root_out.as_posix() + '/*') (root_out / f'{dt:%Y/%j}').mkdir(exist_ok=True,parents=True) idx_gz = idx_now[filetypes[idx_now] == 2] idx_gz = idx_gz[[x.name.endswith('.gz') for x in tarfiles[idx_gz]]] for ii in idx_gz: tar_out = root_out / f'{dt:%Y/%j}' / (tarfiles[ii].stem[:-4] + '.crx') system('gunzip -c '+ tarfiles[ii].as_posix() + '> ' + tar_out.as_posix()) system('./crx2rnx ' + tar_out.as_posix()) system('rm ' + tar_out.as_posix()) gather_data(working_leo, tarfiles[idx_now], filetypes[idx_now], f'C2E{leo}') if True: # 'ucar' not in working_leo: calc_tec(f'C2E{leo}', [dt], nProcess) update_progress(idx_now, 1) if 'repr' in working_leo: log_message(f'-- Processing DCB of {dt:%Y-%m-%d} ...') leos = ['C2E1', 'C2E2', 'C2E4', 'C2E5', 'C2E6'] ants = [1, 2] leos, ants = np.meshgrid(leos, ants, indexing='ij') leos = leos.ravel() ants = ants.ravel() dt0 = [dt] * len(leos) dur = [1] * len(leos) allargs = list(zip(leos, ants, dt0, dur)) with mp.Pool(processes=len(allargs)) as pool: pool.map(process_dcb, allargs) # save_progress() if len(idx)>0: save_progress() new_proc = True log_message(f'- Progress completed.') except Exception as e: log_message('- Unexpected error occurs.') log_message(f"-- {type(e).__name__} >> {str(e)}") for tb in extract_tb(e.__traceback__): log_message(f"-- line {tb.lineno} of {tb.filename}") update_progress(progress == 0, -2) save_progress() new_proc = True if len(working_leo)>0 or not new_proc: if nloop > 1: log_message(f'- No new item, exit program.') lastaccess[0] = datetime(2000, 1, 1) save_progress() break # extract_orb() # extract_att() # gather_crx() # leos_ref = ['GN04', 'GN05'] # ants_ref = [0, 1] # gnsses_ref = ['C06', 'C07', 'C08', 'C09', 'C10', 'C11', 'C12', 'C13', 'C14', 'C16', 'C19', 'C20', 'C21', 'C22', # 'C23', 'C24', 'C25', 'C26', 'C27', 'C28', 'C29', 'C30', 'C32', 'C33', 'C34', 'C35', 'C36', 'C37', # 'C38', 'C39', 'C40', 'C41', 'C42', 'C43', 'C44', 'C45', 'C47', 'C48', 'C49', 'E02', 'E03', 'E04', # 'E05', 'E06', 'E07', 'E08', 'E09', 'E10', 'E11', 'E12', 'E13', 'E14', 'E15', 'E16', 'E18', 'E19', # 'E21', 'E23', 'E24', 'E25', 'E26', 'E27', 'E29', 'E30', 'E31', 'E33', 'E34', 'E36', 'G01', 'G02', # 'G03', 'G04', 'G05', 'G06', 'G07', 'G08', 'G09', 'G10', 'G11', 'G12', 'G13', 'G14', 'G15', 'G16', # 'G17', 'G18', 'G19', 'G20', 'G21', 'G22', 'G23', 'G24', 'G25', 'G26', 'G27', 'G28', 'G29', 'G30', # 'G31', 'G32', 'J02', 'J03', 'J04', 'R02', 'R03', 'R04', 'R05', 'R06', 'R07', 'R08', 'R09', 'R11', # 'R12', 'R14', 'R15', 'R16', 'R17', 'R18', 'R19', 'R20', 'R21', 'R22', 'R24', 'R25', 'R26', 'R27', # 'R28'] # leos_nd, ants_nd, gnss_nd = np.meshgrid(leos_ref,ants_ref,gnsses_ref,indexing='ij') # leos_nd = leos_nd.flatten() # ants_nd = ants_nd.flatten() # gnss_nd = gnss_nd.flatten() # allargs = zip(leos_nd,ants_nd,gnss_nd) # import multiprocessing as mp # with mp.Pool(processes=24) as pool: # pool.map(calc_bore_sight, allargs) # leos_ref = ['GN04', 'GN05'] # ants_ref = [0, 1] # dts_ref = [datetime(2026,1,1) + timedelta(days=x) for x in range(31)] # leos_nd, ants_nd, dts_nd = np.meshgrid(leos_ref,ants_ref,dts_ref,indexing='ij') # leos_nd = leos_nd.flatten() # ants_nd = ants_nd.flatten() # dts_nd = dts_nd.flatten() # idx = np.full(ants_nd.shape, True) # for nn,leo,dt,ant in zip(range(len(dts_nd)),leos_nd,dts_nd,ants_nd): # outfile = f'./p04_teczdcb/teczdcb_{leo}.{dt:%Y-%m-%d}.A{ant:02d}.npz' # if os.path.exists(outfile) and datetime(1970,1,1)+timedelta(seconds=os.path.getmtime(outfile))>datetime(2026,4,10): # idx[nn] = False # leos_nd = leos_nd[idx] # ants_nd = ants_nd[idx] # dts_nd = dts_nd[idx] # allargs = zip(leos_nd, dts_nd, ants_nd) # import multiprocessing as mp # with mp.Pool(processes=10) as pool: # pool.map(calculate_dcb, allargs) # for dd in range(31): # reffiles = glob(f'./p01_podPair/{datetime(2026,1,1) + timedelta(days=dd):%Y-%m-%d}/*.npz') # for reffile in reffiles: # try: # sec_ref_0, tec, coszn, meps = calculate_tec(reffile) # except Exception as e: # print(f"{datetime.now():%H:%M:%S} - Error detected for {os.path.basename(reffile)}", flush=True) # print(f"{datetime.now():%H:%M:%S} -- {type(e).__name__} >> {str(e)}", flush=True) # for tb in extract_tb(e.__traceback__): # print(f"{datetime.now():%H:%M:%S} -- line {tb.lineno} of {tb.filename}", flush=True) # def rename(): # dts = [datetime(2026,1,1)+timedelta(days=dd) for dd in range(31)] # folder = Path('./p05_podTec/') # for dt in dts: # files = (folder/f'{dt:%Y-%m-%d}').glob('podTec*.nc') # for file in files: # ctime = datetime.fromtimestamp(file.stat().st_mtime) # new_name = sub(r'.*_(GN\d\d)\..*\.(\d{14})\.(\d{14})\.(A\d\d)\.([A-Z]\d\d)\.nc',f'ROSW_TEC_\\1_\\2_\\3_\\5_\\4_c{ctime:%Y%m%d%H%M%S}_nes_ops.nc',file.name) # file.rename(file.with_name(new_name))