src/model/point-trajectory.piece.ts

0001 /**
0002  * A data-model piece for "PointTrajectory" typed data (DEX v1).
0003  *
0004  * Track parameters live in per-name columns (parallel arrays, trajectory id
0005  * equals array index). The per-trajectory point lists are ragged, so they
0006  * stay nested: `points[i]` is the point-tuple list of trajectory i, tuple
0007  * layout declared by `pointColumns`.
0008  *
0009  *   {
0010  *     "name": "CentralTrackSegments",
0011  *     "type": "PointTrajectory",
0012  *     "version": "1.0",
0013  *     "count": 2,
0014  *     "columns": { "theta": [...], "phi": [...], "pdg": [...] },
0015  *     "pointColumns": ["x", "y", "z", "t", "dx", "dy", "dz", "dt"],
0016  *     "points": [ [[x,y,z,t,...], ...], [[x,y,z,t,...], ...] ],
0017  *     "refs": { "particle": "McParticles" }        // optional
0018  *   }
0019  *
0020  * A writer declares only the columns it has: a simulation writer has pdg and
0021  * momentum, a reconstruction writer has theta/phi/qOverP. Readers bind to
0022  * declared columns, never to a fixed set.
0023  */
0024
0025 import {
0026   EntityRefLink,
0027   EventPiece,
0028   EventPieceFactory,
0029   readPieceCount,
0030 } from './event-piece';
0031
0032 /** Common PDG codes -> short particle names, for entity labels. */
0033 const PDG_NAMES: Record<number, string> = {
0034   11: 'e-', [-11]: 'e+', 13: 'mu-', [-13]: 'mu+', 22: 'gamma', [-22]: 'opt.photon',
0035   211: 'pi+', [-211]: 'pi-', 321: 'K+', [-321]: 'K-', 111: 'pi0',
0036   2212: 'p', [-2212]: 'p~', 2112: 'n', [-2112]: 'n~',
0037 };
0038
0039 /** One parameter column: numbers for measured values, strings for labels. */
0040 export type ParamColumn = (number | string | null)[];
0041
0042 export class PointTrajectoryPiece extends EventPiece {
0043   static type = 'PointTrajectory';
0044
0045   /** Number of trajectories. */
0046   count = 0;
0047
0048   /**
0049    * Parameter columns by name; each holds `count` values, trajectory id is
0050    * the index. Example keys: "theta", "phi", "q_over_p", "pdg", "charge".
0051    */
0052   columns: Record<string, ParamColumn> = {};
0053
0054   /**
0055    * The point columns define the meaning of each entry in a point tuple.
0056    * Example: ["x","y","z","t","dx","dy","dz","dt"]
0057    */
0058   pointColumns: string[] = [];
0059
0060   /**
0061    * Ragged per-trajectory point lists: `points[i]` belongs to trajectory i
0062    * and holds point tuples laid out per `pointColumns`.
0063    */
0064   points: number[][][] = [];
0065
0066   /** Reference declarations: column name -> target piece name; -1 = no reference. */
0067   refs: Record<string, string> = {};
0068
0069   constructor(name: string, origin?: any) {
0070     super(name, PointTrajectoryPiece.type, origin);
0071   }
0072
0073   /** Reads one parameter of one trajectory; null when the column is absent. */
0074   param(columnName: string, trajectoryIndex: number): number | string | null {
0075     const column = this.columns[columnName];
0076     if (column === undefined) return null;
0077     const value = column[trajectoryIndex];
0078     return value === undefined ? null : value;
0079   }
0080
0081   /**
0082    * Label from the columns painters already use for naming: type or pdg
0083    * (mapped to a particle name), charge as a fallback, momentum when the
0084    * writer declared px/py/pz.
0085    */
0086   override entityLabel(trajectoryIndex: number): string {
0087     let name = '';
0088     const type = this.param('type', trajectoryIndex);
0089     const pdg = this.param('pdg', trajectoryIndex);
0090     const charge = this.param('charge', trajectoryIndex);
0091     if (typeof type === 'string' && type) {
0092       name = type;
0093     } else if (typeof pdg === 'number') {
0094       name = PDG_NAMES[Math.floor(pdg)] ?? `pdg=${Math.floor(pdg)}`;
0095     } else if (typeof charge === 'number') {
0096       name = charge > 0 ? 'q+' : charge < 0 ? 'q-' : 'q0';
0097     }
0098
0099     let momentum = '';
0100     const px = this.param('px', trajectoryIndex);
0101     const py = this.param('py', trajectoryIndex);
0102     const pz = this.param('pz', trajectoryIndex);
0103     if (typeof px === 'number' && typeof py === 'number' && typeof pz === 'number') {
0104       momentum = ` p=${(Math.sqrt(px * px + py * py + pz * pz) / 1000).toFixed(2)}`;
0105     }
0106
0107     return `#${trajectoryIndex}${name ? ' ' + name : ''}${momentum}`;
0108   }
0109
0110   /** Resolves the `refs` declaration for one trajectory (-1 = no reference). */
0111   override entityRefs(trajectoryIndex: number): EntityRefLink[] {
0112     const links: EntityRefLink[] = [];
0113     for (const [column, targetPiece] of Object.entries(this.refs)) {
0114       const value = this.param(column, trajectoryIndex);
0115       if (typeof value === 'number' && value >= 0) {
0116         links.push({ column, targetPiece, targetIndex: Math.floor(value) });
0117       }
0118     }
0119     return links;
0120   }
0121
0122   /** calculate time range */
0123   override get timeRange(): [number, number] | null {
0124     if (this.count === 0) return null;
0125
0126     // Find the index of time column
0127     const timeIndex = this.pointColumns.indexOf("t");
0128     if (timeIndex === -1) return null;
0129
0130     let minTime: number | null = null;
0131     let maxTime: number | null = null;
0132
0133     for (const trajectoryPoints of this.points) {
0134       for (const point of trajectoryPoints) {
0135         const time = point[timeIndex];
0136         if (time == null) continue;
0137         if (minTime === null || time < minTime) minTime = time;
0138         if (maxTime === null || time > maxTime) maxTime = time;
0139       }
0140     }
0141
0142     if (minTime !== null && maxTime !== null) {
0143       return [minTime, maxTime];
0144     }
0145     return null;
0146   }
0147
0148   /**
0149    * Convert this piece to a DEX v1 piece object
0150    */
0151   override toDexObject(): any {
0152     const result: any = {
0153       name: this.name,
0154       type: this.type,
0155       version: '1.0',
0156       count: this.count,
0157       columns: { ...this.columns },
0158       pointColumns: [...this.pointColumns],
0159       points: this.points,
0160     };
0161     if (Object.keys(this.refs).length > 0) result.refs = { ...this.refs };
0162     if (this.origin !== undefined) result.origin = this.origin;
0163     return result;
0164   }
0165 }
0166
0167 /**
0168  * Factory class to deserialize from the DEX piece object to our piece instance.
0169  * Sizes are checked loudly — a malformed file throws instead of pairing
0170  * parameters with the wrong trajectory.
0171  */
0172 export class PointTrajectoryPieceFactory implements EventPieceFactory {
0173   type = PointTrajectoryPiece.type;
0174
0175   fromDexObject(obj: any): PointTrajectoryPiece {
0176     const piece = new PointTrajectoryPiece(obj["name"], obj["origin"]);
0177     piece.count = readPieceCount(obj);
0178
0179     if (Array.isArray(obj["pointColumns"])) {
0180       piece.pointColumns = [...obj["pointColumns"]];
0181     }
0182
0183     const points = obj["points"];
0184     if (!Array.isArray(points) || points.length !== piece.count) {
0185       throw new Error(`DEX piece '${piece.name}': "points" must hold count=${piece.count} ` +
0186         `entries, got ${Array.isArray(points) ? points.length : typeof points}`);
0187     }
0188
0189     const columns = obj["columns"] ?? {};
0190     for (const columnName of Object.keys(columns)) {
0191       const column = columns[columnName];
0192       if (!Array.isArray(column) || column.length !== piece.count) {
0193         throw new Error(`DEX piece '${piece.name}': column "${columnName}" must hold ` +
0194           `count=${piece.count} values (trajectory id must equal array index)`);
0195       }
0196       piece.columns[columnName] = column;
0197     }
0198
0199     if (obj["refs"] && typeof obj["refs"] === 'object') {
0200       piece.refs = { ...obj["refs"] };
0201     }
0202
0203     // Copy point lists sorted by time so partial-track painting can assume
0204     // monotonic time per trajectory
0205     const sortStart = performance.now();
0206     const timeIndex = piece.pointColumns.indexOf("t");
0207     let totalPoints = 0;
0208     piece.points = points.map((trajectoryPoints: number[][]) => {
0209       const copy = Array.isArray(trajectoryPoints) ? [...trajectoryPoints] : [];
0210       totalPoints += copy.length;
0211       if (timeIndex !== -1 && copy.length > 1) {
0212         copy.sort((a, b) => {
0213           if (a.length <= timeIndex || b.length <= timeIndex) return 0;
0214           return a[timeIndex] - b[timeIndex];
0215         });
0216       }
0217       return copy;
0218     });
0219     const sortMs = performance.now() - sortStart;
0220     if (sortMs > 10) {
0221       console.log(`[load-timing] adopt PointTrajectory '${piece.name}': copy+sort ` +
0222         `${sortMs.toFixed(1)} ms (${piece.count} trajectories, ${totalPoints} points)`);
0223     }
0224
0225     return piece;
0226   }
0227 }