Merge dev: HPC/parallelism fixes, PSF photometry, and tpf-dependency removal - #46
Merged
Conversation
…ly after download
multiprocessing.cpu_count() returns the total node CPU count, not the cores allocated by SLURM. With n_jobs=-1 joblib would spawn workers for every node CPU, causing severe oversubscription on jobs with a smaller allocation and making parallel steps appear serial. _available_cores() now checks SLURM_CPUS_PER_TASK first, then os.sched_getaffinity(0) (cgroup-aware on Linux), before falling back to cpu_count(). Also removes a duplicate self.num_cores assignment that was overwriting the resolved value.
helpers.py: - fix_background_anomalies: replace labeled==i component loops with np.bincount/np.isin vectorisation in sep_validate, bad_bkg_mask, and _fit_residual; precompute distance array per SEP object to avoid redundant sqrt calls in annulus search - grad_clip_fill_bkg: replace five Python component loops (size accumulation, size filtering, overlap ratio calculation, ratio filtering, small-region removal) with np.bincount + np.isin - blend_dynamic_background: extract _blend_frame helper and dispatch all T frames via Parallel(n_jobs=n_jobs); add n_jobs parameter - parallel_strap_fit: extract _strap_fit_col helper and dispatch strap columns via Parallel(n_jobs=n_jobs); add n_jobs parameter - regional_stats_mask: extract _clip_region helper and dispatch regions via Parallel(n_jobs=n_jobs); add n_jobs parameter - Add _shift_one and _shift_ref_one module-level helpers for picklable per-frame shift operations tessreduce.py: - shift_images: parallelise both median and normal branches via Parallel(n_jobs=self.num_cores) using the new shift helpers - _bkg_median: replace list comprehension with np.nanmedian(axis=(1,2)) - blend_dynamic_background call: pass n_jobs=self.num_cores background_separator.py: - Replace nested list comprehension for yx_all coordinate array with np.mgrid + np.column_stack (two locations)
helpers.py: - _clip_region: fix image[rx,ry] → image[ry,rx] in the threshold check so it is consistent with the stats computation; previously stats were computed on image[ry,rx] (y,x) but the mask was applied to image[rx,ry] (x,y), clipping the wrong pixels - grad_clip_fill_bkg: fix component size filter to cover all labeled regions 1..n_objects; previous [:n_objects] slice excluded the last real component and included the background (label 0), allowing the final artifact to bypass size filtering tessreduce.py: - _calc_qe: replace zero-only guard with np.isfinite check so infinities produced by division by zero background are also converted to NaN before sigma-clipping; suppress the runtime warning with np.errstate
psf_photom.py: - psf_position: invert finite check so invalid (NaN/inf) ext_shift is zeroed rather than valid shifts being discarded - minimize_psf_flux: convolve a local copy of self.psf with the kernel instead of mutating self.psf in-place; previously each optimizer iteration re-convolved the already-convolved PSF, progressively broadening it through the fit tessreduce.py: - reduce: replace nanmin with nanpercentile(1) for reference baseline removal in difference imaging mode; nanmin is dominated by a single outlier negative pixel (cosmic ray, bad pixel), biasing every difference frame by that amount helpers.py: - _clip_region: fix image[rx,ry] → image[ry,rx] in threshold check so mask is applied to the same pixels as the sigma-clipped stats - grad_clip_fill_bkg: fix component size arrays to cover labels 1..n_objects (exclude background label 0, include last component)
…on HPC On OzSTAR and similar HPC systems, joblib's default loky backend (process-based) can fail silently due to process-spawning restrictions, security policies, or memory limits, causing every Parallel() call to fall back to sequential execution with no error or warning. All heavy functions (Smooth_bkg, shift, inpaint_biharmonic, griddata, Background2D, gaussian_filter, sigma_clipped_stats) are implemented in C extensions that release the GIL, so thread-based parallelism gives genuine speedup without requiring process spawning. Changed 33 Parallel() calls across tessreduce.py, helpers.py, sep_aligner.py, adaptive_background.py, lastpercent.py, background.py, and rescale_straps.py. Also fixes the standalone _fit_residual_surface call which was still using n_jobs=-1 without SLURM awareness.
prefer='threads' caused high system CPU time on OzSTAR because all threads share one address space: simultaneous malloc/free calls from numpy temporary arrays contend on the glibc allocator lock, and cross-socket memory access on NUMA hardware requires OS arbitration — both appearing as system rather than user CPU time. backend='multiprocessing' uses Linux fork() directly. Each worker gets its own address space via copy-on-write, eliminating allocator contention and NUMA cross-talk. Unlike the loky backend (which uses a forkserver that was failing silently on OzSTAR), the multiprocessing backend uses a direct fork() call that works reliably in SLURM environments.
…kling backend='multiprocessing' uses standard pickle to dispatch jobs even with fork, which cannot serialize closures. The affected functions were silently falling back to sequential execution. helpers.py: - Extract _process closure from fix_background_anomalies to module-level _fix_bkg_frame, passing per-frame 2D slices instead of the full 3D cube to avoid large-array pickling - Extract _fit_residual closure to module-level _fit_residual_bkg with explicit res_box and n_sigma parameters - Both Parallel calls now use backend='multiprocessing' explicitly tessreduce.py: - Extract _fit_frame closure from _fit_residual_surface to module-level _fit_bkg_surface_frame with explicit parameters
…orker to avoid pickle failure
…o fix pickle failure Gaussian2DKernel was passed as an argument to each joblib worker. Like SigmaClip and MedianBackground, it contains numpy dtype dispatch wrappers that cannot be pickled by standard pickle, causing silent fallback to serial execution for all 3372 frames (~50s overhead in residual surface rerun). Creating it inside _blend_frame avoids any pickling of the object.
backend='loky' is now a parameter (default) threaded through every Parallel call, since raw 'multiprocessing' breaks in Jupyter/IPython (spawned workers try to re-run the ipykernel launcher). Callers on HPC/SLURM can still pass backend='multiprocessing'. verbose is now a tiered scheme: 0 silent, 1 (default) stage announcements, 2 adds joblib's per-task Parallel output. Also fixes a TypeError in diff_lc where wcs.all_world2pix returned a length-1 array that newer NumPy refuses to cast directly to int/float.
…ibrimbore The hardcoded color-term coefficients in PS1_to_TESS_mag/SM_to_TESS_mag were stale relative to what calibrimbore currently produces for the same tess.dat bandpass and g-r color cuts (i/z terms off by ~20-28%). Re-derived via sauron(band='tess.dat', system='ps1'|'skymapper', gr_lims=[-.5,.8]).
verbose logging, and fix a NumPy 2.x crash in diff_lc scene_photom.py: shared linear scene-fit engine (bucketed PRF cache, 2D polynomial background surface instead of a flat constant, PSF-derivative columns to absorb poor difference-image subtraction residuals, catalog- informed ridge-prior crowding for neighbours, vectorized one-shot solve across an entire time series instead of a per-frame nonlinear fit). Wired in as a new tessreduce.scene_photometry() method, alongside (not replacing) psf_photometry(). field_calibration.py: two additive calibration pathways built on the same scene-fit engine -- field_calibrate_scene(catalog='ps1'|'gaia'). 'ps1' reuses the existing Tonry-locus extinction correction and PS1/SkyMapper synthetic- TESS-mag reconstruction; 'gaia' queries Gaia DR3 via astroquery (matching tessreduce's existing catalog convention) and calibrates directly against Rp with the correct Vega->AB offset. Both use magnitude-binned robust zeropoint combining and formal per-star error gating instead of field_calibrate()'s coarse sanity check; failures warn and fall back to zp_ab=20.6 rather than raising. Neither pathway touches the existing field_calibrate(). verbose is now tiered: 1 (default) prints top-level reduce() stage announcements with timing, 2 adds background()'s internal sub-steps (also timed), 3 adds joblib's per-task Parallel output and core-selection diagnostics. diff_lc: wcs.all_world2pix can return 0-d/length-1 arrays depending on input type, which newer NumPy refuses to cast directly to int/float; flatten via np.ravel(...)[0] first.
Extract flux_raw, quality, column, row, camera, ccd, mjd, wcs at load time (both the direct-tpf path and get_TESS) and close/discard the tpf object afterward so its flux cube isn't kept duplicated in memory alongside self.flux. Smooth_motion and the reference-frame quality check no longer depend on a live tpf object, which also fixes a crash in the flux-array-only init path (self.tpf is None there) since reduce()'s default alignment methods called Smooth_motion(..., self.tpf). Also fixes make_mask(useref=False) incorrectly passing ref=self.ref to Cat_mask, which made it behave like useref=True whenever self.ref was already set.
Removes tpf dependency from catalogue/mask functions (allows raw flux arrays instead of a tpf object), stores tpf-derived values as class attributes, and closes/drops the tpf object after extraction to avoid duplicating the flux cube in memory. Also fixes a useref bug in make_mask and a crash in Smooth_motion for the flux-array-only path.
…r 107 mjd_end values for sectors 56-102 were corrected (previously truncated to the start of the next sector's gap rather than its true end), and new rows were added for sectors 103-107.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
functions for picklability, worker-local SigmaClip/Background2D/Gaussian2D
kernel construction to avoid pickle failures, threading/backend tuning for
Parallel calls, and diagnostic/progress logging throughout the pipeline.
component off-by-one, QE infinity handling, and TPF cache file handle/
deletion ordering (Windows deletion fix).
tiered verbose logging, and a NumPy 2.x crash fix in diff_lc.
coefficients via calibrimbore; joblib backend/verbosity now configurable.
passed directly) so a flux cube can be supplied without a TargetPixelFile;
tpf-derived values (flux_raw, quality, column, row, camera, ccd) are now
stored as class attributes and the tpf object is closed and dropped after
extraction to avoid duplicating the flux cube in memory. Fixes a
make_mask(useref=False) bug and a Smooth_motion crash in the flux-array
path introduced along the way.