Press n or j to go to the next uncovered block, b, p or k for the previous block.
| 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 | 2x 2x 2x 2x 2x 2x 2x 33x 2x 2x 33x 2x 1x 1x 2x 29x 18x 18x 18x 29x 795x 795x 795x 18x 18x 29x 2x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 2x 1x 1x 1x 1x 1x 1x 1x 1x 1x 1x 2x 9x 9x 9x 10x 10x 10x 10x 10x 9x 2x 27x 27x 27x 27x 27x 27x 2x 1x 1x 2x 2x 2x 2x 1x 1x 1x 2x 2x 1x 1x 1x 1x 2x 1x 1x 2x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 26x 2x 36x 36x 36x 36x 36x 36x 36x 36x 2x 27x 1x 1x 1x 1x 1x 1x 27x 2x | /**
* The imported pixel table, as the rest of the reader wants it.
*
* The third {@link PixelStore}, after `CoolerStore` over HDF5 and `HicStore`
* over a `.hic` — and the easiest of the three, because the difficulty the
* interface exists to absorb is *where the pixels are*, and here they are in
* two `Int32Array`s that this object owns.
*
* ## Why this one answers questions the others refuse
*
* `HicStore.binOffsets()` throws by name: a `.hic` has no per-bin pixel index
* and inventing one would mean reading the whole file to open it. That argument
* does not apply here for the plainest of reasons — the whole file *was* read,
* and the index is a by-product of having read it. `bin1` comes out of the
* importer already ascending, so the CSR row pointer is one linear pass and is
* exact.
*
* So a pairs-backed map is a *fuller* citizen than a `.hic`: `nnz` is a real
* number rather than `null`, the stored-pixel density column works, and nothing
* downstream has a branch for the case where it does not. That is worth stating
* because it is the payoff of binning up front rather than on demand — the cost
* was paid once, at open, and every question afterwards is answered from
* memory.
*
* ## Semantics are the cooler path's, deliberately
*
* `forEachWindow` honours `plan.span` and ignores the column range, which is
* exactly what `CoolerStore` does and exactly what `scatterDense` is written
* against — the columns are discarded by the scatter. Narrowing by column here
* would be free and is still not done: two stores that emit different pixel
* sets for the same query are two stores that can disagree, and the whole point
* of the interface is that nothing above it can tell which file it is looking
* at.
*/
import { yieldToEventLoop } from '../async'
import type { Extent, PixelChunk } from '../types'
import type { ReadPlan, SpanScanOptions } from '../metacooler/matrix'
import { PIXEL_BAND_GAP } from '../metacooler/matrix'
import type { PixelStore, ScanAllOptions } from '../metacooler/store'
import { ViewerError } from '../errors'
import type { BinnedPixels } from './bin'
/** Pixels per chunk handed to a callback; the cooler path's own default. */
const CHUNK = 1 << 20
export class PairsStore implements PixelStore {
private offsets: Float64Array | null = null
constructor(private readonly pixels: BinnedPixels) {}
get nnz(): number {
return this.pixels.nnz
}
/**
* Counts are integers and are stated as such.
*
* A pairs file is a list of observations, so a bin's value is a count of
* them: there is no balancing, no expected value and no float anywhere in the
* path. Saying `'i'` is what lets the float-preservation guarantee stay off.
*/
readonly countKind = 'i' as const
scoreColumns(): string[] {
return ['count']
}
/**
* The CSR row pointer, built once.
*
* `bin1` is ascending out of the importer — the merge pass guarantees it even
* for a file that was not sorted — so this is a counting pass and nothing
* more.
*/
binOffsets(): Float64Array {
if (this.offsets !== null) return this.offsets
const { bin1, nnz, nBins } = this.pixels
const offsets = new Float64Array(nBins + 1)
let at = 0
for (let bin = 0; bin < nBins; bin++) {
while (at < nnz && bin1[at] === bin) at++
offsets[bin + 1] = at
}
this.offsets = offsets
return offsets
}
async forEachWindow(
_rows: Extent,
_cols: Extent,
plan: ReadPlan,
onChunk: (chunk: PixelChunk) => void,
options: SpanScanOptions = {},
): Promise<boolean> {
this.refuseOtherColumns(options.value)
const range = this.pixelRange(plan.span[0], plan.span[1])
return this.emit(range[0], range[1], onChunk, options)
}
async forEachBands(
bands: readonly Extent[],
onChunk: (chunk: PixelChunk) => void,
options: SpanScanOptions = {},
): Promise<boolean> {
this.refuseOtherColumns(options.value)
for (const [lo, hi] of this.mergedBands(bands)) {
if (!(await this.emit(lo, hi, onChunk, options))) return false
}
return true
}
async scanAll(options: ScanAllOptions): Promise<void> {
const width = Math.max(1, Math.trunc(options.chunkSize ?? CHUNK))
const total = this.pixels.nnz
for (let at = 0; at < total; at += width) {
const end = Math.min(at + width, total)
if (options.onChunk(this.chunk(at, end)) === false) return
options.onProgress?.(end, total)
await yieldToEventLoop()
}
}
/** Pixel index range for a half-open bin range, from the CSR pointer. */
private pixelRange(lo: number, hi: number): Extent {
const offsets = this.binOffsets()
const nBins = this.pixels.nBins
const from = offsets[Math.max(0, Math.min(nBins, lo))]
const to = offsets[Math.max(0, Math.min(nBins, hi))]
return [from, Math.max(from, to)]
}
/**
* Bin ranges as merged pixel ranges.
*
* The same shape, and the same gap constant, as `CoolerStore.pixelBandsOf`.
* The gap is about the *source's* 1 MiB reads there and about nothing here —
* but merging still matters, because each range costs a chunk boundary and a
* yield, and a selection of four hundred contigs would otherwise be four
* hundred of both.
*/
private mergedBands(extents: readonly Extent[]): Extent[] {
const bands: Extent[] = []
for (const [lo, hi] of extents) {
if (hi <= lo) continue
const [p0, p1] = this.pixelRange(lo, hi)
if (p1 > p0) bands.push([p0, p1])
}
bands.sort((a, b) => a[0] - b[0])
const merged: Extent[] = []
for (const band of bands) {
const last = merged[merged.length - 1]
if (last !== undefined && band[0] <= last[1] + PIXEL_BAND_GAP) {
last[1] = Math.max(last[1], band[1])
} else {
merged.push([band[0], band[1]])
}
}
return merged
}
/** Hand `[from, to)` over in chunks, yielding so a cancel is seen. */
private async emit(
from: number,
to: number,
onChunk: (chunk: PixelChunk) => void,
options: SpanScanOptions,
): Promise<boolean> {
const width = Math.max(1, Math.trunc(options.chunkSize ?? CHUNK))
const total = to - from
for (let at = from; at < to; at += width) {
if (options.cancelled?.() === true) return false
const end = Math.min(at + width, to)
onChunk(this.chunk(at, end))
options.onProgress?.(end - from, total)
await yieldToEventLoop()
}
return true
}
/**
* One chunk, as a `PixelChunk` — **copied**, never a view.
*
* A `subarray` would be free and is the wrong call, for a reason particular
* to this store: it is the only one whose arrays are the *session's own
* data*. `CoolerStore` hands back what it just read out of HDF5 and
* `HicStore` what it just decoded from a block, so a caller that keeps or
* neuters one of their chunks damages nothing. Here a view is a window onto
* the imported table itself, and `transfersFor` — which walks a result
* payload and transfers every buffer it finds — would *detach the whole
* pixel table* the first time a chunk reached one. The map would go blank
* with no error anywhere.
*
* So the three arrays are materialised per chunk, which is exactly what the
* other two stores pay anyway, and the failure mode stops existing rather
* than being documented as a rule for future callers to remember.
*/
private chunk(from: number, to: number): PixelChunk {
const { bin1, bin2, count } = this.pixels
return {
bin1: bin1.slice(from, to),
bin2: bin2.slice(from, to),
count: Float64Array.from(count.subarray(from, to)),
length: to - from,
}
}
/** A column that is not `count` is a question this container cannot answer. */
private refuseOtherColumns(value: string | undefined): void {
if (value === undefined || value === 'count') return
throw new ViewerError(
'unsupported',
`A map binned from a .pairs file stores one value per pixel — the number of contacts ` +
`— so there is no '${value}' column to draw. Extra per-pixel columns come from a ` +
'cooler that was written with them.',
{ value },
)
}
}
|