#!/usr/bin/env python3 """OPE1 — tải toàn bộ trường EC cần thiết cho một lần khởi tạo. Mười nhóm yêu cầu dưới đây được lấy NGUYÊN VĂN từ quy trình đã chạy của người dùng (`~/Data/QG.25.196/Data_Download/EC_S2S/get_*.py`), chỉ khác ở chỗ chúng được gom vào một chương trình có tham số ngày, và hỗ trợ thêm chế độ THỜI GIAN THỰC. HAI CHẾ ĐỘ reforecast -- tái dự báo, cần cả `date` (phiên bản mô hình) lẫn `hdate` (năm tái dự báo). Đây là chế độ đã dùng để dựng bộ 2080 mẫu huấn luyện. realtime -- dự báo thật, chỉ cần ngày khởi tạo. HAI ĐƯỜNG TRUY CẬP ecds -- ECMWF Data Store qua `cdsapi` (đường mới, người dùng đã thử ở Data_Download/Test_EC_New). webapi -- ECMWF Web API qua `ecmwfapi` (đường cũ, đang được ECMWF ngừng dần). RÀNG BUỘC KHÔNG ĐƯỢC PHÁ: mọi yêu cầu đều là `type = cf`, tức CONTROL FORECAST. OPE1 được huấn luyện trên control; nếu tải trung bình tổ hợp thì z-score sẽ co về 0 và hệ dự báo thiếu một cách hệ thống mà không hề báo lỗi. """ from __future__ import annotations import argparse import datetime as dt import os import sys # Hai miền khác nhau, KHÔNG được dùng chung: # AREA_VN -- lưới dự báo 0.25 độ trên Việt Nam, dùng cho 40 biến đầu vào; # AREA_LS -- miền quy mô lớn châu Á - Thái Bình Dương ở 1 độ, chỉ dùng cho H500, # T850, UV200 vì chúng nuôi phần phân tích EOF/PCA. Tải nhóm này ở # miền nhỏ là sai hẳn cơ sở của các thành phần chính. AREA_VN = [24, 101, 8, 111] GRID_VN = "0.25/0.25" AREA_LS = [45, 80, -15, 150] GRID_LS = "1.0/1.0" LARGE_SCALE = {"h500", "t850", "uv200"} def area_grid(name): return ((AREA_LS, GRID_LS) if name in LARGE_SCALE else (AREA_VN, GRID_VN)) NLT = 46 # 47 giá trị, KHÔNG phải 46: các trường tích luỹ của nhóm Srf_24hour được lấy hiệu # hai bước liên tiếp (procc_Srf_24hour.sh: sub -seltimestep,2/47 -seltimestep,1/46) # nên cần đủ 47 bước mới ra được 46 ngày. Các nhóm mực khí áp chỉ dùng 1/46 nhưng # vẫn tải 47 cho khớp đúng bộ số liệu huấn luyện. STEP_INST = "/".join(str(24 * k) for k in range(NLT + 1)) # 0,24,...,1104 STEP_ACC = "/".join("%d-%d" % (24 * k, 24 * (k + 1)) for k in range(NLT)) STEP_6H = "/".join(str(6 * k) for k in range(1, 4 * NLT + 1)) # 6,12,...,1104 # Ánh xạ mã param của MARS (dùng trong các get_*.py cũ) sang TÊN BIẾN của giao diện # ECDS mới. Giao diện mới không nhận mã param nữa, nên đây là phần bắt buộc phải # dịch; sai một tên là thiếu một biến trong 40 biến đầu vào của OPE1. PARAM2NAME = { "167": "2_m_temperature", "168": "2_m_dewpoint_temperature", "235": "skin_temperature", "59": "convective_available_potential_energy", "136": "total_column_water", "228164": "total_cloud_cover", "228086": "soil_moisture_top_20_cm", "228087": "soil_moisture_top_100_cm", "228095": "soil_temperature_top_20_cm", "228096": "soil_temperature_top_100_cm", "146": "surface_sensible_heat_flux", "147": "surface_latent_heat_flux", "151": "mean_sea_level_pressure", "165": "10_m_u_component_of_wind", "166": "10_m_v_component_of_wind", "169": "surface_solar_radiation_downwards", "175": "surface_thermal_radiation_downwards", "176": "surface_net_solar_radiation", "177": "surface_net_thermal_radiation", "179": "top_net_thermal_radiation", "228143": "convective_precipitation", "228228": "total_precipitation", "121": "maximum_2_m_temperature_in_the_last_6_hours", "122": "minimum_2_m_temperature_in_the_last_6_hours", "134": "surface_pressure", "156": "geopotential_height", "130": "temperature", "131": "u_component_of_wind", "132": "v_component_of_wind", "133": "specific_humidity", "135": "vertical_velocity", } # (tên nhóm, thư mục đích, tiền tố file, param, levtype, levelist, step) GROUPS = [ ("srf_daily", "Srf_daily_ave", "Srf_daily_cf", "167/168/235/59/136/228164/228086/228087/228095/228096", "sfc", None, STEP_ACC), ("srf_24hour", "Srf_24hour", "Srf_24hour_cf", "146/147/151/165/166/169/175/176/177/179/228143/228228", "sfc", None, STEP_INST), # tên tệp phải là Tx_tm_6hours_cf: procc_Tx_tm.sh đọc đúng tiền tố này ("tx_tm", "Tx_tm", "Tx_tm_6hours_cf", "121/122", "sfc", None, STEP_6H), # P_srf và QUV_plev sau đó được ghép thành pquv_cf_*.nc cho calc_IVT.py ("p_srf", "QUV_plev", "P_srf_cf", "134", "sfc", None, STEP_INST), ("h_plev", "H_plev", "H_plev_cf", "156", "pl", "500/700/850/925/1000", STEP_INST), ("t_plev", "T_plev", "T_plev_cf", "130", "pl", "500/700/850/925", STEP_INST), ("quv_plev", "QUV_plev", "QUV_plev_cf", "131/132/133", "pl", "300/500/700/850/925/1000", STEP_INST), ("w_plev", "W_plev", "W_plev_cf", "135", "pl", "700", STEP_INST), ("h500", "H500", "H500_cf", "156", "pl", "500", STEP_INST), ("t850", "T850", "T850_cf", "130", "pl", "850", STEP_INST), ("uv200", "UV200", "UV200_cf", "131/132", "pl", "200", STEP_INST), ] def req_webapi(g, init, hdate): """Yêu cầu theo ECMWF Web API — đúng dạng các get_*.py hiện có.""" _, _, _, param, levtype, levelist, step = g area, grid = area_grid(g[0]) r = {"class": "s2", "dataset": "s2s", "expver": "prod", "model": "glob", "origin": "ecmf", "levtype": levtype, "param": param, "step": step, "time": "00:00:00", "area": "%d/%d/%d/%d" % tuple(area), "grid": grid, "type": "cf", "format": "netcdf", "date": init.isoformat(), "stream": "enfh" if hdate else "enfo"} if hdate: r["hdate"] = hdate.isoformat() if levelist: r["levelist"] = levelist return r def req_ecds(g, init, hdate, realtime): """Yêu cầu theo ECMWF Data Store. Giao diện của bộ `s2s-forecasts` (thời gian thực) KHÁC bộ `s2s-reforecasts`: tên biến mô tả thay cho mã param, `level_type` nhận "pressure" chứ không phải "pressure_level", mức nằm ở `level_value` dạng "850_hpa", và định dạng chỉ có grib. Bản tái dự báo giữ nguyên dạng cũ để còn dựng lại được bộ huấn luyện. """ _, _, _, param, levtype, levelist, step = g area, grid = area_grid(g[0]) lead = [s_.replace("-", "_") for s_ in step.split("/")] names = [PARAM2NAME[p] for p in param.split("/")] r = {"origin": "ecmwf", "year": "%04d" % init.year, "month": "%02d" % init.month, "day": "%02d" % init.day, "time": "00:00", "forecast_type": "control_forecast", "leadtime_hour": lead, "area": area, "grid": grid} if realtime: r["level_type"] = "single_level" if levtype == "sfc" else "pressure" r["variable"] = names r["data_format"] = "grib" if levelist: r["level_value"] = ["%s_hpa" % L for L in levelist.split("/")] else: r["level_type"] = ("single_level" if levtype == "sfc" else "pressure_level") r["variable"] = param.split("/") r["data_format"] = "grib2" r["hyear"] = ["%04d" % hdate.year] r["hmonth"] = ["%02d" % hdate.month] r["hday"] = ["%02d" % hdate.day] if levelist: r["pressure_level"] = levelist.split("/") return r def fetch(route, dataset, request, target): if os.path.exists(target) and os.path.getsize(target) > 0: print(" đã có, bỏ qua: %s" % os.path.basename(target), flush=True) return True tmp = target + ".part" if os.path.exists(tmp): os.remove(tmp) try: if route == "ecds": import cdsapi cdsapi.Client(url="https://ecds.ecmwf.int/api").retrieve( dataset, request).download(tmp) else: from ecmwfapi import ECMWFDataServer request = dict(request, target=tmp) ECMWFDataServer().retrieve(request) except Exception as e: print(" LỖI: %s" % e, flush=True) return False if not os.path.exists(tmp) or os.path.getsize(tmp) == 0: print(" LỖI: tệp rỗng", flush=True) return False os.replace(tmp, target) print(" xong: %s (%.1f MB)" % (os.path.basename(target), os.path.getsize(target) / 1e6), flush=True) return True def main(): ap = argparse.ArgumentParser() ap.add_argument("--init", required=True, help="ngày khởi tạo YYYY-MM-DD") ap.add_argument("--hdate", default=None, help="ngày tái dự báo YYYY-MM-DD; có thì chạy chế độ " "reforecast, không có thì chạy chế độ realtime") ap.add_argument("--route", default="ecds", choices=["ecds", "webapi"]) ap.add_argument("--dataset", default=None, help="tên bộ dữ liệu trên ECDS; mặc định " "'s2s-reforecasts' cho reforecast và 's2s-forecasts' " "cho thời gian thực") ap.add_argument("--root", default="/home/tanpv/Data/QG.25.196/EC") ap.add_argument("--only", nargs="*", default=None, help="chỉ tải một số nhóm, ví dụ: srf_daily tx_tm") args = ap.parse_args() init = dt.date.fromisoformat(args.init) hdate = dt.date.fromisoformat(args.hdate) if args.hdate else None mode = "reforecast" if hdate else "realtime" dataset = args.dataset or ("s2s-reforecasts" if hdate else "s2s-forecasts") groups = [g for g in GROUPS if args.only is None or g[0] in args.only] print("[dl] chế độ %s | khởi tạo %s%s | đường %s | bộ %s" % (mode, init, (" | hdate %s" % hdate) if hdate else "", args.route, dataset), flush=True) print("[dl] %d nhóm trường, type=cf (CONTROL FORECAST); miền Việt Nam %s @ %s, " "miền quy mô lớn %s @ %s" % (len(groups), AREA_VN, GRID_VN, AREA_LS, GRID_LS), flush=True) tag = (hdate or init).isoformat() ok = 0 for g in groups: name, sub, prefix = g[0], g[1], g[2] out_dir = os.path.join(args.root, sub, str((hdate or init).year)) os.makedirs(out_dir, exist_ok=True) ext = "grib" if (args.route == "ecds" and hdate is None) else "nc" target = os.path.join(out_dir, "%s_%s.%s" % (prefix, tag, ext)) print("[dl] %-11s -> %s" % (name, target), flush=True) r = (req_ecds(g, init, hdate, hdate is None) if args.route == "ecds" else req_webapi(g, init, hdate)) ok += bool(fetch(args.route, dataset, r, target)) print("[dl] %d/%d nhóm tải xong" % (ok, len(groups)), flush=True) return 0 if ok == len(groups) else 1 if __name__ == "__main__": sys.exit(main())