.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/traj_hk.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_traj_hk.py: TrajGRU (Hong Kong) ======================================================== This example demonstrates how to use pre-trained TrajGRU model to predict the next two-hour forecast using thirty-minute radar echo maps. .. GENERATED FROM PYTHON SOURCE LINES 10-13 Definitions -------------------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 15-16 Import all required modules and methods: .. GENERATED FROM PYTHON SOURCE LINES 16-181 .. code-block:: default # Python package to allow system command line functions import os # Python package to manage warning message import warnings # Python package for datetime function import datetime # Python package for xarrays to read and handle netcdf data import xarray as xr # Python package for deep learning import mxnet as mx import mxnet.ndarray as nd # Python package for projection description from pyresample.geometry import AreaDefinition from pyresample import get_area_def # Python package for projection import cartopy.crs as ccrs # Python package for land/sea features import cartopy.feature as cfeature # Python package for reading map shape file import cartopy.io.shapereader as shpreader # Python package for creating plots from matplotlib import pyplot as plt # Python package for colorbars from matplotlib.colors import BoundaryNorm, ListedColormap # swirlspy regrid function from swirlspy.core.resample import grid_resample # swirlspy standardize data function from swirlspy.utils import standardize_attr, FrameType # swirlspy cuda checking function from swirlspy.utils import locate_cuda # swirlspy deep learning function from swirlspy.qpf.dl.config import cfg, cfg_from_file from swirlspy.qpf.dl.utils import parse_ctx from swirlspy.qpf.dl.hko_benchmark import * from swirlspy.qpf.dl.encoder_forecaster import EncoderForecasterBaseFactory, encoder_forecaster_build_networks, EncoderForecasterStates, load_encoder_forecaster_params from swirlspy.qpf.dl.hko_evaluation import rainfall_to_pixel, pixel_to_dBZ from swirlspy.qpf.dl.operators import * from swirlspy.qpf.dl.ops import * # directory constants from swirlspy.tests.samples import DATA_DIR from swirlspy.tests.outputs import OUTPUT_DIR from swirlspy.tests.models import DL_MODEL_DIR warnings.filterwarnings("ignore") def get_loss_weight_symbol(data, mask, seq_len): if cfg.MODEL.USE_BALANCED_LOSS: balancing_weights = cfg.HKO.EVALUATION.BALANCING_WEIGHTS weights = mx.sym.ones_like(data) * balancing_weights[0] thresholds = [rainfall_to_pixel(ele) for ele in cfg.HKO.EVALUATION.THRESHOLDS] for i, threshold in enumerate(thresholds): weights = weights + \ (balancing_weights[i + 1] - balancing_weights[i]) * (data >= threshold) weights = weights * mask else: weights = mask if cfg.MODEL.TEMPORAL_WEIGHT_TYPE == "same": return weights elif cfg.MODEL.TEMPORAL_WEIGHT_TYPE == "linear": upper = cfg.MODEL.TEMPORAL_WEIGHT_UPPER assert upper >= 1.0 temporal_mult = 1 + \ mx.sym.arange(start=0, stop=seq_len) * \ (upper - 1.0) / (seq_len - 1.0) temporal_mult = mx.sym.reshape( temporal_mult, shape=(seq_len, 1, 1, 1, 1)) weights = mx.sym.broadcast_mul(weights, temporal_mult) return weights elif cfg.MODEL.TEMPORAL_WEIGHT_TYPE == "exponential": upper = cfg.MODEL.TEMPORAL_WEIGHT_UPPER assert upper >= 1.0 base_factor = np.log(upper) / (seq_len - 1.0) temporal_mult = mx.sym.exp(mx.sym.arange( start=0, stop=seq_len) * base_factor) temporal_mult = mx.sym.reshape( temporal_mult, shape=(seq_len, 1, 1, 1, 1)) weights = mx.sym.broadcast_mul(weights, temporal_mult) return weights else: raise NotImplementedError class HKONowcastingFactory(EncoderForecasterBaseFactory): def __init__(self, batch_size, in_seq_len, out_seq_len, ctx_num=1, name="hko_nowcasting"): super(HKONowcastingFactory, self).__init__(batch_size=batch_size, in_seq_len=in_seq_len, out_seq_len=out_seq_len, ctx_num=ctx_num, height=cfg.HKO.ITERATOR.HEIGHT, width=cfg.HKO.ITERATOR.WIDTH, name=name) self._central_region = cfg.HKO.EVALUATION.CENTRAL_REGION def _slice_central(self, data): """Slice the central region in the given symbol Parameters ---------- data : mx.sym.Symbol Returns ------- ret : mx.sym.Symbol """ x_begin, y_begin, x_end, y_end = self._central_region return mx.sym.slice(data, begin=(0, 0, 0, y_begin, x_begin), end=(None, None, None, y_end, x_end)) def _concat_month_code(self): # TODO raise NotImplementedError def loss_sym(self, pred=mx.sym.Variable('pred'), mask=mx.sym.Variable('mask'), target=mx.sym.Variable('target')): """Construct loss symbol. Optional args: pred: Shape (out_seq_len, batch_size, C, H, W) mask: Shape (out_seq_len, batch_size, C, H, W) target: Shape (out_seq_len, batch_size, C, H, W) """ self.reset_all() weights = get_loss_weight_symbol( data=target, mask=mask, seq_len=self._out_seq_len) mse = weighted_mse(pred=pred, gt=target, weight=weights) mae = weighted_mae(pred=pred, gt=target, weight=weights) gdl = masked_gdl_loss(pred=pred, gt=target, mask=mask) avg_mse = mx.sym.mean(mse) avg_mae = mx.sym.mean(mae) avg_gdl = mx.sym.mean(gdl) global_grad_scale = cfg.MODEL.NORMAL_LOSS_GLOBAL_SCALE if cfg.MODEL.L2_LAMBDA > 0: avg_mse = mx.sym.MakeLoss(avg_mse, grad_scale=global_grad_scale * cfg.MODEL.L2_LAMBDA, name="mse") else: avg_mse = mx.sym.BlockGrad(avg_mse, name="mse") if cfg.MODEL.L1_LAMBDA > 0: avg_mae = mx.sym.MakeLoss(avg_mae, grad_scale=global_grad_scale * cfg.MODEL.L1_LAMBDA, name="mae") else: avg_mae = mx.sym.BlockGrad(avg_mae, name="mae") if cfg.MODEL.GDL_LAMBDA > 0: avg_gdl = mx.sym.MakeLoss(avg_gdl, grad_scale=global_grad_scale * cfg.MODEL.GDL_LAMBDA, name="gdl") else: avg_gdl = mx.sym.BlockGrad(avg_gdl, name="gdl") loss = mx.sym.Group([avg_mse, avg_mae, avg_gdl]) return loss .. GENERATED FROM PYTHON SOURCE LINES 182-188 Initialising --------------------------------------------------- This section demonstrates the settings to run the pre-trained model. .. GENERATED FROM PYTHON SOURCE LINES 188-246 .. code-block:: default # set the directory containing the pre-trained model configuration # set the path to '[DL_MODEL_DIR]/hko/convgru' when using ConvGRU base_dir = os.path.join(DL_MODEL_DIR, 'hko/trajgru') cfg_from_file(os.path.join(base_dir, 'cfg0.yml'), target=cfg.MODEL) cfg.MODEL.LOAD_DIR = base_dir cfg.MODEL.SAVE_DIR = base_dir # set the load_iter = 49999 when using ConvGRU cfg.MODEL.LOAD_ITER = 79999 # no training on inference data: mode = 'fixed' and finetune = 0 # train on inference data: mode = 'online' and finetune = 1 mode = 'fixed' finetune = 0 # use GPU if CUDA installed; otherwise use CPU if locate_cuda(): ctx = 'gpu' else: ctx = 'cpu' ctx = parse_ctx(ctx) # Initialize the model hko_nowcasting_online = HKONowcastingFactory(batch_size=1, in_seq_len=cfg.MODEL.IN_LEN, out_seq_len=cfg.MODEL.OUT_LEN) # Initialize the encoder network and forecaster network t_encoder_net, t_forecaster_net, t_loss_net =\ encoder_forecaster_build_networks( factory=hko_nowcasting_online, context=ctx, for_finetune=True) t_encoder_net.summary() t_forecaster_net.summary() t_loss_net.summary() # Load the pre-trained model params to the networks created load_encoder_forecaster_params(load_dir=cfg.MODEL.LOAD_DIR, load_iter=cfg.MODEL.LOAD_ITER, encoder_net=t_encoder_net, forecaster_net=t_forecaster_net) # 2019-04-20 10:30 to 2019-04-20 15:24 pd_path = cfg.HKO_PD.RAINY_EXAMPLE # Set the save directory save_dir = os.path.join(base_dir, "iter%d_%s_finetune%d" % (cfg.MODEL.LOAD_ITER + 1, 'example', finetune)) # Create the environment to run the model env = HKOBenchmarkEnv(pd_path=pd_path, save_dir=save_dir, mode=mode) states = EncoderForecasterStates(factory=hko_nowcasting_online, ctx=ctx[0]) .. rst-class:: sphx-glr-script-out .. code-block:: none Loading YAML config file from <_io.TextIOWrapper name='/tmp/build/docs/swirlspy/swirlspy/tests/models/hko/trajgru/cfg0.yml' mode='r' encoding='UTF-8'> ebrnn1_0_ 96 96 ebrnn2_0_ 32 32 ebrnn3_0_ 16 16 fbrnn1_0_ 96 96 fbrnn2_0_ 32 32 fbrnn3_0_ 16 16 dbrnn1_0_ 96 96 dbrnn2_0_ 32 32 dbrnn3_0_ 16 16 .. GENERATED FROM PYTHON SOURCE LINES 247-255 Inference and Evaluation --------------------------------------------------- The pre-trained model loops through the radar echo maps and generate the inference over the next two hours. Evaluation of the inference is to compare the inference with the truth forecast data. .. GENERATED FROM PYTHON SOURCE LINES 255-289 .. code-block:: default stored_data = [] stored_prediction = [] counter = 0 finetune_iter = 0 while not env.done: states.reset_all() # Get radar data in_frame_dat, in_datetime_clips, out_datetime_clips, begin_new_episode, need_upload_prediction =\ env.get_observation(batch_size=1) in_frame_nd = nd.array(in_frame_dat, ctx=ctx[0]) t_encoder_net.forward(is_train=False, data_batch=mx.io.DataBatch(data=[in_frame_nd] + states.get_encoder_states())) outputs = t_encoder_net.get_outputs() states.update(states_nd=outputs) if need_upload_prediction: counter += 1 t_forecaster_net.forward(is_train=False, data_batch=mx.io.DataBatch( data=states.get_forecaster_state())) # Get inference from model pred_nd = t_forecaster_net.get_outputs()[0] pred_nd = nd.clip(pred_nd, a_min=0, a_max=1) # Generate prediction movies and inference evaluation env.upload_prediction(prediction=pred_nd.asnumpy()) fingerprint = env._fingerprint print("Saving videos to %s" % os.path.join(save_dir, fingerprint)) .. rst-class:: sphx-glr-script-out .. code-block:: none Saving videos to /tmp/build/docs/swirlspy/swirlspy/tests/models/hko/trajgru/iter80000_example_finetune0/hko7_rainy_example_in5_out20_stride5_fixed .. GENERATED FROM PYTHON SOURCE LINES 290-293 Printing of Evaluation -------------------------------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 293-296 .. code-block:: default env.save_eval(print_out=True) .. rst-class:: sphx-glr-script-out .. code-block:: none Saving evaluation result to /tmp/build/docs/swirlspy/swirlspy/tests/models/hko/trajgru/iter80000_example_finetune0/hko7_rainy_example_in5_out20_stride5_fixed Total Sequence Num: 6, Out Seq Len: 20, Use Central: 0 Threshold = 0.5: POD: [0.9949990623241858, 0.9918333506062804, 0.9879483471286296, 0.9835357957698712, 0.9815608919382505, 0.9803458002719322, 0.9763011105374515, 0.9738750651119926, 0.9709277727410094, 0.9684372362233209, 0.9651381657989587, 0.9614454714596993, 0.9560960295113169, 0.9549910786063018, 0.9553564472146474, 0.9534676883380898, 0.9518710296291221, 0.9511550963547714, 0.9490770764569803, 0.9474938728434812] FAR: [0.046496071246717716, 0.0891856138113502, 0.1348375159430183, 0.1700057240984545, 0.19915031863051355, 0.22666186416936807, 0.25931046135649716, 0.287582896713443, 0.31619759778082057, 0.33715963878209315, 0.35671624859028517, 0.3777800672238606, 0.39334781715523626, 0.409948573193853, 0.42273719642447316, 0.4375580141713126, 0.4535575010395478, 0.4638890252081595, 0.4768551123631474, 0.48599792484605797] CSI: [0.9489561700732321, 0.9040345358536194, 0.8561270653903243, 0.8186204064873742, 0.7889799837506463, 0.7615313071442084, 0.7276074360212458, 0.6990573296038771, 0.6700825228206536, 0.6488238699861315, 0.6286757752633123, 0.6070726298397382, 0.5902104000872125, 0.5740865157944366, 0.5620999673384255, 0.5474159699231244, 0.5317505016978439, 0.5217467443563537, 0.5088614295574735, 0.4997668005028998] GSS: [0.943182998181027, 0.8931192648410953, 0.8395793763158353, 0.797794985421462, 0.7645093750142621, 0.7336281527411268, 0.6956896462756536, 0.6633676695748449, 0.6311404476579361, 0.6071824659777206, 0.5844684792499397, 0.5600591866095282, 0.5405972815568172, 0.5224201325828103, 0.5085016657873864, 0.49220862854570174, 0.47473985519984613, 0.46285356382837206, 0.44843375962538284, 0.4377804348867319] HSS: [0.9707608589246828, 0.9435425241589963, 0.9127949433715427, 0.8875260993504582, 0.8665404512323237, 0.8463500683017294, 0.8205388855250004, 0.7976200111481078, 0.7738640146704173, 0.7555862247518292, 0.7377470576462551, 0.7179973573011715, 0.7018022010405319, 0.6863021861074788, 0.6741811127161943, 0.6597048417089044, 0.6438286095353581, 0.6328091550285562, 0.6191981602822678, 0.608967022035211] POD stat: avg 0.967793/final 0.947494 FAR stat: avg 0.317249/final 0.485998 CSI stat: avg 0.669775/final 0.499767 GSS stat: avg 0.630063/final 0.43778 HSS stat: avg 0.762883/final 0.608967 Threshold = 2: POD: [0.9852889083369798, 0.9801533837703389, 0.9740197134270798, 0.9700417257610326, 0.9662766478804214, 0.9624266899364815, 0.9592937913161438, 0.9548222943462756, 0.9502367281152413, 0.9449535753044737, 0.9402282739287144, 0.9343528311460894, 0.9282525729198348, 0.9234230701358797, 0.9242890533696922, 0.9219379178232143, 0.9228487808697865, 0.9212618992338054, 0.9190034168013705, 0.9154884905569706] FAR: [0.017536204782979077, 0.04274378023846637, 0.07551523057602633, 0.10435170424352842, 0.12914534947918585, 0.15433925928567607, 0.18342022491274937, 0.20794070233168802, 0.2340173286019375, 0.25689269599437364, 0.27999693322088476, 0.30362466074421224, 0.32165159847143254, 0.3397055730809674, 0.35321938719877877, 0.36967344020420845, 0.38615144215118175, 0.39719706033156726, 0.4101549395524521, 0.4183332620895595] CSI: [0.9682604689259752, 0.9390545411218015, 0.9022364945750281, 0.8715408267851043, 0.845167397458452, 0.8186337548008458, 0.7892325887428581, 0.7634479223400507, 0.7364410683196839, 0.7122741904942246, 0.6884897360703812, 0.6638930520600628, 0.6445534436950746, 0.6260160440246975, 0.6142385606109637, 0.5983899821109123, 0.5838845124851367, 0.573268077728816, 0.5606966114306636, 0.5520254998551145] GSS: [0.9647095082312246, 0.9322123859498771, 0.8911644447912733, 0.8570114906524168, 0.8275104320050067, 0.797813166628722, 0.7649794848564586, 0.7359488590529116, 0.7059869485631705, 0.6789478278566798, 0.6522955135234061, 0.6247045047258717, 0.6026429343115548, 0.5819447182225841, 0.5683267745936315, 0.550749142348133, 0.5344983834646916, 0.5220620384952848, 0.5080226858850864, 0.4980031333854099] HSS: [0.9820378067999749, 0.9649170999300895, 0.9424505068776613, 0.9230007406699742, 0.9056150022597953, 0.8875373497511865, 0.8668423530358172, 0.8478923272594864, 0.8276580886597898, 0.8087777554391486, 0.7895627727420635, 0.7690069214540338, 0.7520613873612857, 0.7357333180092869, 0.7247555596197616, 0.7103007537559463, 0.6966424848983755, 0.6859931136728131, 0.6737600045942524, 0.664889308021604] POD stat: avg 0.94493/final 0.915488 FAR stat: avg 0.249281/final 0.418333 CSI stat: avg 0.722587/final 0.552025 GSS stat: avg 0.689977/final 0.498003 HSS stat: avg 0.807972/final 0.664889 Threshold = 5: POD: [0.9204261130838569, 0.9123557911660547, 0.9132637148726112, 0.9172258786079924, 0.9174234311415721, 0.9138993623361781, 0.9105009521175307, 0.9039197293394745, 0.8963952529034714, 0.8900604958118284, 0.8847186403416057, 0.8787300121710664, 0.871170845179633, 0.8640429153556507, 0.8584284189788212, 0.8601423854725547, 0.8611418123973515, 0.8609210660675993, 0.8632869058058201, 0.8648629456043175] FAR: [0.04149552276166528, 0.05882221101585701, 0.0825329256566224, 0.10195421904942405, 0.12197396740801476, 0.14529453484771296, 0.16247027624741905, 0.18338840063605116, 0.19980549155871133, 0.21883731075900606, 0.24013689427912635, 0.2612749292147022, 0.2809379960063489, 0.29966386273558876, 0.31812142897896306, 0.33362908963014476, 0.34851954019395626, 0.35710418557125534, 0.36617757149273844, 0.37231818114885656] CSI: [0.885155297776657, 0.8631388549171146, 0.8439307507304565, 0.8307214305889686, 0.8137175748710593, 0.7910102803047112, 0.7738237315106145, 0.7513908480011214, 0.7324529182475822, 0.7124221183800623, 0.6914054433568744, 0.6703807762809275, 0.6499494781984928, 0.6308210253590715, 0.6129493134827015, 0.6012274805636891, 0.5895481224807297, 0.5824082453146912, 0.5760062045881297, 0.571619184272812] GSS: [0.873810091293894, 0.8495994789775003, 0.8282223533958314, 0.8135058492532931, 0.7945724989755995, 0.7694482850387966, 0.7503913955363465, 0.7254491098110649, 0.7047207217697593, 0.6826409884804459, 0.6594526797338668, 0.6360795807695043, 0.613112562936445, 0.5919862380023466, 0.5720597607969308, 0.5591386930419584, 0.5461380194770304, 0.537709106448511, 0.5303740256169345, 0.525163341804218] HSS: [0.9326559776295312, 0.9186848165065205, 0.906041162725584, 0.8971637445649842, 0.8855284469467434, 0.8697041801613612, 0.8573984052365787, 0.8408814907215675, 0.8267873004302453, 0.8113923209453291, 0.7947833497001263, 0.7775655759609619, 0.7601609175002141, 0.7437077329828953, 0.7277837332429833, 0.7172404809620246, 0.7064544207531454, 0.6993638838367845, 0.6931299365240162, 0.6886650464375403] POD stat: avg 0.888146/final 0.864863 FAR stat: avg 0.224723/final 0.372318 CSI stat: avg 0.708704/final 0.571619 GSS stat: avg 0.678179/final 0.525163 HSS stat: avg 0.802755/final 0.688665 Threshold = 10: POD: [0.8797434737276966, 0.8623155925399022, 0.8467330816826867, 0.8365058516790496, 0.8319829676503637, 0.8227749100578758, 0.8187597809076682, 0.8155337892434126, 0.8100978753879208, 0.8045177411294353, 0.8002547821413643, 0.8084692235882972, 0.8052185672805378, 0.8048080160003168, 0.8043121352981901, 0.8016070933776669, 0.8075337965636012, 0.8012230588422772, 0.8012173191312768, 0.7959853005097404] FAR: [0.2281404259725914, 0.25932373619015736, 0.29862196792650497, 0.3310916729370454, 0.36039039767216297, 0.3810763505861242, 0.3909608137013809, 0.4146462217444196, 0.4307004054242975, 0.44239854801391026, 0.45772862152684113, 0.4679911122729055, 0.4877335931272695, 0.5020400161731012, 0.5153804898396996, 0.5255877150855657, 0.5348623853211009, 0.5484265714638555, 0.5500210872122705, 0.5545383177983436] CSI: [0.6981936269535214, 0.6623456328119153, 0.6223649209118732, 0.5915682360117562, 0.5664434214234904, 0.546117607361167, 0.536685131045802, 0.5169135833312514, 0.5022694755538014, 0.4910683781541315, 0.4776239976239976, 0.47246178310679837, 0.4557868311738019, 0.4443011128845355, 0.43350604490500866, 0.42456261727307987, 0.4187184490074019, 0.40607963493922383, 0.40478829085172874, 0.39981343098435] GSS: [0.6826728570275701, 0.6449739655564909, 0.6030081781669192, 0.5707881963747381, 0.5444019313141256, 0.5228770528406542, 0.512971314394788, 0.4922912653057141, 0.47710156116471986, 0.4654397384838655, 0.45115436217901367, 0.44544060861785756, 0.42791364201846643, 0.41559700528456817, 0.4041143638109575, 0.39468033534166413, 0.38853463057562554, 0.3753452114872746, 0.3737072857551623, 0.3685510230943937] HSS: [0.811414832272872, 0.784175286735675, 0.7523457289612514, 0.7267538649603742, 0.7050003244309545, 0.6866963447447326, 0.6780978720670385, 0.6597790615692738, 0.6459969628473172, 0.6352219422757109, 0.621786866976126, 0.6163388602231006, 0.5993550722207225, 0.5871685285192073, 0.5756145998166932, 0.5659796375417836, 0.5596326112724407, 0.5458196361935674, 0.5440857592157645, 0.5386003398851333] POD stat: avg 0.81798/final 0.795985 FAR stat: avg 0.434083/final 0.554538 CSI stat: avg 0.503581/final 0.399813 GSS stat: avg 0.478078/final 0.368551 HSS stat: avg 0.641993/final 0.5386 Threshold = 20: POD: [0.847693200833248, 0.827469052744887, 0.8135791728014992, 0.8020088984764729, 0.778738933567615, 0.7748988404910061, 0.7686053572633426, 0.7588077695914267, 0.7512246207898382, 0.7267514843087363, 0.7134427124194154, 0.7064948278785823, 0.7019424653851324, 0.6901295041266859, 0.6805022554366121, 0.6754052770806485, 0.66625678860539, 0.6494918458993146, 0.6403772827613887, 0.6252173331550087] FAR: [0.27861086893344955, 0.3236272649784157, 0.368755841728113, 0.39816875758802106, 0.42820989951618904, 0.45098653304102726, 0.46439570256482593, 0.4782868984572876, 0.486881879225604, 0.5042582735477898, 0.5171522230943257, 0.5253383767032767, 0.5289992333017635, 0.5397798460712367, 0.551276359600444, 0.5595293389069516, 0.5706739446229696, 0.577527892471229, 0.5759937991363083, 0.5764186204553177] CSI: [0.6386158991510162, 0.5927755554484554, 0.551477759883866, 0.5239815018718343, 0.4918788551426803, 0.4734983066342538, 0.4612316786784355, 0.44750355506399114, 0.4385909548135145, 0.4178565855181024, 0.4044161720064193, 0.39647887323943665, 0.392502114065583, 0.38140621523399837, 0.3706386021525091, 0.36351854556074764, 0.3533375600036229, 0.344034482141899, 0.3424615879943478, 0.3378075658488998] GSS: [0.6283683233911987, 0.5811386413816828, 0.5385984144320107, 0.5104097828146703, 0.47755553634259684, 0.45867847128078054, 0.4459448528913043, 0.43175107519283673, 0.4227397333870319, 0.4016801544266556, 0.3880386148354192, 0.3798340281179515, 0.3756414119179976, 0.36428825770304335, 0.3533324680904233, 0.34614719143839545, 0.3358669253746832, 0.3263262264167268, 0.3246830898193664, 0.32013460589787246] HSS: [0.771776648273991, 0.7350887849706248, 0.7001156499057485, 0.6758560340671448, 0.6464129768342831, 0.6288959223180167, 0.6168213843005079, 0.6031091335268401, 0.5942615131449805, 0.5731409596663073, 0.559117895839561, 0.5505503131214016, 0.5461327474785118, 0.534034146590648, 0.5221665428436544, 0.5142783696165166, 0.5028448852126186, 0.49207535810906344, 0.49020492873301513, 0.48500297540512854] POD stat: avg 0.729952/final 0.625217 FAR stat: avg 0.485244/final 0.576419 CSI stat: avg 0.436201/final 0.337808 GSS stat: avg 0.420558/final 0.320135 HSS stat: avg 0.587094/final 0.485003 Threshold = 30: POD: [0.8079593740283967, 0.7816115283133753, 0.7634045555331587, 0.748710879665084, 0.7305992393873985, 0.7230745147150908, 0.7058912922344359, 0.6756948311648864, 0.6609632769804382, 0.6214223542061338, 0.5952331172904035, 0.5930957570841994, 0.5780462237308441, 0.5638950331791676, 0.5385965792245574, 0.5241141232612448, 0.5171144278606965, 0.49129574678536103, 0.4766483516483517, 0.45601523958384216] FAR: [0.3426644182124789, 0.3981872167526176, 0.4473386596621794, 0.48342667934763467, 0.5111584883600977, 0.5370506481357744, 0.5439073086112657, 0.5657817684950225, 0.5811160743186379, 0.6018567117176223, 0.6215102366391917, 0.6106748302941568, 0.6237816764132553, 0.6361188607020049, 0.648706120980142, 0.6559010765590108, 0.6479712795502269, 0.658531555066685, 0.6543086885362556, 0.6504810752124592] CSI: [0.5685116312987676, 0.5151839464882944, 0.4718432691708715, 0.4402449641259644, 0.41418290941933983, 0.39322909276653706, 0.38326136582414255, 0.35933125134611243, 0.34479910476393477, 0.320423464808068, 0.3010150137449778, 0.3072559678562988, 0.2951596869981903, 0.28396749449381026, 0.27002967359050445, 0.2621835081392691, 0.2649367862969005, 0.25227284270404793, 0.25058028575849794, 0.2466708941027267] GSS: [0.5607126857067353, 0.5063021181523711, 0.4620997694274454, 0.43003435479207536, 0.40359421342590696, 0.3823774915111035, 0.3720010191427285, 0.34779242605426863, 0.33313522363364567, 0.30859673892925893, 0.2892277733360422, 0.29529991078692164, 0.28288534221918565, 0.27164217440832705, 0.25762087891414104, 0.2498922576171006, 0.25259318984640045, 0.23989758258652077, 0.23831065830950024, 0.2346185053837437] HSS: [0.7185341553789301, 0.6722451121205364, 0.6321042914990713, 0.6014322010531022, 0.5750867445382382, 0.5532171839590925, 0.5422751352986123, 0.5160919728157974, 0.49977709346789173, 0.4716452819250702, 0.44868374591035726, 0.45595604280945373, 0.4410142245913254, 0.4272305210932744, 0.4096956137315037, 0.3998620778618409, 0.40331241123444833, 0.3869635459504262, 0.3848964017395019, 0.3800664000428532] POD stat: avg 0.627669/final 0.456015 FAR stat: avg 0.571024/final 0.650481 CSI stat: avg 0.347254/final 0.246671 GSS stat: avg 0.335932/final 0.234619 HSS stat: avg 0.496005/final 0.380066 MSE: [314.61374, 561.646, 862.4934, 1132.2228, 1381.6515, 1651.9943, 1965.568, 2278.7437, 2560.9277, 2812.3176, 3080.173, 3402.1716, 3704.8875, 3981.3147, 4229.9194, 4482.1465, 4771.235, 5021.6626, 5250.041, 5443.3755] MAE: [1497.3224, 2071.152, 2699.64, 3233.034, 3736.8633, 4283.821, 4910.3013, 5534.777, 6078.5493, 6570.014, 7081.8594, 7663.645, 8218.234, 8702.816, 9153.38, 9616.606, 10146.479, 10650.319, 11085.315, 11487.362] Balanced MSE: [1070.8331, 1634.8705, 2263.3179, 2824.0442, 3341.3748, 3771.3533, 4360.8774, 4870.9004, 5281.554, 5947.5796, 6573.6094, 7365.2114, 8027.2095, 8434.544, 8809.261, 9131.528, 9834.229, 10723.3, 10968.574, 11710.198] Balanced MAE: [10147.587, 12380.348, 14265.581, 15578.719, 16676.086, 17684.6, 19159.518, 20519.938, 21459.555, 22900.639, 24149.828, 25751.506, 27186.436, 28011.447, 28775.604, 29444.352, 30696.69, 32317.549, 33036.81, 34442.59] GDL: [1246.7345, 1350.3164, 1375.5074, 1386.1107, 1412.33, 1445.2308, 1488.2333, 1488.806, 1480.857, 1481.2399, 1491.5039, 1517.3276, 1509.305, 1503.5643, 1508.6593, 1511.4779, 1517.622, 1521.2721, 1523.7344, 1523.9825] MSE stat: avg 2944.46/final 5443.38 MAE stat: avg 6721.07/final 11487.4 Balanced MSE stat: avg 6347.22/final 11710.2 Balanced MAE stat: avg 23229.3/final 34442.6 GDL stat: avg 1464.19/final 1523.98 .. GENERATED FROM PYTHON SOURCE LINES 297-300 Generating radar reflectivity maps ----------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 302-304 Step 1: Retrieve observation and prediction for one-hour and two-hour forecasts. .. GENERATED FROM PYTHON SOURCE LINES 304-329 .. code-block:: default # Ground truth at one-hour forecast gt_10 = env._out_frame_dat[9:10, 0, 0] # Ground truth at two-hour forecast gt_20 = env._out_frame_dat[19:, 0, 0] # translate pixel to dBZ gt_10 = pixel_to_dBZ(gt_10) gt_20 = pixel_to_dBZ(gt_20) pred_nd = pixel_to_dBZ(pred_nd.asnumpy()) # Prediction at one-hour forecast pd_10 = pred_nd[9:10, 0, 0] # Prediction at two-hour forecast pd_20 = pred_nd[19:, 0, 0] # Save observation and forecast in list reflec_np = [] reflec_np.append(gt_10) reflec_np.append(pd_10) reflec_np.append(gt_20) reflec_np.append(pd_20) .. GENERATED FROM PYTHON SOURCE LINES 330-331 Step 2: Generate xarrays for one-hour and two-hour forecast. .. GENERATED FROM PYTHON SOURCE LINES 331-406 .. code-block:: default # Get forecast base time basetime = in_datetime_clips[0][0] # Create a datetime list datetime_list = [] curr = basetime for i in range(len(reflec_np)): datetime_list.append(curr) curr = curr + datetime.timedelta(hours=1) y = [] y_init = 255466.666667 for i in range(480): y.append(y_init) y_init -= 1066.666667 x = [] x_init = -255466.666667 for i in range(480): x.append(x_init) x_init += 1066.666667 # Define source grid area_id = 'aeqd' description = 'Azimuthal Equidistant Projection centered at the radar site extending up to 256000.000000m in the x direction and256000.000000m in the y directionwith a 480.000000x480.000000 grid resolution' proj_id = 'aeqd' projection = {'datum': 'WGS84', 'ellps': 'WGS84', 'lat_0': '22.411667', 'lon_0': '114.123333', 'no_defs': 'True', 'proj': 'aeqd', 'units': 'm'} width = 480 height = 480 area_extent = (-256000.00000000003, -256000.00000000003, 256000.00000000003, 256000.00000000003) area_def = AreaDefinition(area_id, description, proj_id, projection, width, height, area_extent) reflec_xr_list = [] for i in range(len(reflec_np)): reflec_xr = xr.DataArray( reflec_np[i], dims=('time', 'y', 'x'), coords={ 'time': [datetime_list[i]], 'y': y, 'x': x }, attrs={ 'long_name': 'Reflectivity', 'units': 'dBZ', 'projection': 'Centered Azimuthal', 'site': (114.12333341315389, 22.41166664287448, 948), 'radar_height': 20, 'area_def': area_def } ) reflec_xr_list.append(reflec_xr) # Defining target grid area_id = "hk1980_250km" description = ("A 1km resolution rectangular grid " "centred at HKO and extending to 250 km " "in each direction in HK1980 easting/northing coordinates") proj_id = 'hk1980' projection = ('+proj=tmerc +lat_0=22.31213333333334 ' '+lon_0=114.1785555555556 +k=1 +x_0=836694.05 ' '+y_0=819069.8 +ellps=intl +towgs84=-162.619,-276.959,' '-161.764,0.067753,-2.24365,-1.15883,-1.09425 +units=m ' '+no_defs') x_size = 500 y_size = 500 area_extent = (587000, 569000, 1087000, 1069000) area_def_tgt = get_area_def( area_id, description, proj_id, projection, x_size, y_size, area_extent ) .. GENERATED FROM PYTHON SOURCE LINES 407-409 Step 3: Reproject the radar data from read_iris_grid() from Centered Azimuthal (source) projection to HK 1980 (target) projection. .. GENERATED FROM PYTHON SOURCE LINES 409-422 .. code-block:: default reproj_xr_list = [] for i in range(len(reflec_xr_list)): reproj_xr = grid_resample( reflec_xr_list[i], area_def, area_def_tgt, coord_label=['easting', 'northing'] ) reproj_xr_list.append(reproj_xr) # Combining the DataArrays reflec_concat = xr.concat(reproj_xr_list, dim='time') standardize_attr(reflec_concat, frame_type=FrameType.dBZ, zero_value=9999.) .. GENERATED FROM PYTHON SOURCE LINES 423-424 Step 4: Generating radar reflectivity maps .. GENERATED FROM PYTHON SOURCE LINES 424-531 .. code-block:: default # Defining colour scale and format levels = [ -32768, 10, 15, 20, 24, 28, 32, 34, 38, 41, 44, 47, 50, 53, 56, 58, 60, 62 ] cmap = ListedColormap([ '#FFFFFF', '#08C5F5', '#0091F3', '#3898FF', '#008243', '#00A433', '#00D100', '#01F508', '#77FF00', '#E0D100', '#FFDC01', '#EEB200', '#F08100', '#F00101', '#E20200', '#B40466', '#ED02F0' ]) norm = BoundaryNorm(levels, ncolors=cmap.N, clip=True) # Defining the crs crs = area_def_tgt.to_cartopy_crs() # Get one-hour forecast date time datetime_one_hr = out_datetime_clips[0][9] # Get two-hour forecast date time datetime_two_hr = out_datetime_clips[0][-1] timelist = [ datetime_one_hr, datetime_one_hr, datetime_two_hr, datetime_two_hr ] # Obtaining the slice of the xarray to be plotted da_plot = reflec_concat # Defining coastlines map_shape_file = os.path.join(DATA_DIR, "shape/rsmc") ocean_color = np.array([[[178, 208, 254]]], dtype=np.uint8) land_color = cfeature.COLORS['land'] coastline = cfeature.ShapelyFeature( list(shpreader.Reader(map_shape_file).geometries()), ccrs.PlateCarree() ) # Plotting p = da_plot.plot( col='time', col_wrap=2, subplot_kws={'projection': crs}, cbar_kwargs={ 'extend': 'max', 'ticks': levels[1:], 'format': '%.3g' }, cmap=cmap, norm=norm ) for idx, ax in enumerate(p.axes.flat): # ocean ax.imshow(np.tile(ocean_color, [2, 2, 1]), origin='upper', transform=ccrs.PlateCarree(), extent=[-180, 180, -180, 180], zorder=-1) # coastline, color ax.add_feature(coastline, facecolor=land_color, edgecolor='none', zorder=0) # overlay coastline without color ax.add_feature(coastline, facecolor='none', edgecolor='gray', linewidth=0.5, zorder=3) ax.gridlines() if idx == 0 or idx == 2: ax.set_title( "Ground Truth\n" f"Based @ {timelist[idx].strftime('%H:%MH')}", loc='left', fontsize=9 ) ax.set_title( '' ) ax.set_title( f"{basetime.strftime('%Y-%m-%d')} \n" f"Valid @ {timelist[idx].strftime('%H:%MH')} ", loc='right', fontsize=9 ) else: ax.set_title( "Reflectivity\n" f"Based @ {basetime.strftime('%H:%MH')}", loc='left', fontsize=9 ) ax.set_title( '' ) ax.set_title( f"{basetime.strftime('%Y-%m-%d')} \n" f"Valid @ {timelist[idx].strftime('%H:%MH')} ", loc='right', fontsize=9 ) plt.savefig( os.path.join(OUTPUT_DIR, 'traj-output-map-hk.png'), dpi=300 ) .. image-sg:: /auto_examples/images/sphx_glr_traj_hk_001.png :alt: Ground Truth Based @ 14:24H, 2019-04-20 Valid @ 14:24H , Reflectivity Based @ 13:00H, 2019-04-20 Valid @ 14:24H , Ground Truth Based @ 15:24H, 2019-04-20 Valid @ 15:24H , Reflectivity Based @ 13:00H, 2019-04-20 Valid @ 15:24H :srcset: /auto_examples/images/sphx_glr_traj_hk_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 1 minutes 53.377 seconds) .. _sphx_glr_download_auto_examples_traj_hk.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: traj_hk.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: traj_hk.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_