guychuk commited on
Commit
7c60349
·
verified ·
1 Parent(s): 4994b2f

Upload train_v2_prod.py with huggingface_hub

Browse files
Files changed (1) hide show
  1. train_v2_prod.py +751 -0
train_v2_prod.py ADDED
@@ -0,0 +1,751 @@
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
1
+ # /// script
2
+ # requires-python = ">=3.10"
3
+ # dependencies = [
4
+ # "torch>=2.1",
5
+ # "numpy",
6
+ # "pandas",
7
+ # "scikit-learn",
8
+ # "huggingface-hub",
9
+ # "trackio",
10
+ # ]
11
+ # ///
12
+ """
13
+ Flight-JEPA v2 — bundled training script for HF Jobs.
14
+
15
+ Self-contained: downloads the dataset from HF, trains either the supervised
16
+ baseline (`--lambda-jepa 0`) or the JEPA-augmented model, runs blindspot
17
+ scoring + extrapolation eval, and pushes the result to a hub repo.
18
+
19
+ Usage (HF Jobs):
20
+ python train_v2_prod.py --tag baseline --lambda-jepa 0.0 \
21
+ --hub-model-id guychuk/flight-jepa-v2 --push-to-hub
22
+
23
+ python train_v2_prod.py --tag jepa --lambda-jepa 0.5 \
24
+ --hub-model-id guychuk/flight-jepa-v2 --push-to-hub
25
+ """
26
+ from __future__ import annotations
27
+ import argparse
28
+ import copy
29
+ import json
30
+ import math
31
+ import os
32
+ import shutil
33
+ import sys
34
+ import time
35
+
36
+ import numpy as np
37
+ import pandas as pd
38
+ import torch
39
+ import torch.nn as nn
40
+ import torch.nn.functional as F
41
+ from torch.utils.data import Dataset, DataLoader
42
+
43
+ try:
44
+ import trackio
45
+ HAS_TRACKIO = True
46
+ except ImportError:
47
+ HAS_TRACKIO = False
48
+
49
+
50
+ # ============================================================================
51
+ # DATA UTILITIES (inlined from flight_jepa.data)
52
+ # ============================================================================
53
+
54
+ def load_atfm(dset_name, mode, path):
55
+ variables = ["X", "Y", "Z"]
56
+ data, labels = [], None
57
+ for var in variables:
58
+ df = pd.read_csv(os.path.join(path, f"{dset_name}_{mode}_{var}.tsv"),
59
+ sep="\t", header=None, na_values="NaN")
60
+ if labels is None:
61
+ labels = df.values[:, 0]
62
+ data.append(df.values[:, 1:])
63
+ return np.stack(data, axis=-1), labels.astype(int)
64
+
65
+
66
+ def compute_features(traj_xyz: np.ndarray) -> np.ndarray:
67
+ if traj_xyz.shape[0] < 2:
68
+ T = traj_xyz.shape[0]
69
+ return np.concatenate([
70
+ traj_xyz, np.zeros((T, 3), dtype=traj_xyz.dtype),
71
+ np.zeros((T, 3), dtype=traj_xyz.dtype)
72
+ ], axis=1)
73
+ x, y, z = traj_xyz[:, 0], traj_xyz[:, 1], traj_xyz[:, 2]
74
+ diff = np.diff(traj_xyz, axis=0)
75
+ norms = np.maximum(np.linalg.norm(diff, axis=1, keepdims=True), 1e-8)
76
+ u = diff / norms
77
+ u = np.vstack([u, u[-1:]])
78
+ r = np.sqrt(x ** 2 + y ** 2)
79
+ theta = np.arctan2(y, x)
80
+ return np.column_stack([
81
+ traj_xyz, u,
82
+ r[:, None], np.sin(theta)[:, None], np.cos(theta)[:, None]
83
+ ]).astype(np.float32)
84
+
85
+
86
+ def ensure_data(airport: str, data_dir: str = "data"):
87
+ target = os.path.join(data_dir, airport)
88
+ if os.path.isdir(target) and any(f.endswith(".tsv") for f in os.listdir(target)):
89
+ return target
90
+ print(f"[data] downloading {airport} from HF ...")
91
+ from huggingface_hub import snapshot_download
92
+ snap = snapshot_download(
93
+ "petchthwr/ATFMTraj",
94
+ repo_type="dataset",
95
+ allow_patterns=[f"{airport}/*"],
96
+ )
97
+ os.makedirs(data_dir, exist_ok=True)
98
+ src = os.path.join(snap, airport)
99
+ if not os.path.isdir(target):
100
+ shutil.copytree(src, target)
101
+ return target
102
+
103
+
104
+ # ============================================================================
105
+ # DATASET — variable-length blindspot
106
+ # ============================================================================
107
+
108
+ PAD_VALUE = 0.0
109
+
110
+
111
+ class BlindspotDataset(Dataset):
112
+ def __init__(self, airport, mode, data_dir,
113
+ past_max=256, past_min=60,
114
+ delta_min=30, delta_max=120,
115
+ seed=0, epoch_multiplier=4):
116
+ ensure_data(airport, data_dir)
117
+ airport_dir = os.path.join(data_dir, airport)
118
+ raw, labels = load_atfm(airport, mode, airport_dir)
119
+
120
+ self.past_max = past_max
121
+ self.past_min = past_min
122
+ self.delta_min = delta_min
123
+ self.delta_max = delta_max
124
+ self.epoch_multiplier = epoch_multiplier
125
+ self.rng_seed = seed
126
+
127
+ lengths = np.array(
128
+ [int(np.sum(~np.isnan(raw[i, :, 0]))) for i in range(raw.shape[0])],
129
+ dtype=np.int64,
130
+ )
131
+ min_required = past_min + delta_max + 1
132
+ keep = lengths >= min_required
133
+ if keep.sum() == 0:
134
+ raise RuntimeError(
135
+ f"No trajectories of length >= {min_required} in {airport}/{mode}"
136
+ )
137
+ raw = raw[keep]
138
+ lengths = lengths[keep]
139
+ self.labels = labels[keep].astype(np.int64)
140
+
141
+ self.positions = []
142
+ for i in range(raw.shape[0]):
143
+ L = int(lengths[i])
144
+ self.positions.append(np.nan_to_num(raw[i, :L], nan=0.0).astype(np.float32))
145
+ del raw
146
+
147
+ self.n_traj = len(self.positions)
148
+ print(f"[data] {airport}/{mode}: {self.n_traj} trajectories "
149
+ f"(after filtering for L >= {min_required})")
150
+
151
+ def __len__(self):
152
+ return self.n_traj * self.epoch_multiplier
153
+
154
+ def __getitem__(self, idx):
155
+ traj_idx = idx % self.n_traj
156
+ rng = np.random.default_rng(self.rng_seed + idx * 9173)
157
+ positions = self.positions[traj_idx]
158
+ L = positions.shape[0]
159
+ delta = int(rng.integers(self.delta_min, self.delta_max + 1))
160
+ t_in_max = L - delta - 1
161
+ t_in_min = self.past_min
162
+ t_in = int(rng.integers(t_in_min, t_in_max + 1))
163
+
164
+ past_start = max(0, t_in - self.past_max)
165
+ past_pos = positions[past_start:t_in]
166
+ target_pos = positions[t_in:t_in + delta]
167
+
168
+ past_features = compute_features(past_pos)
169
+ T_past = past_features.shape[0]
170
+ feat_pad = np.full((self.past_max, 9), PAD_VALUE, dtype=np.float32)
171
+ feat_pad[:T_past] = past_features
172
+ tgt_pad = np.zeros((self.delta_max, 3), dtype=np.float32)
173
+ tgt_pad[:delta] = target_pos
174
+ return {
175
+ "past_features": torch.from_numpy(feat_pad),
176
+ "past_length": torch.tensor(T_past, dtype=torch.long),
177
+ "target_pos": torch.from_numpy(tgt_pad),
178
+ "delta": torch.tensor(delta, dtype=torch.long),
179
+ "label": torch.tensor(int(self.labels[traj_idx]), dtype=torch.long),
180
+ }
181
+
182
+
183
+ # ============================================================================
184
+ # MODEL
185
+ # ============================================================================
186
+
187
+ def sinusoidal_embedding(values, dim):
188
+ half = dim // 2
189
+ device = values.device
190
+ freqs = torch.exp(-math.log(10000.0)
191
+ * torch.arange(half, device=device) / half)
192
+ angles = values.float().unsqueeze(-1) * freqs
193
+ emb = torch.cat([torch.sin(angles), torch.cos(angles)], dim=-1)
194
+ if dim % 2 == 1:
195
+ emb = F.pad(emb, (0, 1))
196
+ return emb
197
+
198
+
199
+ class LearnablePosEnc(nn.Module):
200
+ def __init__(self, max_len, d_model):
201
+ super().__init__()
202
+ self.pe = nn.Parameter(torch.randn(1, max_len, d_model) * 0.02)
203
+ def forward(self, x):
204
+ return x + self.pe[:, :x.size(1)]
205
+
206
+
207
+ class PatchTokenizer(nn.Module):
208
+ def __init__(self, in_channels=9, d_model=256, patch_size=8, max_patches=64):
209
+ super().__init__()
210
+ self.patch_size = patch_size
211
+ self.d_model = d_model
212
+ self.embed = nn.Sequential(
213
+ nn.Conv1d(in_channels, d_model // 2, 5, padding=2),
214
+ nn.GELU(),
215
+ nn.Conv1d(d_model // 2, d_model, 3, padding=1),
216
+ nn.GELU(),
217
+ )
218
+ self.pos_enc = LearnablePosEnc(max_patches, d_model)
219
+ self.norm = nn.LayerNorm(d_model)
220
+
221
+ def forward(self, features, lengths):
222
+ B, T, C = features.shape
223
+ h = self.embed(features.transpose(1, 2))
224
+ N = max(1, T // self.patch_size)
225
+ h = h[:, :, :N * self.patch_size]
226
+ h = h.reshape(B, self.d_model, N, self.patch_size).mean(-1)
227
+ h = h.transpose(1, 2)
228
+ h = self.norm(self.pos_enc(h))
229
+ patch_lengths = (lengths.float() / self.patch_size).clamp(min=1).long()
230
+ patch_lengths = patch_lengths.clamp(max=N)
231
+ return h, patch_lengths
232
+
233
+
234
+ class CausalEncoder(nn.Module):
235
+ def __init__(self, d_model=256, n_heads=8, n_layers=4, d_ff=1024, dropout=0.1):
236
+ super().__init__()
237
+ layer = nn.TransformerEncoderLayer(
238
+ d_model=d_model, nhead=n_heads, dim_feedforward=d_ff,
239
+ dropout=dropout, activation="gelu", batch_first=True,
240
+ norm_first=True,
241
+ )
242
+ self.tf = nn.TransformerEncoder(layer, num_layers=n_layers)
243
+ self.norm = nn.LayerNorm(d_model)
244
+
245
+ def forward(self, x, key_padding_mask):
246
+ N = x.size(1)
247
+ causal_mask = torch.triu(
248
+ torch.ones(N, N, dtype=torch.bool, device=x.device), diagonal=1
249
+ )
250
+ return self.norm(
251
+ self.tf(x, mask=causal_mask, src_key_padding_mask=key_padding_mask)
252
+ )
253
+
254
+
255
+ def last_valid_token(encoded, patch_lengths):
256
+ B, N, D = encoded.shape
257
+ idx = (patch_lengths - 1).clamp(min=0).view(B, 1, 1).expand(-1, 1, D)
258
+ return encoded.gather(1, idx).squeeze(1)
259
+
260
+
261
+ class DeltaEmbedding(nn.Module):
262
+ def __init__(self, d_model=256, d_freq=64):
263
+ super().__init__()
264
+ self.d_freq = d_freq
265
+ self.proj = nn.Sequential(
266
+ nn.Linear(d_freq * 2, d_model),
267
+ nn.GELU(),
268
+ nn.Linear(d_model, d_model),
269
+ )
270
+ def forward(self, delta, t_past):
271
+ d_emb = sinusoidal_embedding(delta.float(), self.d_freq)
272
+ rel = delta.float() / t_past.float().clamp(min=1.0)
273
+ rel_emb = sinusoidal_embedding(rel * 100.0, self.d_freq)
274
+ return self.proj(torch.cat([d_emb, rel_emb], dim=-1))
275
+
276
+
277
+ class GaussianHead(nn.Module):
278
+ def __init__(self, d_model=256, d_hidden=256):
279
+ super().__init__()
280
+ self.net = nn.Sequential(
281
+ nn.Linear(d_model, d_hidden), nn.GELU(),
282
+ nn.Linear(d_hidden, d_hidden), nn.GELU(),
283
+ )
284
+ self.mu_head = nn.Linear(d_hidden, 3)
285
+ self.log_sigma_head = nn.Linear(d_hidden, 3)
286
+ self.rho_head = nn.Linear(d_hidden, 1)
287
+
288
+ def forward(self, h):
289
+ z = self.net(h)
290
+ delta_mu = self.mu_head(z)
291
+ log_sigma = self.log_sigma_head(z).clamp(min=-7.0, max=2.0)
292
+ rho = torch.tanh(self.rho_head(z)).squeeze(-1) * 0.99
293
+ return delta_mu, log_sigma, rho
294
+
295
+
296
+ def gaussian_nll_xyz(true_delta, mu, log_sigma, rho):
297
+ sx = log_sigma[:, 0].exp()
298
+ sy = log_sigma[:, 1].exp()
299
+ sz = log_sigma[:, 2].exp()
300
+ dx = true_delta[:, 0] - mu[:, 0]
301
+ dy = true_delta[:, 1] - mu[:, 1]
302
+ dz = true_delta[:, 2] - mu[:, 2]
303
+ omr2 = (1.0 - rho * rho).clamp(min=1e-6)
304
+ z2 = (((dx / sx) ** 2)
305
+ - 2.0 * rho * (dx / sx) * (dy / sy)
306
+ + ((dy / sy) ** 2)) / omr2
307
+ log_det = 2.0 * (log_sigma[:, 0] + log_sigma[:, 1]) + torch.log(omr2)
308
+ nll_xy = 0.5 * (z2 + log_det + 2.0 * math.log(2.0 * math.pi))
309
+ nll_z = 0.5 * ((dz / sz) ** 2 + 2.0 * log_sigma[:, 2]
310
+ + math.log(2.0 * math.pi))
311
+ return nll_xy + nll_z
312
+
313
+
314
+ class FuturePredictor(nn.Module):
315
+ def __init__(self, d_model=256, pred_dim=128, dropout=0.1):
316
+ super().__init__()
317
+ self.proj_in = nn.Linear(d_model * 2, pred_dim)
318
+ layer = nn.TransformerEncoderLayer(
319
+ d_model=pred_dim, nhead=4, dim_feedforward=pred_dim * 2,
320
+ dropout=dropout, activation="gelu", batch_first=True, norm_first=True,
321
+ )
322
+ self.tf = nn.TransformerEncoder(layer, num_layers=2)
323
+ self.proj_out = nn.Linear(pred_dim, d_model)
324
+ self.norm = nn.LayerNorm(d_model)
325
+
326
+ def forward(self, z_in, delta_emb):
327
+ h = self.proj_in(torch.cat([z_in, delta_emb], dim=-1)).unsqueeze(1)
328
+ h = self.tf(h)
329
+ return self.norm(self.proj_out(h.squeeze(1)))
330
+
331
+
332
+ class FlightJEPAv2(nn.Module):
333
+ def __init__(self, cfg):
334
+ super().__init__()
335
+ self.cfg = cfg
336
+ d = cfg.get("d_model", 256)
337
+ h_ = cfg.get("n_heads", 8)
338
+ n_l = cfg.get("n_layers", 4)
339
+ d_ff = cfg.get("d_ff", 1024)
340
+ dr = cfg.get("dropout", 0.1)
341
+ ps = cfg.get("patch_size", 8)
342
+ past_max = cfg.get("past_max", 256)
343
+ max_patches = past_max // ps
344
+ self.lambda_jepa = cfg.get("lambda_jepa", 0.0)
345
+ self.ema_decay = cfg.get("ema_decay", 0.998)
346
+
347
+ self.tokenizer = PatchTokenizer(9, d, ps, max_patches)
348
+ self.encoder = CausalEncoder(d, h_, n_l, d_ff, dr)
349
+ self.delta_emb = DeltaEmbedding(d, 64)
350
+ self.head = GaussianHead(d, d)
351
+ self.fuse_in = nn.Sequential(
352
+ nn.Linear(d * 2, d), nn.GELU(),
353
+ nn.Linear(d, d),
354
+ )
355
+ self.step_cell = nn.GRUCell(input_size=3, hidden_size=d)
356
+
357
+ self.target_tokenizer = copy.deepcopy(self.tokenizer)
358
+ self.target_encoder = copy.deepcopy(self.encoder)
359
+ for p in self.target_tokenizer.parameters():
360
+ p.requires_grad = False
361
+ for p in self.target_encoder.parameters():
362
+ p.requires_grad = False
363
+ self.predictor = FuturePredictor(d, d // 2, dr)
364
+
365
+ @torch.no_grad()
366
+ def update_ema(self):
367
+ m = self.ema_decay
368
+ for online, target in [(self.tokenizer, self.target_tokenizer),
369
+ (self.encoder, self.target_encoder)]:
370
+ for po, pt in zip(online.parameters(), target.parameters()):
371
+ pt.data.mul_(m).add_(po.data, alpha=1.0 - m)
372
+
373
+ def encode_past(self, past_features, past_length):
374
+ patches, patch_lens = self.tokenizer(past_features, past_length)
375
+ N = patches.size(1)
376
+ pad_mask = (torch.arange(N, device=patches.device).unsqueeze(0)
377
+ >= patch_lens.unsqueeze(1))
378
+ encoded = self.encoder(patches, key_padding_mask=pad_mask)
379
+ z_in = last_valid_token(encoded, patch_lens)
380
+ return z_in, encoded, patch_lens
381
+
382
+ @torch.no_grad()
383
+ def encode_future_target(self, target_features, target_length):
384
+ patches, patch_lens = self.target_tokenizer(target_features, target_length)
385
+ N = patches.size(1)
386
+ pad_mask = (torch.arange(N, device=patches.device).unsqueeze(0)
387
+ >= patch_lens.unsqueeze(1))
388
+ encoded = self.target_encoder(patches, key_padding_mask=pad_mask)
389
+ return last_valid_token(encoded, patch_lens)
390
+
391
+ def forward(self, past_features, past_length, target_pos, delta, last_pos):
392
+ B = past_features.size(0)
393
+ device = past_features.device
394
+ delta_max = target_pos.size(1)
395
+
396
+ z_in, _, _ = self.encode_past(past_features, past_length)
397
+ delta_e = self.delta_emb(delta, past_length)
398
+ h = self.fuse_in(torch.cat([z_in, delta_e], dim=-1))
399
+
400
+ prev_pos = last_pos
401
+ nll_total = torch.zeros(B, device=device)
402
+ valid_steps = torch.zeros(B, device=device)
403
+ ade_total = torch.zeros(B, device=device)
404
+
405
+ for t in range(delta_max):
406
+ delta_mu, log_sigma, rho = self.head(h)
407
+ true_pos_t = target_pos[:, t]
408
+ true_delta = true_pos_t - prev_pos
409
+ nll = gaussian_nll_xyz(true_delta, delta_mu, log_sigma, rho)
410
+ mask = (t < delta).float()
411
+ nll_total = nll_total + nll * mask
412
+ ade_total = (ade_total
413
+ + (true_delta - delta_mu).pow(2).sum(-1).sqrt() * mask)
414
+ valid_steps = valid_steps + mask
415
+ h = self.step_cell(true_delta, h)
416
+ prev_pos = true_pos_t
417
+
418
+ nll_loss = (nll_total / valid_steps.clamp(min=1.0)).mean()
419
+ ade_train = (ade_total / valid_steps.clamp(min=1.0)).mean().detach()
420
+
421
+ losses = {"nll": nll_loss, "ade_train": ade_train, "total": nll_loss}
422
+
423
+ if self.lambda_jepa > 0.0:
424
+ tgt_feat = torch.zeros(B, delta_max, 9, device=device)
425
+ tgt_feat[..., :3] = target_pos
426
+ z_target = self.encode_future_target(tgt_feat, delta)
427
+ z_pred = self.predictor(z_in, delta_e)
428
+ jepa_loss = F.l1_loss(z_pred, z_target.detach())
429
+ losses["jepa"] = jepa_loss
430
+ losses["total"] = nll_loss + self.lambda_jepa * jepa_loss
431
+
432
+ return losses
433
+
434
+ @torch.no_grad()
435
+ def rollout(self, past_features, past_length, delta, last_pos, delta_max):
436
+ B = past_features.size(0)
437
+ device = past_features.device
438
+ z_in, _, _ = self.encode_past(past_features, past_length)
439
+ delta_e = self.delta_emb(delta, past_length)
440
+ h = self.fuse_in(torch.cat([z_in, delta_e], dim=-1))
441
+ prev_pos = last_pos
442
+ mu_pos = torch.zeros(B, delta_max, 3, device=device)
443
+ sigma = torch.zeros(B, delta_max, 3, device=device)
444
+ rho_out = torch.zeros(B, delta_max, device=device)
445
+ for t in range(delta_max):
446
+ delta_mu, log_sigma, rho = self.head(h)
447
+ cur_pos = prev_pos + delta_mu
448
+ mu_pos[:, t] = cur_pos
449
+ sigma[:, t] = log_sigma.exp()
450
+ rho_out[:, t] = rho
451
+ h = self.step_cell(delta_mu, h)
452
+ prev_pos = cur_pos
453
+ return mu_pos, sigma, rho_out
454
+
455
+
456
+ # ============================================================================
457
+ # TRAIN + SCORE
458
+ # ============================================================================
459
+
460
+ RMAX_KM = 120.0
461
+ DELTA_BUCKETS = [(30, 60), (60, 90), (90, 120)]
462
+ EXTRAP_DELTAS = [180, 300]
463
+ THRESH_M = [500.0, 1000.0, 2000.0]
464
+
465
+
466
+ def get_last_pos(past_features, past_length):
467
+ B = past_features.size(0)
468
+ idx = (past_length - 1).clamp(min=0)
469
+ return past_features[torch.arange(B, device=past_features.device), idx, :3]
470
+
471
+
472
+ def train_one_epoch(model, loader, optimizer, device, grad_clip=1.0):
473
+ model.train()
474
+ sums = {"nll": 0.0, "ade": 0.0, "jepa": 0.0, "total": 0.0, "n": 0}
475
+ for batch in loader:
476
+ past_f = batch["past_features"].to(device)
477
+ past_l = batch["past_length"].to(device)
478
+ target = batch["target_pos"].to(device)
479
+ delta = batch["delta"].to(device)
480
+ last_pos = get_last_pos(past_f, past_l)
481
+ losses = model(past_f, past_l, target, delta, last_pos)
482
+ optimizer.zero_grad()
483
+ losses["total"].backward()
484
+ torch.nn.utils.clip_grad_norm_(model.parameters(), grad_clip)
485
+ optimizer.step()
486
+ if model.lambda_jepa > 0.0:
487
+ model.update_ema()
488
+ bs = past_f.size(0)
489
+ sums["nll"] += losses["nll"].item() * bs
490
+ sums["ade"] += losses["ade_train"].item() * bs
491
+ if "jepa" in losses:
492
+ sums["jepa"] += losses["jepa"].item() * bs
493
+ sums["total"] += losses["total"].item() * bs
494
+ sums["n"] += bs
495
+ n = max(sums["n"], 1)
496
+ return {k: v / n for k, v in sums.items() if k != "n"} | {
497
+ "ade_train": sums["ade"] / n
498
+ }
499
+
500
+
501
+ @torch.no_grad()
502
+ def score_loader(model, loader, device, extrap_delta=None):
503
+ model.train(False)
504
+ delta_max_dataset = loader.dataset.delta_max
505
+ per_sample = []
506
+ for batch in loader:
507
+ past_f = batch["past_features"].to(device)
508
+ past_l = batch["past_length"].to(device)
509
+ target = batch["target_pos"].to(device)
510
+ delta = batch["delta"].to(device)
511
+ last_pos = get_last_pos(past_f, past_l)
512
+ if extrap_delta is not None:
513
+ forced = torch.full_like(delta, extrap_delta)
514
+ roll_len = extrap_delta
515
+ else:
516
+ forced = delta
517
+ roll_len = int(delta.max().item())
518
+ if roll_len > delta_max_dataset:
519
+ continue
520
+ mu_pos, sigma, rho = model.rollout(past_f, past_l, forced, last_pos, roll_len)
521
+ active_len = torch.minimum(forced, delta).clamp(min=1)
522
+ for i in range(past_f.size(0)):
523
+ L = int(active_len[i].item())
524
+ per_sample.append({
525
+ "mu": mu_pos[i, :L].cpu().numpy(),
526
+ "sigma": sigma[i, :L].cpu().numpy(),
527
+ "rho": rho[i, :L].cpu().numpy(),
528
+ "target": target[i, :L].cpu().numpy(),
529
+ "delta_orig": int(delta[i].item()),
530
+ })
531
+
532
+ if not per_sample:
533
+ return {}
534
+ ades, fdes = [], []
535
+ in_circle = {t: [] for t in THRESH_M}
536
+ nlls, coverage95, delta_orig = [], [], []
537
+ for s in per_sample:
538
+ diff = s["target"] - s["mu"]
539
+ per_step_l2 = np.linalg.norm(diff, axis=1) * RMAX_KM * 1000.0
540
+ ades.append(per_step_l2.mean())
541
+ fdes.append(per_step_l2[-1])
542
+ for t in THRESH_M:
543
+ in_circle[t].append(per_step_l2[-1] <= t)
544
+ sx = max(s["sigma"][-1, 0], 1e-9)
545
+ sy = max(s["sigma"][-1, 1], 1e-9)
546
+ sz = max(s["sigma"][-1, 2], 1e-9)
547
+ rho_xy = s["rho"][-1]
548
+ dx = diff[-1, 0]; dy = diff[-1, 1]; dz = diff[-1, 2]
549
+ omr2 = max(1.0 - rho_xy * rho_xy, 1e-6)
550
+ z2 = ((dx / sx) ** 2 - 2 * rho_xy * (dx / sx) * (dy / sy)
551
+ + (dy / sy) ** 2) / omr2
552
+ coverage95.append(z2 <= 5.991)
553
+ log_det = 2 * (math.log(sx) + math.log(sy)) + math.log(omr2)
554
+ nll_xy = 0.5 * (z2 + log_det + 2 * math.log(2 * math.pi))
555
+ nll_z = 0.5 * ((dz / sz) ** 2 + 2 * math.log(sz) + math.log(2 * math.pi))
556
+ nlls.append(nll_xy + nll_z)
557
+ delta_orig.append(s["delta_orig"])
558
+ ades = np.array(ades); fdes = np.array(fdes)
559
+ nlls = np.array(nlls); coverage95 = np.array(coverage95, dtype=float)
560
+ delta_orig = np.array(delta_orig)
561
+ out = {
562
+ "ade_m": float(ades.mean()),
563
+ "fde_m": float(fdes.mean()),
564
+ "fde_median_m": float(np.median(fdes)),
565
+ "nll_xy_z": float(nlls.mean()),
566
+ "coverage_95": float(coverage95.mean()),
567
+ "n": len(ades),
568
+ }
569
+ for t in THRESH_M:
570
+ out[f"miss_rate_{int(t)}m"] = float(1.0 - np.mean(in_circle[t]))
571
+ if extrap_delta is None:
572
+ per_bucket = {}
573
+ for lo, hi in DELTA_BUCKETS:
574
+ mask = (delta_orig >= lo) & (delta_orig <= hi)
575
+ if mask.sum() == 0:
576
+ continue
577
+ per_bucket[f"delta_{lo}_{hi}"] = {
578
+ "ade_m": float(ades[mask].mean()),
579
+ "fde_m": float(fdes[mask].mean()),
580
+ "coverage_95": float(coverage95[mask].mean()),
581
+ "n": int(mask.sum()),
582
+ }
583
+ out["per_bucket"] = per_bucket
584
+ return out
585
+
586
+
587
+ def main():
588
+ p = argparse.ArgumentParser()
589
+ p.add_argument("--airport", default="RKSIa")
590
+ p.add_argument("--data-dir", default="data")
591
+ p.add_argument("--tag", default="run")
592
+ p.add_argument("--out-dir", default="runs")
593
+ p.add_argument("--epochs", type=int, default=30)
594
+ p.add_argument("--batch-size", type=int, default=64)
595
+ p.add_argument("--lr", type=float, default=1e-4)
596
+ p.add_argument("--weight-decay", type=float, default=1e-4)
597
+ p.add_argument("--past-max", type=int, default=256)
598
+ p.add_argument("--past-min", type=int, default=60)
599
+ p.add_argument("--delta-min", type=int, default=30)
600
+ p.add_argument("--delta-max", type=int, default=120)
601
+ p.add_argument("--extrap-delta-max", type=int, default=300)
602
+ p.add_argument("--epoch-multiplier", type=int, default=4)
603
+ p.add_argument("--lambda-jepa", type=float, default=0.0)
604
+ p.add_argument("--ema-decay", type=float, default=0.998)
605
+ p.add_argument("--d-model", type=int, default=256)
606
+ p.add_argument("--n-layers", type=int, default=4)
607
+ p.add_argument("--n-heads", type=int, default=8)
608
+ p.add_argument("--patch-size", type=int, default=8)
609
+ p.add_argument("--seed", type=int, default=0)
610
+ p.add_argument("--num-workers", type=int, default=2)
611
+ p.add_argument("--push-to-hub", action="store_true")
612
+ p.add_argument("--hub-model-id", default=None)
613
+ p.add_argument("--trackio-name", default=None)
614
+ args = p.parse_args()
615
+
616
+ torch.manual_seed(args.seed)
617
+ np.random.seed(args.seed)
618
+
619
+ device = "cuda" if torch.cuda.is_available() else "cpu"
620
+ print(f"[v2] device={device} tag={args.tag} lambda_jepa={args.lambda_jepa}")
621
+
622
+ if HAS_TRACKIO and args.trackio_name:
623
+ trackio.init(project="flight-jepa-v2", name=args.trackio_name,
624
+ config=vars(args))
625
+
626
+ train_ds = BlindspotDataset(
627
+ airport=args.airport, mode="TRAIN", data_dir=args.data_dir,
628
+ past_max=args.past_max, past_min=args.past_min,
629
+ delta_min=args.delta_min, delta_max=args.delta_max,
630
+ seed=args.seed, epoch_multiplier=args.epoch_multiplier,
631
+ )
632
+ test_ds = BlindspotDataset(
633
+ airport=args.airport, mode="TEST", data_dir=args.data_dir,
634
+ past_max=args.past_max, past_min=args.past_min,
635
+ delta_min=args.delta_min, delta_max=args.delta_max,
636
+ seed=args.seed + 1, epoch_multiplier=1,
637
+ )
638
+ extrap_ds = BlindspotDataset(
639
+ airport=args.airport, mode="TEST", data_dir=args.data_dir,
640
+ past_max=args.past_max, past_min=args.past_min,
641
+ delta_min=args.delta_min, delta_max=args.extrap_delta_max,
642
+ seed=args.seed + 99, epoch_multiplier=1,
643
+ )
644
+
645
+ train_dl = DataLoader(train_ds, batch_size=args.batch_size, shuffle=True,
646
+ num_workers=args.num_workers, pin_memory=True,
647
+ drop_last=True)
648
+ test_dl = DataLoader(test_ds, batch_size=args.batch_size, shuffle=False,
649
+ num_workers=args.num_workers, pin_memory=True)
650
+ extrap_dl = DataLoader(extrap_ds, batch_size=args.batch_size, shuffle=False,
651
+ num_workers=args.num_workers, pin_memory=True)
652
+
653
+ cfg = {
654
+ "d_model": args.d_model, "n_heads": args.n_heads,
655
+ "n_layers": args.n_layers, "d_ff": args.d_model * 4,
656
+ "dropout": 0.1, "patch_size": args.patch_size,
657
+ "past_max": args.past_max, "lambda_jepa": args.lambda_jepa,
658
+ "ema_decay": args.ema_decay,
659
+ }
660
+ model = FlightJEPAv2(cfg).to(device)
661
+ n_params = sum(p.numel() for p in model.parameters())
662
+ print(f"[v2] params={n_params/1e6:.2f}M")
663
+
664
+ optimizer = torch.optim.AdamW(model.parameters(), lr=args.lr,
665
+ weight_decay=args.weight_decay)
666
+ scheduler = torch.optim.lr_scheduler.CosineAnnealingLR(optimizer, T_max=args.epochs)
667
+
668
+ os.makedirs(args.out_dir, exist_ok=True)
669
+ history = []
670
+ best_fde = float("inf")
671
+ best_state = None
672
+
673
+ for epoch in range(args.epochs):
674
+ t0 = time.time()
675
+ train_stats = train_one_epoch(model, train_dl, optimizer, device)
676
+ scheduler.step()
677
+
678
+ score_stats = None
679
+ if (epoch + 1) % 5 == 0 or epoch == args.epochs - 1:
680
+ score_stats = score_loader(model, test_dl, device)
681
+ if score_stats and score_stats["fde_m"] < best_fde:
682
+ best_fde = score_stats["fde_m"]
683
+ best_state = {k: v.detach().cpu().clone()
684
+ for k, v in model.state_dict().items()}
685
+
686
+ elapsed = time.time() - t0
687
+ log = {
688
+ "epoch": epoch + 1, "elapsed_s": elapsed,
689
+ "lr": optimizer.param_groups[0]["lr"],
690
+ "train": train_stats, "score": score_stats,
691
+ }
692
+ history.append(log)
693
+ msg = (f"[v2] ep {epoch+1:03d} | loss={train_stats['total']:.4f} "
694
+ f"nll={train_stats['nll']:.4f} ade_t={train_stats['ade_train']:.4f} "
695
+ f"jepa={train_stats['jepa']:.4f}")
696
+ if score_stats:
697
+ msg += f" | fde={score_stats['fde_m']:.0f}m ade={score_stats['ade_m']:.0f}m"
698
+ msg += f" | {elapsed:.0f}s"
699
+ print(msg, flush=True)
700
+
701
+ if HAS_TRACKIO and args.trackio_name:
702
+ tlog = {f"train/{k}": v for k, v in train_stats.items()}
703
+ if score_stats:
704
+ tlog.update({f"test/{k}": v for k, v in score_stats.items()
705
+ if isinstance(v, (int, float))})
706
+ trackio.log(tlog, step=epoch + 1)
707
+
708
+ final = {"in_distribution": score_loader(model, test_dl, device)}
709
+ for d in EXTRAP_DELTAS:
710
+ final[f"extrap_delta_{d}"] = score_loader(model, extrap_dl, device, extrap_delta=d)
711
+
712
+ if best_state is not None:
713
+ model.load_state_dict(best_state)
714
+
715
+ out_path = os.path.join(args.out_dir, f"{args.tag}.pt")
716
+ torch.save({
717
+ "state_dict": model.state_dict(),
718
+ "config": cfg, "args": vars(args),
719
+ "history": history, "final": final,
720
+ "best_fde_m": best_fde,
721
+ }, out_path)
722
+ print(f"[v2] saved {out_path}")
723
+
724
+ summary_path = os.path.join(args.out_dir, f"{args.tag}_summary.json")
725
+ with open(summary_path, "w") as f:
726
+ json.dump({
727
+ "tag": args.tag, "lambda_jepa": args.lambda_jepa,
728
+ "n_params": n_params, "best_fde_m": best_fde,
729
+ "final": final, "args": vars(args),
730
+ }, f, indent=2, default=float)
731
+ print(f"[v2] summary -> {summary_path}", flush=True)
732
+
733
+ if args.push_to_hub and args.hub_model_id:
734
+ try:
735
+ from huggingface_hub import HfApi
736
+ api = HfApi()
737
+ api.create_repo(args.hub_model_id, exist_ok=True)
738
+ for path, fname in [(out_path, f"{args.tag}.pt"),
739
+ (summary_path, f"{args.tag}_summary.json")]:
740
+ api.upload_file(path_or_fileobj=path, path_in_repo=fname,
741
+ repo_id=args.hub_model_id)
742
+ print(f"[v2] uploaded to {args.hub_model_id}")
743
+ except Exception as e:
744
+ print(f"[v2] hub upload failed: {e}")
745
+
746
+ if HAS_TRACKIO and args.trackio_name:
747
+ trackio.finish()
748
+
749
+
750
+ if __name__ == "__main__":
751
+ main()