Cretaceous sequence stratigraphy of the northern Baltimore Canyon Trough: Implications for basin evolution and carbon storage

Sequence stratigraphy uses stratal patterns to identify genetically related strata bounded by unconformities as sequences (e.g., Mitchum et al., 1977a, 1977b). Sequence boundaries can be objectively identified on seismic profiles by reflection terminations (onlap, downlap, erosional truncation, and toplap; e.g., Mitchum et al., 1977a, 1977b). Sequence boundaries can also be recognized in well logs and cores by the interpretation of cyclical successions of vertical changes in composition, grain size, and other physical criteria known as stacking patterns (e.g., Van Wagoner et al., 1987; Neal and Abreu, 2009; Miller et al., 2018) and in outcrops and cores by the observation of erosional surfaces, facies changes, and/or evidence for hiatuses (Van Wagoner et al., 1987; Miller et al., 2013a, 2018). Geological samples (wells, cores, and outcrops) provide important lithologic and biostratigraphic ground truth for sequence interpretations and environments of deposition. Biostratigraphy is essential in not only evaluating the fidelity of physical correlation but also evaluating temporal relationships of strata and depositional environments. An integrated application of these stratigraphic techniques allows mapping of sequences formed as the basin responds to variations in sediment supply, sea level, and basin geodynamics. We have used this comprehensive approach to determine the distribution of Cretaceous lithologies and depositional environments across a well-studied region rich in public domain data to evaluate processes molding the stratigraphic record and suitability of this region for carbon storage.

The offshore U.S. Middle Atlantic region is a passive continental margin encompassing the Baltimore Canyon Trough (BCT) basin, an elongate, asymmetric offshore basin with 2–16 km of Jurassic and younger sedimentary rocks formed following Late Triassic–Early Jurassic rifting (Grow et al., 1988) (see Supplemental Material¹ for additional background information). Accommodation in the BCT was driven by thermal subsidence (~5 km maximum), flexural loading (~2/3 of total), and compaction (Watts, 1981). The region has also been imprinted by changes in mantle dynamic topography (MDT) that caused differential uplift and subsidence (Moucha et al., 2008; Müller et al., 2008; Schmelz et al., 2021). The greatest rate of creation of accommodation occurred in the Jurassic (accounting for ~70% of all accommodation), which was associated initially with the deposition of terrestrial sandstones, coals, and evaporites and later with marine shelf limestones. A carbonate barrier reef prograded seaward from the Jurassic into the Barremian (Jansa, 1981; Mountain and Tucholke, 1985; Poag, 1985). Lower Cretaceous siliciclastic sediments overstepped and buried the carbonate reef (Mountain and Tucholke, 1985; Poag, 1985). Cretaceous formations in the BCT have been previously identified as Missisauga (Berriasian to Barremian heterolithic sand-prone strata), Logan Canyon (Middle Cretaceous sandstones and subordinate mudstones), and Dawson Canyon (Upper Cretaceous mudstones) Formations (Libby-French, 1984; Table 1). The strata become progressively more marine downdip toward the deep western North Atlantic basin (Mountain and Tucholke, 1985; Poag, 1985). Deformation of Upper Jurassic and Cretaceous sediments in the northern BCT occurred in part due to the Great Stone Dome (GSD, also known as the Schlee Dome), the largest postrift igneous structure along the U.S. Atlantic Margin (Grow, 1980).

The BCT was studied with seismic and well data collected for reconnaissance purposes between 1975 and 1981 by the U.S. Geological Survey (USGS) and hydrocarbon exploration companies. Industry data sets on the outer continental shelf offshore New Jersey include >10,000 km of seismic-reflection profiles and 32 industry and two Continental Offshore Stratigraphic Test (COST) wells (Fig. 1). These well and seismic data were interpreted and analyzed by Libby-French (1984), Poag (1985), Poag and Sevon (1989), Prather (1991), and Klitgord et al. (1994) in building toward a generalized understanding of Mesozoic sedimentation in the basin, providing general insights into the Jurassic and Cretaceous stratigraphy on the outer continental shelf (OCS) and continental slope offshore New Jersey. In addition, several generations of academic seismic grids and core holes were collected on the New Jersey shelf and slope to study the effect of sea-level change, primarily in the Cenozoic (e.g., Mountain et al., 2010). Integration of multichannel seismic profiles (MCS) and cores facilitated a well-developed understanding of changes in Oligocene to Miocene sea level and sedimentation on this margin (e.g., Miller et al., 2013a, 2013b; Hodgson et al., 2017). However, the Mesozoic strata offshore New Jersey have been relatively sparsely studied since the 1980s.

The BCT was a target for oil and gas exploration due to the presence of potentially suitable geologic structures, reservoirs, and seals. Industry exploration efforts focused on normal and associated antithetic faults, downthrown simple closures, fault closures, and traps over and around the GSD, and structural closures associated with the carbonate shelf-edge reef (Prather, 1991). Although these targets were promising, very few encountered hydrocarbons, and no economically viable shows were found (Prather, 1991). The GSD was thought to trap large amounts of hydrocarbons as it appeared to have excellent closure, and bids within blocks around the GSD reached $650 million (Amato and Giordano, 1985). Despite its potential, no hydrocarbons were encountered in seven wells drilled on structure (Prather, 1991). In general, the scarcity of hydrocarbons in wells drilled in the BCT was attributed to poor migration or a lack of thermally mature source rocks, because reservoirs, traps, and seals appear to be adequate (Prather, 1991).

The poor hydrocarbon potential for the BCT coupled with the presence of excellent reservoirs and traps make it an attractive target for geological storage of carbon dioxide, carbon capture, and storage and sequestration (CCS). In applying CCS, carbon dioxide is stored in subsurface reservoirs as a supercritical fluid, relying on initial structural trapping in a porous reservoir (Metz et al., 2005). Carbon dioxide attains a supercritical state at burial pressures that require sub-surface depths greater than ~800 m (Bachu, 2002). Since carbon injections drive off ambient pore fluids, it is not ideal to inject into natural gas or oil reservoirs unless hydro carbon recovery is targeted. Storage requires good seals that avoid leakage through overburden along faults or into existing exploration wells.

The recent release of industry MCS seismic data (Triezenberg et al., 2016) and recently digitized logs from 29 shelf wells (Schmelz et al., 2020b; this study) offer the opportunity to comprehensively assess the Cretaceous geologic history of the BCT basin. Using these legacy data, we present an integrated seismic, well-log, and biostratigraphic analysis of the NBCT from nearshore to continental slope settings, building out from an initial sequence stratigraphic framework and CCS evaluation (Miller et al., 2018; Fukai et al., 2020; Schmelz et al., 2020b). These sequences are identified in geophysical well logs using stacking patterns, guided by biostratigraphy (Figs. 2A, 2B, and 4; Fig. S1, see ^{footnote 1}), and seismic data using reflection terminations (e.g., Fig. 3B). We evaluate these data interpreting the structure, morphology, and sedimentary characteristics of each sequence and their distributions in space and time. Our goal is to understand the tectonic, paleogeographic, and long-term sea-level evolution of this passive margin basin through the Cretaceous as a crucial step toward predicting the spatial continuity of porous reservoirs bound by mudstoneprone seals, which were originally targeted for hydrocarbons (Prather, 1991) but are now targeted for CCS.

Guided by recently released, reprocessed, and cataloged seismic (Figs. 3 and 5–8), well-log (Figs. 2–4; Fig. S1, see ^{footnote 1}), biostratigraphic (Fig. S1), and limited core data, we develop a sequence stratigraphic framework for the Cretaceous fluvial-deltaic sedimentary sequences of the BCT and map these sequences throughout the northern BCT (Fig. 9). We divided the stratigraphy of the northern Baltimore Canyon Trough (NBCT) into sequences by integrating interpretations from reflection terminations observed on MCS profiles and stacking patterns observed in geophysical well logs, supplemented by biostratigraphy, lithologies described from cuttings, and limited core data (see Miller et al., 2018, for core-log integration at COST B-2).

Jurassic and Cretaceous lithologic units of the BCT were named using formational names derived from lithologically similar units on the Scotian Shelf (McIver, 1972) and were interpreted in BCT commercial wells by Libby-French (1984). Relying primarily on lithologic correlations, Libby-French (1984) described nine Upper Jurassic and Cretaceous stratigraphic units (Table 1). We have revisited this previously defined lithostratigraphic framework of the BCT (Libby-French, 1984) and have expanded upon recent sequence stratigraphic work by Miller et al. (2018) at COST B-2 and the GSD and by Schmelz et al. (2020b) in the southern BCT. We focus on Barremian to Cenomanian sandstones and mudstones in sequences that we identified independently in wells and in MCS profiles. We then tied well logs to MCS profiles to integrate facies, systems tracts, and sequence boundaries and lastly relied on biostratigraphy from well reports to constrain sequence ages and test correlations.

We began locating surfaces in the COST B-2 well with a complete log suite and conventional core within the Barremian to Cenomanian interval. Previous studies integrated limited core data and photographs with gamma-ray logs providing ground truth (see Miller et al., 2018, for core-log integration at COST B-2). Successions of fining or coarsening strata in the cores generally correspond to patterns of increasing or decreasing gamma response, respectively. We used these stacking patterns to identify important stratal surfaces such as sequence boundaries (SBs), maximum flooding surface (MFS), and transgressive surface (TS) in the COST B-2 well logs.

We extended this work to 29 of the 34 exploration and stratigraphic test wells drilled in the NBCT between 1978 and 1984 (Fig. S1, see ^{footnote 1}). We did not interpret sequences in five of the wells for several reasons including: (1) lack of gamma-ray or spontaneous potential logs; (2) lack of gamma-ray or spontaneous potential (SP) logs associated with the target sequence depths; or (3) well location on the modern slope in an environment where carbonate deposition and down-slope processes have been dominant, thus complicating the correlation of well, log, and seismic information.

Glauconite causes high gamma-log values that can be mistakenly interpreted as clay-prone mudstone (e.g., Lanci et al., 2002). Glauconite is common in the Cretaceous onshore, though offshore well reports note that it is more common in the post-Cenomanian section, is found only in a few concentrated zones in the mid-Cretaceous, and is absent in the Missisauga and older sections (e.g., Scholle, 1977). We caution that the presence of glauconite in the sections may complicate our interpretations, though the inferred fining- and coarsening-upward patterns (e.g., Fig. S1) appear regionally traceable and likely reflect grain-size trends and not compositional changes.

Age control for the sequences comes from the highest occurrence (HO) of specific species of planktonic foraminifera, calcareous nannofossils, and dinoflagellates identified in industry paleontology reports, revised using quantitative biostratigraphic methods to ensure consistency among wells (Jordan, 2019; Jordan et al., 2022). The ages of the sequences were temporally constrained based on chronostratigraphically significant biostratigraphic markers that we identified (Table 1):

1. The late Cenomanian and younger Dcx sequence is identified by planktonic foraminiferal markers that tie its oldest parts to the late Cenomanian (Rotalipora cushmani and T. greenhornensis Zone).

2. The early Cenomanian (possibly latest Albian) LC1 sequence lacks age diagnostic taxa but can be placed in a temporal framework by late Cenomanian markers in the sequence above and Albian marker in the sequence below.

3. The Albian LC2 sequence is constrained in age by nannofossil, planktonic foraminiferal, and dinocyst markers (Braarudosphaera africana, Planomalina buxtorfi, and Spinidinium vestitum).

4. The Aptian LC3 sequence is tied to the dino cyst Cyclonephelium tabulatum Zone.

5. The Barremian and older Missisauga sequence(s) are correlated to the time scale by dinocysts (Aptea anaphrissa, Pseudoceratium pelliferum, and Muderongia simplex).

The ages of highest occurrence were determined using GTS2012 (Gradstein et al., 2012). Only shortlived, common (in more than five wells) species were used in this analysis. Obtaining accurate age estimates below the Aptian (126–113 Ma) proved difficult due to a paucity of Jurassic and Lower Cretaceous biostratigraphic data in paleontology reports; therefore, most interpretations of this study focus on Aptian and younger sequences.

We objectively derived sequence boundaries and systems tracts in each well based on stacking patterns in gamma-ray log cross sections and guided by biostratigraphy (e.g., Miller et al., 2018). Using COST B-2 as an example for ground-truth to gamma response (Fig. 4), arrows or triangular columns pointing in the fining direction placed next to the gamma-ray log indicate a pattern of increasing or decreasing gamma response, interpreted as sediment fining or coarsening upward, respectively. Groupings of these stacking patterns were used to indicate surfaces based on the sequence stratigraphic model (Fig. 4; see Miller et al., 2018, for additional examples from the COST B-2 well). Sequence boundaries were also guided by bio-stratigraphic data (i.e., chosen so as not to cross biostratigraphic markers).

Sequences generally follow a repetitive regressive-transgressive-regressive pattern that defines three packages: the LST, the TST, and the HST (Van Wagoner et al., 1987). The LST is regressive, has coarsening-upward parasequences in deltaic depositional environments, and is characterized by a progradational parasequence stacking pattern that typically grades to a more aggradational stacking pattern moving up through the section (e.g., Vail, 1987; Fig. 2). We recognize LSTs based on aggradational-progradational log stacking patterns above sequence boundaries (Fig. 2). Above the LST, the TS marks a change in stacking pattern from aggradational to retrogradational that is generally well expressed in Logan Canyon sequences (Fig. 2). The TST is transgressive and identified by generally fining-upward parasequences, with an overall retrogradational stacking pattern (Fig. 2). The MFS is placed at the maximum transgression or the transition between transgression to regression at the boundary between the TST and the HST (Fig. 2). MFS is readily identifiable by a change in gamma-ray stacking patterns (from increasing to decreasing upsection) and supported by paleontologic evidence for maximum water depths (Jordan, 2019; Jordan et al., 2022). The HST is regressive, associated with aggradational to progradational to degradational parasequence sets, is placed above the MFS, and is identified by an overall coarsening-upward stacking pattern (e.g., Vail, 1987; Neal and Abreu, 2009).

Sequence boundaries were placed at the top of the HST, following Miller et al. (2018). In most wells for the LC3, LC2, and LC1 sequence boundaries, this was located at the top of a blocky sandstone that capped a series of genetically related coarsening-upward parasequences (Figs. 2A and 2B). The contact between the HST below and LST above is generally sharp and regionally traceable in wells. In our wells (Figs. 2A and 2B) and in general (e.g., Vail, 1987), Blocky sandstones that typify the top of the HST may be classified as falling stage systems tracts (FSST; Posamentier et al., 1992; Plint and Nummedal, 2000), though we did not make this distinction (see Miller et al., 2018, for discussion of FSST).

Biostratigraphy helped to guide sequence stratigraphic interpretations, particularly in downdip sections where log character becomes less clear. A lack of clear coarsening and fining-upward parasequences within the LC1 sequence of some wells (i.e., Fig. 2B) made placing the DCx sequence boundary difficult. In wells where this was the case, we used biostratigraphy, placing the DCx sequence boundary in all wells below the foraminifera Rotalipora cushmani. In general, there was good correlation between our well-log sequence boundaries and our seismic sequences boundaries (±30 m; Figs. 2A and 2B).

We used seismic-reflection terminations to identify sequences independent of well data (Fig. 3). Diagnostic criteria included erosional truncation, toplap, onlap, and downlap that highlighted major regional unconformities formed by erosion or non-deposition and bounding a depositional sequence (Mitchum et al., 1977a, 1977b). The majority of interpreted seismic data include >10,000 km of MCS reflection profiles collected by industry for exploration purposes in the 1970s and 1980s, covering ~30,000 km² of the outer continental shelf and slope with ~2 km line spacing (Fig. 1). These data were recently released by the National Archive of Marine Seismic Surveys (NAMSS; Triezenberg et al., 2016). Within the interval of interest for this study (between ~1–3 s two-way travel time, TWTT), the vertical seismic resolution decreases with depth and ranges from ~25 m to 75 m, respectively, based on checkshot calibration (see Seismic to Well Ties section).

This data set was supplemented by two higher-resolution seismic surveys (Fig. 1). The first includes ~500 km of USGS seismic lines collected between 1973 and 1978 and recently reprocessed in 2016–2017 by Absolute Imaging for the Mid-Atlantic U.S. Offshore Carbon Storage Resource Assessment Project (Gupta, 2019; https://www.osti.gov/servlets/purl/1566748). The second is a 2-D MCS academic survey collected on cruise R/V Ewing 9009 in 1990, covering the continental shelf offshore New Jersey. The Ewing 9009 sound source achieved slightly better seismic resolution (~10 m) than the earlier data sets; yet the useful acoustic penetration was limited to roughly 1 s TWTT (~1 km). The Ewing 9009 seismic profiles were useful in correlating seismic horizons landward of the GSD, although the decrease in data quality with depth made tracing seismic reflections difficult using this data set on the OCS where the strata we examined lie ~1 km below seafloor.

We were able to clearly identify seismic sequence boundaries using reflection terminations (Mitchum et al., 1977a, 1977b; e.g., red arrows in Fig. 3) and trace them as significant regional unconformities throughout the NBCT. These surfaces could be a loop correlated across 10,000 km of seismic lines, encompassing the entirety of the NBCT within the restrictions of seismic resolution. We used seismic facies derived from the geometric arrangement of reflections within sequences wherever possible to infer depositional environments and lithology within each sequence (Mitchum et al., 1977a, 1977b). This proved difficult in many cases due to the combination of thinner sequences (<300 m) and low seismic resolution (25–75 m; Figs. 3, 5–7), preventing identification of internal configurations as progradational, retrogradational, or aggradational. As a result, identification of internal surfaces and systems tracts largely relied on stacking patterns in logs tied to the seismic profiles. This limitation made it particularly difficult to identify seismic facies within terrestrial strata where reflection continuity was especially poor and stratal thickness was especially thin.

Well logs were tied to seismic profiles to integrate interpreted facies, systems tracts, sequence boundaries, and biostratigraphy observed in wells with sequence boundaries independently determined by seismic sequence analysis. Checkshot surveys collected by industry for each well were the primary means of converting from time on profiles to depth in wells. Synthetic seismograms created by integration of checkshots, sonic logs (DTs) and bulk density (RHOB) provided a way to check the well-seismic tie by comparison of the generated synthetic trace and the nearby observed seismic trace. Following this comparison, the original depth-time chart was shifted (±10 ms) to best fit the correlation between the synthetic and seismic traces. Error for the seismic-well tie increased for wells that were drilled on the GSD that encountered significant topography and presumed difference in the regional time-depth function in Cretaceous strata (e.g., Conoco 590-1), as well as for wells that were drilled at a significant distance from a seismic line (e.g., Exxon 684-1).

To generate structural contour and isopach maps in meters, seismic surfaces were converted to depth following the procedure described by Schmelz et al. (2020b). This included supplementing velocity-depth grids generated by Klitgord et al. (1994) with a higher-resolution time-depth profile for the top 1 s (~800 m) of data generated by Mountain and Monteverde (2012). The accuracy of the conversion of the seismic data (TWTT) to depth was checked against values obtained from OCS well velocity surveys at the well locations (checkshot surveys and velocity logs; Schmelz et al., 2020b).

We identified, traced, and mapped three mid-Cretaceous sequences (LC3, LC2, and LC1; Table 1) using geophysical logs at 29 of the 34 wells on the New Jersey continental shelf and slope, including seven wells on the GSD (Fig. S1, see ^{footnote 1}). In addition, the sequence boundaries constraining the Missisauga (below LC3 to LC1) and the DCx (above LC3 to LC1) composite sequences were delineated, although individual systems tracts within each sequence were not identified. As a result, six distinct seismic reflectors (from oldest to youngest: LK1, MK3, MK2, MK1, UK2, UK1) associated with the sequences identified in well logs (Table 1), were traced across the northern Baltimore Canyon Trough region (i.e., Fig. 3). Each of these six reflectors represents a seismic sequence boundary and correlates within ±25 m to sequence boundaries identified based on well logs (e.g., thick dashed lines mark seismic sequence boundaries; red lines mark well-log sequence boundaries; Figs. 2A and 2B). We describe the evolution of the Middle Cretaceous sequences from oldest to youngest.

The LK2 seismic reflector correlates to the base of the Missisauga (hereafter called “Miss”) composite sequence boundary identified in well logs and can be correlated to the onshore Wastegate Formation (Table 1). As discussed below, seismic, well-log, and cutting samples all indicate that this section, which includes both the Missisauga and Naskapi Formations of Libby-French (1984), are heterolithic sandstones and mudstones deposited in terrestrial or possibly marginal marine environments. The LK2/Miss could likely be divided into two or more sequences, but the paucity of biostratigraphic data and chaotic seismic facies of the LK2/Miss made it a challenge to subdivide. Thus, we focused on the Aptian and younger and the excellent reservoirs of the overlying Logan Canyon sequences.

The LK1 seismic reflector correlates to the basal LC3 sequence boundary identified in well logs (Table 1). LC3 is the oldest of three sequences within the Logan Canyon Formation and is Aptian based primarily on the occurrence of dinoflagellate Cyclonephelium tabulatum. The LC3 sequence boundary is placed at a change from blocky sandstones of the underlying Miss (Barremian) sequence interpreted as the Miss/LK2 HST, to an overlying coarsening-upward package (LST) of the LC3 (Aptian) sequence (Figs. 2A, 2B, and 4). This log pattern is clear in several wells (e.g., COST B-2, HOM 855, Conoco 590; Figs. 2 and 4). Placement of the basal sequence boundary becomes progressively more difficult to identify in well logs in the southwestern portion of our study area, as the sequence become finer grained in the southern wells likely due to greater distance from source (Fig. 2A). In these wells, the placement of the basal sequence boundary is primarily guided by biostratigraphy. In general, the LST of LC3 is identified in gamma logs by a thin package of coarsening-upward parasequences that thicken to the south of the GSD (Fig. 2A) and thin downdip (Fig. 2B; Miller et al., 2018). The LST sandstones are maintained downdip, until reaching the modern shelf break, at which point they transition to mudstones (i.e., COST B-3, Fig. 2B). Similarly, the sandstone content of the LC3 lowstand remains consistent along strike until the Mobil 17-2 well, at which point the systems tract is more mudstone-prone to the southwest (Fig. 2A), again suggesting greater distance from source. In well logs, several fining-upward parasequences overlying this LST are interpreted as the LC3 TST, capped by a MFS. A thick, coarsening upward HST occurs above the MFS that becomes finer grained and thins downdip (to the east; Fig. 2B) and in the southwestern portion of our study area (Fig. 2A).

Gamma-ray log facies interpretations were calibrated with core data at COST B-2 for the LC3 sequence (Miller et al., 2018). Interpretations of gamma-ray logs incorporated with observations from core data at COST B-2 show that blocky, low gamma-ray–value sandstones within the HST of LC3 thicken and increase in percent sandstone upsection and have sharp upper contacts (Fig. 4). Thus, based on the stacking patterns observed in the well logs that were calibrated to depositional environments using core data, Miller et al. (2018) concluded that the HST of the LC3 sequence is deltaic; the blocky sandstones are provisionally interpreted as delta front (mouth bars, distributary channel fills, or possibly shorefaces). The permeability and porosity values of the LC3 HST blocky sandstones are high (>50 mD; 10%–40%), while the permeability and porosity of the transgressive mudstones are low (<1 mD; <10%; Fig. 4).

The LK1 seismic reflector associated with the basal LC3 sequence boundary is identified above the middle Barremian structural unconformity produced by the GSD igneous intrusion event (Figs. 5 and 6). LK1 is continuous throughout the basin and can be loop-correlated along seismic lines in the northern and southern BCT; salt diapirism is minimal and limited to one identified dome (Figs. 7 and 10), and the seismic sequence can be generally correlated across faults discussed below (see also Schmelz et al., 2020b). Seismic sequence boundary LK1 is associated with erosional truncation below and downlap and onlap onto the surface (i.e., Fig. 3). In the north, reflector LK1 is continuous and high amplitude; farther south, reflector LK1 is difficult to trace and develops a discontinuous, wavy character with varying amplitude and discontinuous nature (Fig. 7), suggesting facies changes to more seismically incoherent terrestrial deposits.

There are few internal reflections within the LK1 seismic sequence in the north (bracketed by LK1 at the base and reflector MK3 at the top), and those that are present are low amplitude, discontinuous, and appear to onlap onto the LK1 surface (i.e., Fig. 8). The LK1 seismic sequence is thinnest in the northern BCT (on average ~200 m) and thickens slightly to the south within our data (~300 m) (Fig. 9). Schmelz et al. (2020b) reported that the LK1 seismic sequence is significantly thicker offshore Maryland, implying a predominant southern source of sediment, although, as we note, a secondary northern source is indicated by progressive fining toward the southwest, away from COST B-2.

The MK3 seismic sequence boundary correlates to the basal LC2 sequence boundary recognized in the geophysical well logs (Table 1). The LC2 sequence contains thick, porous, and permeable sandstones that correlate to the upper Logan Canyon Formation defined by Libby-French (1984) (Table 1; Miller et al., 2018). Planktonic foraminiferal, nannofossil, and dinocyst biostratigraphy indicate the age of the LC2 sequence is Albian. In the NBCT, the lower part of the LC2 sequence is composed of small, blocky sandstone bodies that thicken down dip and are interpreted as LST. The upper part of the LC2 sequence is associated with a gradual upward decrease in gamma response culminating in several thick (>10 m), blocky sands, which we assign to the HST. These blocky sandstones at the top of the HST (and the sequence) could possibly be identified as falling stage systems tract (FSST). However, we do not make that distinction because the differentiation of HST and FSST, or the identification of the basal surface of forced regression that identifies the transition from normal to forced regression, can be both difficult and imprecise (Miller et al., 2018). In between, generally fining-upward parasequences are interpreted as TST; the TSTs are generally thinner (25–50 m; Fig. 2) than the HSTs. At the top of the sequence, the LC2/LC1 boundary is represented by a sharp change from low gamma-ray values of the HST to higher gamma-ray values of the overlying shalier sequence (Figs. 2A and 2B). In the northern portion of our study area, sandstone is preserved both updip and downdip in the HST (Fig. 2B). To the south of the GSD, HST and LST sandstones of the LC2 sequence give way to mudstones (Fig. 2A) and a uniform high gamma response indicating claystone is preserved in both the LC2 and LC3 sequences (e.g., Shell 272-1 well in Fig. 2A; Schmelz et al., 2020b). This suggests a secondary source that originates north of our study area (Fig. 9).

Like the LC3 sequence, LC2 was interpreted by Miller et al. (2018) to represent deltaic deposition. Paleodepths estimated by Sikora and Olsson (1991) place the LC2 at COST B-2 in middle to outer neritic (~60–130 m) paleodepths (Fig. 4). Peak permeabilities (>1000 mD) were originally reported at COST B-2 within the blocky highstand sandstones of the LC2. These values were remeasured by Core Laboratories for an unpublished contract report (P. McLaughlin and M. KunleDare, 2018, written commun.) and confirm high permeabilities (hundreds of mD) at ~8240 ftbkb; although the very high permeabilities (>1000 mD) are likely attributable to a measurement artifact (Fig. 4). Porosities at this location range from 10%–30%. The MK3 seismic sequence boundary is associated with the basal LC2 sequence boundary defined in the geophysical well logs (Table 1) and is a high-amplitude, continuous reflection in the north (i.e., Fig. 5) that is irregular and discontinuous to the south as the section becomes shalier (i.e., Fig. 7). The surface is associated with reflection terminations that indicate erosional truncation of underlying reflections and downlapping reflections that terminate onto it. Internal reflections are limited in the northern BCT. Isochron maps (Fig. 9) show approximately uniform thickness of the LC2/LK3 sequence throughout the northern BCT (on average 200 m thick), with a slight thickening occurring in the southern part of the study area (~300 m thick).

The MK2 seismic reflector correlates with the base of the LC1 sequence (Table 1). The LC1 sequence lacks age-diagnostic taxa but is interpreted as early Cenomanian or possibly latest Albian because late Cenomanian markers are present in the DCx sequence above and Albian markers are identified in the LC2 sequence below. The lower Cenomanian LC1 sequence is less sandstone-prone than the previously described LC3 and LC2 sequences (Figs. 2A, 2B, and 4) and exhibits two distinct gamma-ray log patterns. Updip, LC1 is relatively thick (~200 m). The sequence, which is characterized by low-value gamma-ray response, contains both thick coarsening-upward sandstones that we interpret to be the HST and thinner coarsening-upward sandstones below interpreted as the LST (Fig. 2B). Only updip on the structural high of the GSD are thick HST sandstones present (Fig. S1, see ^{footnote 1}). Downdip, off the structural high of the GSD, the sandstones of the HST grade to mudstones (Fig. 2B), the sequence thins, and the gamma-ray log within the sequence has an hourglass shape. This trend is present in all wells downdip and to the southwest (Fig. 2A).

Like the HSTs of the LC2 and LC3 sequences, the sandstones of the LC1 HST are permeable and porous in wells on the GSD (Fig. 4). Off structure, the sequence is mudstone-prone, and permeability is low, providing a probable seal for the LC2 sequence. Paleodepths at COST B-2 during deposition of LC1 are deeper than those of LC2 (outer neritic, ~130–180 m versus middle to outer neritic, ~60–130 m; Sikora and Olsson, 1991), reflecting a Cenomanian relative sea-level rise.

The LC1/MK2 sequence is the thickest of the Logan Canyon sequences in the NBCT near the GSD (on average 350 m; Fig. 9), and it greatly thins downdip (<100 m). Unlike the two underlying sequences, internal reflections that appear to be clinoforms can be observed in the thickest section of LC1/MK2 (Fig. 8). Reflections around and to the south of the GSD generally downlap onto MK2 (Figs. 5, 6, and 8). The sequence depocenter lies offshore Atlantic City, New Jersey, although another, thinner depocenter lies in the south (Fig. 9).

DCx is a composite sequence consisting of several embedded sequences with the MK1 seismic reflector at its base. It is identified in well logs by a thick (up to 1 km), homogeneous, high gamma-ray, mudstone-prone package, associated with the foraminiferal species Rotalipora cushmani (93.9–95.6 Ma) near its base (i.e., Figs. 2A and 2B). In the upper part of the composite sequence, there is a thin package of laterally extensive low gamma-ray sandstone that can be correlated across all wells. Individual systems tracts within this sequence are difficult to identify because of the homogeneous nature of the log data. However, within DCx, there is a fining-upward trend leading to peak gamma-ray values in all wells across the outer continental shelf (i.e., Fig. 2A). This peak occurs above the MK1 seismic reflector and above the highest common occurrence of R. cushmani in wells. In addition, this peak in gamma-ray corresponds to an interval of maximum paleodepths at COST B-2 (greater than ~180 m at B-2; Sikora and Olsson, 1991). Foraminiferal biostratigraphy correlates this mudstone with the deepest water indicators with the Cenomanian–Turonian boundary and OAE 2 (see Discussion; Kennedy et al., 2005). The low permeability and porosity values (Fig. 4) and a very thick package of mudstone make the DCx composite sequence a good seal for the LC sequences below.

The corresponding MK1 seismic reflector associated with the base of the DCx composite sequence is a prominent seismic reflection that is regionally associated with downlap onto its surface and represents the base of a thick package (~500 m) of prograding Turonian clinoforms (Figs. 5 and 8). The MK1 sequence generally thickens downdip toward the modern shelf and slope break and is thickest in the northern part of the study area (Fig. 9), though it is thinner overlying the GSD. An anomalous ENE-WSW–striking band of thicker sediments is located in the southern part of the study area (Fig. 9) and correlates to the Delaware High (DEH), a local uplift observed in the structure map (Fig. 10).

The UK2 seismic sequence lies above the seismic MK1 sequence. Although it is not defined in the gamma-ray logs, the base is a prominent reflector in the seismic data and corresponds to the upper bounding surface of the prograding clinoforms that make up MK1 (Figs. 5 and 8). The reflection is associated with toplapping and downlapping internal reflections and generally correlates to mid-Turonian biostratigraphic markers in the wells. The sequence is bound at its top by the top Campanian seismic reflector UK1 (equivalent to UK1 of Prather, 1991, and Miller et al., 2018), which was originally defined by Greenlee and Moore (1988) and Greenlee et al. (1992) as their top Cretaceous “Red” reflector. This surface correlates with the Campanian/Maastrichtian boundary and global sea-level fall (Miller et al., 2018). UK1 is regionally identified as a strong reflection that truncates seismic reflections below (i.e., Fig. 3). Like the MK1 sequence, the UK2 sequence thickens downdip toward the slope. Within this interval, there are likely additional sequences that cannot be resolved with available seismic data.

The previously described sequences are locally influenced by the large, igneous intrusion of the GSD. While the igneous dike swarm of the GSD is not clearly imaged on seismic profiles, the deformation of sediments due to intrusion is observed on several strike and dip seismic lines. Deformed sediments above the GSD are represented by mottled, discontinuous, and incoherent seismic reflections that form a concentric pattern around the dome (Figs. 5 and 6). Jurassic and Lower Cretaceous reflections on the flanks of the intrusion are high angle and truncated by a major angular unconformity (Figs. 5 and 6; Lippert, 1983; Amato and Giordano, 1985; Jansa and Pe-Piper, 1988; Prather, 1991). This unconformity occurs in the Barremian according to biostratigraphy from the well-seismic tie and agrees with the estimated age of the lamprophyre stock (ca. 116–127 ± 5 Ma on conventional K-Ar ages) interpreted by Jansa and Pe-Piper (1988). In addition, reflections are domed over the GSD for the entire Cretaceous section (i.e., Fig. 10). This doming reflects differential compaction around the GSD (see Discussion).

Sedimentation in a basin is controlled by the balance between the rate that space becomes available (accommodation) and the rate that sediment is supplied to the basin. To assess the contribution of various processes affecting Cretaceous sedimentation in the BCT, we show (Fig. 11) the tectonic and structural features of the basin, sediment accumulation rates, the potential location of sediment dispersal routes and depocenters inferred from isopach maps (Fig. 9), and the inferred regional depositional processes as the basin evolved from terrestrial, to deltaic, to shelf deposition (Fig. 11).

Even in a passive tectonic setting such as the Mid-Atlantic continental margin, tectonic and structural processes contribute to the deposition and subsequent modification of sediments and sedimentary rocks. The effects of these processes are clearly observed in the NBCT. The postrift tectonic history of the BCT is characterized by simple passive thermal subsidence and flexural accommodation from sediment loading (Steckler and Watts, 1978; Grow and Sheridan, 1988; Kominz et al., 1998; Kominz et al., 2002). Analysis of basin subsidence at the COST B-2 well (Fig. 11) by Watts and Steckler (1979) using backstripping (progressively accounting for the effects of compaction, Airy loading, and thermal subsidence), showed the most rapid thermal subsidence of the basin postrifting coincided with high sedimentation rates in the Jurassic. This was followed by a decelerated subsidence rate in the Cretaceous at an exponential rate related to thermal subsidence (Fig. 11; Watts and Steckler, 1979). During the Late Jurassic to mid-Cretaceous, there was a 60% reduction in average sedimentation rate, coincident with shifting from thermal to loaddriven subsidence and stabilization of erosion in the highland source terrains (Poag and Sevon, 1989). Additional postrift structural features contributed to local accommodation in the NBCT. These include two main structural uplifts, the GSD, previously targeted for hydrocarbon exploration (Amato and Giordano, 1985), and the newly identified Delaware High (DEH) (Fig. 10), as well as growth (Fig. 8) and interpreted basinal faults (Fig. 10).

Intrusion by the GSD deformed and impacted the morphology of sediments deposited both pre- and post-emplacement. Intrusion of the igneous dike swarm caused alteration of Upper Jurassic and Lower Cretaceous sediments (Jansa and Pe-Piper, 1988); this alteration is observed on seismic profiles as chaotic, discontinuous seismic reflections (Figs. 5 and 6). Subsequent uplift, hypothesized to be on the order of 2 km on what was then a coastal plain (Lippert, 1983), created an area as wide as 30 km (Figs. 5 and 6). Over the course of the Barremian (ca. 126–131 Ma), exposure of sediments due to uplift caused increased erosion that eventually eroded the structure flat. This is recognized on seismic profiles by a major angular unconformity that sits above the deformed seismic reflections (GSD reflection; Figs. 5 and 6).

Although intrusion and subsequent erosion ceased sometime in the Barremian–Aptian (Jansa and Pe-Piper, 1988), the emplaced structure continued to influence sequence morphology throughout the Cretaceous. These influences included decreased accommodation on the structural high and increased accommodation due to differential compaction in the moat of the structure. As a result, mid-Cretaceous sequences were thinned on top of the Dome (Fig. 9), and strata were arched around the structure for the entire Cretaceous (Fig. 10). Differential subsidence occurred because the greater thickness of sediments on the flanks of the GSD experienced greater porosity reduction than the thinner sediments directly above the igneous intrusion. Thus, sediments in the moat compacted and subsided more readily to form a symmetrical domal structure around the intrusion. Today, differential compaction accounts for ~0.2 s (~200 m) of relief in the LC3/LK1, LC2/MK3, and LC1/MK2 sequences. This relief diminishes upsection.

A lower-relief, yet similarly extensive (~20-km-wide) structural feature in the NBCT is the Delaware High (DEH), located offshore Delaware and so named here. Unlike the GSD, the DEH is not associated with an igneous intrusion event; rather, it is bracketed by two antithetic normal faults (Fig. 10). The anomalous high is likely caused by the displacement of sediment due to these faults. There is ~0.1 s (~100 m) of relief on the Cretaceous sequences of this study (i.e., Fig. 10). Upsection, the DEH broadens but maintains relief throughout the entirety of the Cretaceous. Thickening of the LC2/MK3 sequence (Fig. 9) and an anomalous thick depocenter of the DCx/MK1 sequence appear to be related to this fault-bounded structure.

Normal growth faults, termed the Gemini fault system (Poag, 1987), parallel the trend of the outer continental shelf, offset sediment by as much as 0.05 s (~50 m; Figs. 8 and 10), and were active from at least the Late Jurassic to the late Middle Miocene (Poag, 1987). Campanian to Miocene sequences interpreted by Poag (1987) are thickest toward the modern shelf-edge and upper slope in association with this fault system. Sedimentation and loading over the outer shelf and slope likely contributed to the development of the Gemini growth faults. The depocenter of the LC3/LK1 sequence is proximal and parallels several growth faults in the NBCT (Fig. 9). We would expect thickening on the down-thrown side of the normal faults, if the faults were providing accommodation for sedimentation in the basin; however, we find the opposite. Thus, we hypothesize that movement along faults during the mid-Cretaceous had little effect on accommodation during this time. It could be coincidental that thickening occurs on the up-thrown side of the fault, or it is possible that there was a period of shortening during the mid-Cretaceous. However, we note that the DCx/MK1 and UK2 sequences both thicken downdip toward the slope (see also Poag, 1987), with depocenters near the trend of the Gemini fault system. Thus, while it is possible that faults created accommodation at this time, as we discuss later, a regional regression may have driven sediment depocenters seaward. Additional basinal normal faults (some with growth) underlie the BCT between the GSD and the COST B-3 well and parallel the trend of the continental shelf and slope margin, although they appear to have had little effect on the distribution or thickness of the sequences in this study.

Although we have shown that active faulting and intrusion tectonics impact deposition of Cretaceous sequences in the NBCT, these influences are generally local in effect. This is demonstrated by the fact that zones of faulting are largely restricted to the regions beneath the modern outermost continental shelf (Fig. 9), and intrusion is limited to the GSD (Fig. 9) and one lone salt diapir (Figs. 7 and 10). As a result, the sedimentation patterns we observe can be expected to generally correspond to variations in relative sea level. Here, we discuss sedimentation patterns for the Cretaceous sequences of this study and tie them to established records of relative (Kominz et al., 2008) and inferred global (Haq, 2014) sea-level change across the margin to evaluate this hypothesis and other potential forcing mechanisms on deposition within the NBCT. In doing so, we show that the NBCT evolved in three depositional phases that imply progressively deepening water depths through the Cretaceous: from Early Cretaceous fluvial, to mid-Cretaceous deltaic, to Late Cretaceous marine deposition (e.g., Table 2).

Early Cretaceous (Berriasian to Barremian) paleoenvironments preserved in the NBCT were primarily terrestrial. The early Cretaceous LK2/Miss sequence is characterized by serrated gamma-ray (e.g., Fig. 2), heterolithic, interbedded, calcareous sandstones and mudstones containing abundant terrestrial microfossils (Poag, 1985), and variable amplitude, discontinuous chaotic seismic facies (e.g., Fig. 8; Table 2), suggesting non-marine, terrestrial depositional environments. Coeval onshore sequences (Wastegate Formation; Fig. 11) are fluvial sandstones and interbedded varicolored mudstones deposited as floodplain deposits (Hansen, 1982; Andreasen et al., 2015; Miller, 2017). Considering the back reef setting (i.e., landward of the fringing reef of the Jurassic to Early Cretaceous) for these heterolithic, interbedded, calcareous sandstones and mudstones, they are likely fluvial to shallow marginal marine (e.g., lagoonal) depositional environments. Depocenters for the LK2/Miss sequence are located on the modern outer continental shelf and slope. This depocenter is thickest offshore Maryland (Schmelz et al., 2020b), indicating a southern source of sediment. Poag and Sevon (1989) noted a significant change in sediment provenance that took place during the Barremian, shifting from thin, northern sourced deposition to thick southern source deposition. They attributed southern sediment depocenters to high production of the ancient James and Potomac Rivers draining the central Appalachians (Fig. 11) versus sources from the northern Appalachians.

During the Aptian, the depositional environment of the NBCT transitioned from predominantly fluvial to deltaic influenced, as also observed in the southern BCT (Fig. 11; Schmelz et al., 2020b). Blocky, funnel-shaped gamma-log HST sandstones and transgressive-regressive parasequences characterize the Aptian to early Cenomanian sequences (Fig. 2). This log signature has been tied to limited cores at COST B-2 and interpreted as deltaically influenced (e.g., Miller et al., 2018). HST sandstone intervals tend to thicken and coarsen upsection in the Logan Canyon sequences; there are sharp contacts at the upper sequence boundaries including termination or reset of stacking patterns and sharp shifts in well-log values that correlate to seismic sequence boundaries (Figs. 2A and 2B). At COST B-2, Aptian sediments consist of porous, clean sandstones with calcareous mudstones containing lignite and muddy limestones (e.g., Scholle, 1977). Sedimentological features described in cores from COST B-2 in the mid-Cretaceous sequences indicate a marine deltaic setting (Miller et al., 2018). Seismic profiles show clear reflection terminations and high-amplitude, continuous seismic reflections, and variably defined clinoform geometries (Fig. 3 with Table 2).

The Aptian LC3 sequence that overlies the LK1 reflector marks not only a shift from fluvial to deltaic environments; it also marks the beginning of a landward shift in sediment depocenters. The depocenter parallels the modern shelf and slope margin and is thickest offshore southern New Jersey on the modern outer continental shelf (Fig. 9), though the main depocenter for this sequence is offshore Maryland (Schmelz et al., 2020b). Sandstones and mudstones that are present in the northern study area become finer grained downdip and to the southwest. The main sediment source for the Aptian sequence in the NBCT likely follows the trend of the thickest sediments through the area, coming from the WNW (i.e., a central Appalachian source; Fig. 9).

The Albian LC2 sequence that overlies the MK3 reflector sedimentary depocenter is thin and difficult to distinguish in the NBCT, implying a distal sediment source or a low sediment supply. Two sources for LC2 appear likely (Figs. 9 and 11): from the south-southwest and northwest (Fig. 9). During the Albian, a northern source of sediment began to develop. Like the Aptian sequence, the Albian sequence is thickest in the southern BCT, indicating a proximity to a southern source (Schmelz et al., 2020b) draining the central Appalachians via an inferred paleo–Susquehanna River (Poag and Sevon, 1989). Sandstones and mudstones that are present updip and to the northeast become progressively muddier to the southwest. Poag and Sevon (1989) attributed latitudinal depocenter shifts during the Barremian and Aptian–Albian to differential tectonic uplift among the source terrains. Miller et al. (2017) partially attributed variations in sediment thickness on the New Jersey Coastal Plain during this time interval to changing sediment sources but indicated that accommodation changes by mantle dynamic topography can explain the onshore shift of depocenters from the south during the Early Cretaceous to the north in the Late Cretaceous.

A clear spatial shift in depocenters is observed for the early Cenomanian LC1/MK2 sequence. This sequence is significantly thicker (0–0.3 s; on average 400 m), broader, and has a generally higher sedimentation rate than the older sequences in the NBCT. The higher sedimentation rates help resolve ~2 m.y. scale sequences not resolvable in the older section. Sandstones that are present on the GSD become finer grained downdip and to the southwest. The sediment depocenter for this sequence is on the inner continental shelf ~20 km southwest of the GSD and thins significantly downdip. This landward shift in depocenter, coupled with a shift from terrestrial to marine deposition (Sikora and Olsson, 1991; Miller et al., 2018) indicates a long-term sea-level rise from the Early to mid-Cretaceous that is seen in regional backstripped estimates (Kominz et al., 2008) and in global estimates (Haq, 2014; Fig. 11). The landward depocenter shift correlates with peak gamma-ray values in all NBCT well logs and to a maximum increase in regional (Kominz et al., 2008) and inferred peak global mean sea level (Haq, 2014) near the Cenomanian/Turonian boundary.

The Cenomanian/Turonian boundary peak in sea level is global (e.g., Gale et al., 2008) and has been previously correlated with Ocean Anoxic Event (OAE) 2 near the Cenomanian/Turonian boundary (Arthur et al., 1988). We note that this interval consists of mudstones with peak gamma-log values in NBCT wells, likely indicating black claystones. Miller et al. (2018) showed COST B-2 cores with black shaley claystones occur at the MFS of the LC3 sequence. We admittedly lack core ground truth for the gamma-log peak observed at the Cenomanian/Turonian boundary but argue it is a reasonable inference from stacking patterns and regional correlations (e.g., high organic carbon claystones also occur at the Cenomanian/Turonian boundary onshore at ODP 174AX Bass River Site; Sugarman et al., 1999). In addition, regional analysis shows that the sedimentary depocenter, and likely sediment source, began to shift north during this time (Schmelz et al., 2020b). Similar patterns are observed onshore, where Aptian–Albian sequences (Potomac I and II) thicken to the south, and during the Cenomanian (Potomac III) to Campanian, there was a gradual shift in depocenters from the central to northern New Jersey coastal plain (Kulpecz et al., 2008; Miller et al., 2017).

The Upper Cretaceous (DCx/MK1) composite sequence is predominantly mudstone and is the thickest sequence in the NBCT (up to 1 km). It is marked by a series of high-angle, prograding clinoforms that transition upsection to flatter, continuous reflections above seismic sequence boundary MK1. We interpret these seismic facies as primarily shelfal marine deposits (Fig. 8; Table 2), because the delta shifted landward beneath the modern onshore coastal plain in the Late Cretaceous (e.g., Kulpecz et al., 2008; Miller et al., 2017). The depocenter shifted from the inner continental shelf during MK2/LC1 time (early Cenomanian; Fig. 9) to the outer continental shelf during the later Late Cretaceous (Fig. 1). The MK1/DCx sequence boundary (ca. middle Cenomanian; Fig. 10) marks the transition from shelf depocenters of the Aptian–lower Cenomanian sequences to outer shelf and slope depocenters of the Upper Cretaceous sequences. Our isopach and environmental interpretations (Fig. 9) suggest that the sediment source in the study area is from the NE at MK2/LC1 time (Cenomanian) and remains the primary source of sediment for the BCT in the Late Cretaceous.

Homogeneous, high gamma-ray mudstones of the Upper Cretaceous provide a thick confining unit for the Logan Canyon Sandstone. These confining mudstones are interrupted by a thin, blocky, low gamma-ray Campanian to Coniacian sandstone interval at the top of the DCx sequence. The lithology for this interval consists of dominantly dark-gray mudstone at COST B-3, and pyritic claystone with small amounts of calcareous siltstone and sandstone at COST B-2. Upsection at COST B-2, claystones transition to interbedded calcareous sandstones, siltstones, and claystones with abundant foraminifera (Poag, 1985). A thick package of steeply dipping, progradational seismic reflections within the MK1 seismic sequence transitions to thinner, more shingled reflections within the UK2 seismic sequence that is capped by the high-amplitude, regionally observed UK1 reflection (Fig. 8). We interpret the Upper Cretaceous interval to be deposited in a predominantly marine environment, with deltaic influence during the Coniacian–Campanian as observed onshore (Kulpecz et al., 2008; Miller et al., 2017). The interpreted generalized depositional environments for the NBCT indicate a general transgression throughout the mid-Cretaceous, with shoaling after the OAE2 event and sandstone deposition in the Coniacian–Campanian.

The sequences identified in this study are primarily 10 m.y. (stage-level) in duration, except for the LC1 and DCx sequences that are ~1–2 m.y. in duration (Fig. 11). Unconformities and inferred hiatuses occurred across the Barremian/Aptian boundary, Aptian/Albian boundary, Albian/Cenomanian boundary, and in the mid-Cenomanian and mid-Turonian (Fig. 11). The stage-level sequences (e.g., LC3 Aptian and LC2 Albian) reflect longer-term tectonic and global “first- and secondorder” sea-level variations (Vail et al., 1977). We suggest that the breaks across the stage-level sequence boundaries (Barremian/Aptian and Aptian/Albian) are global and associated with sea-level falls (Ray et al., 2019), though we cannot evaluate if the mechanism that causes these stage-level, sea-level falls was the same as the 1–2 m.y. scale falls.

The 1–2 m.y. sequences are global and likely caused by small changes in ice volume (Miller et al., 2005; Ray et al., 2019; Plint et al., 2022). The Albian–Cenomanian, mid-Cenomanian, and mid-Turonian breaks can be identified onshore New Jersey, Europe, India, North Africa, and the U.S. Gulf Coast (e.g., Miller et al., 2005; Gale et al., 2019) and, therefore, are highly likely to reflect global sea-level lowerings (Ray et al., 2019). Backstripping of boreholes onshore New Jersey suggest that relative sea-level lowerings of ~10–30 m coincide temporally with the Albian–Cenomanian, mid-Cenomanian, and mid-Turonian hiatuses we observe herein (Fig. 11; Kominz et al., 2008). Changes on that scale (“third-order’) and higher (e.g., ~405 k.y.; Plint et al., 2022) are likely global and forced by long tilt (1.2 m.y.) and very long eccentricity (2.4 m.y.) cycles (Boulila et al., 2011) due to growth and decay of small continental ice sheets (<30 m sea level) associated with Albian–Cenomanian, mid-Cenomanian, and mid-Turonian events. We do not observe Albian and older third-order events; resolution of Albian–Cenomanian, mid-Cenomanian, and mid-Turonian events are possible due to the high sedimentation rates at those times.

Our results show that although long-term, global sea level is the primary control on accommodation (Fig. 11; e.g., the longer-term rise and fall in the middle to Late Cretaceous), proximity to sediment input and tectonic controls play an important role in determining the composition, thickness, and continuity of stratal packages in the BCT (Figs. 9 and 11). Though the m.y. scale hiatuses from ca. 100–90 Ma (Fig. 11) may be explained by moderate (tens of m) global sea-level changes driven by ice volume changes (Ray et al., 2019), the sea-level processes associated with the major sequence boundaries we observe on the 10 m.y. scale (e.g., Barremian/Aptian and Albian/Aptian hiatuses) are uncertain. These sequence boundaries coincide with standard stage boundaries and the second-order sequences of Haq (2014). Thus, it seems likely the stage-level unconformities are global, but it is entirely possible that they were driven by regional tectonics since their timing is not well constrained.

Subsurface storage of supercritical CO₂ requires a porous and permeable reservoir with a strong seal and cap and burial depths of 800 m or more below the seafloor (Bachu, 2002). The Logan Canyon sandstones offshore New Jersey are porous (20%–40%), permeable (tens to hundreds of mD), below 800 m burial, and capped by thick (hundreds of m) mudstones. In addition, the extensive seismic data set we used has led to the characterization of these reservoirs and seals. Using a sequence stratigraphic approach, we have identified HST and LST sandstones throughout 29 wells in the northern BCT and used seismic profiles to trace the sequences throughout the basin. Although seismic resolution is too low to distinguish between individual systems tracts in the available seismic data sets, seismic sequences of the Logan Canyon Forma tion can be traced across the BCT to provide a means to calculate sequence thickness where wells are unavailable.

Sandstones within the Cretaceous Logan Canyon sandstone sequences are excellent reservoirs that provide thick, porous (<20%) and permeable (>10 mD) reservoirs with multiple confining units. Blocky sandstones of the HST are the thickest and most laterally continuous, with thickness of individual sandstone bodies from 20 m to over 100 m thick. Porosity and permeability measurements suggest excellent storage capability, particularly in these blocky HST sands. As shown in Figure 12, there are multiple caps and seals for these target reservoirs, with the DCx composite sequence providing a regionally extensive, tight seal for the LC sequences below. We cannot fully evaluate the quality of the confining nature of the mudstones within the Logan Canyon Formation from the well logs, and some of the mudstones may in fact be leaky siltstones and not fully confining claystones. However, the DCx upper Cenomanian to lower Turonian sequence appears to be primarily claystone with peak organic content in OAE2 and provides a thick, continuous cap rock and seal for the reservoirs below. Much of the overburden (Fig. 12; especially Maastrichtian to Eocene rocks) are also fine-grained clays and silty clays (equivalent to the Composite Confining Unit; Zapecza, 1989) onshore, providing greater confidence that once injected CO₂ remains buried.

The flanks of the GSD provide the ideal future CCS target (Fig. 12). Differential compaction, forming an anticline around the GSD, has created a structural trap that provides excellent closure for storage. The HST and possibly LST sandstones within the mid-Cretaceous sequences in the vicinity of the GSD appear to be the optimum storage location for supercritical CO₂ in the NBCT and a world-class target for carbon sequestration. We illustrate (Fig. 12) an example of seabed pumping of supercritical CO₂ into the LC2 highstand blocky sandstones on the GSD (Fig. 2A).

We not only see the offshore GSD as a target but also the nearshore region particularly in southern New Jersey to Maryland (Schmelz et al., 2020a) as viable and substantial carbon storage reservoirs. The nearshore zone has the advantage of shorter distance to pipe or transport via ship supercritical CO₂, though the GSD has proven reservoirs, seals, and has the advantage of proximity to mid-Atlantic Bight point sources. The offshore location may be more acceptable to local stakeholders.

Modeling the total cost of storage (including capture and transportation costs) for these reservoirs on the GSD or offshore Maryland suggests their utilization for carbon storage is economically feasible, in addition to politically preferable to storing carbon dioxide at onshore sites (Schmelz et al., 2020a). It should be recognized, however, that firm estimates of storage volume are not possible at present, as published estimates attest. Monteverde et al. (2011) estimated large amounts (14–69 Mt) of CO₂ could be stored in the Logan Canyon Formation in a 50 km × 30 km area of the outer continental shelf and slope, not including the GSD. Fukai et al. (2020) used storage efficiencies of 1% to 13%, to estimate much higher regional-scale prospective storage resources of 37–403 (Gt) of CO₂. They also estimated that 30–51 Mt of CO₂ can be stored using a dynamic injection simulation into a Middle Cretaceous sequence on the eastern flank of the GSD.

Our work has provided better constrained estimates of potential storage volume in the NBCT. We created regional isopach and structural contour maps for the LK1, MK3, MK2, and MK1 sequence boundaries and the correlative Potomac III, Potomac II, Potomac I, and Bass River sequence boundaries delineated onshore by Miller et al. (2017). We calculated volumes by multiplying the fraction of the isopach thickness that is a sandstone reservoir by the study area, the reservoir porosity, the CO₂ density at formation temperature and pressure, and the formation storage efficiency. The fraction of the isopach that is a sandstone reservoir was derived by determining the presence of sandstone within each sequence using gamma-ray cutoffs in normalized well logs (Fukai et al., 2020). The reservoir and sandstone porosity was taken from well logs where available or assumed to be 20% where well data were not available. We calculated CO₂ density at formation temperature and pressure using the map of those variables defined by Bachu (2002) and assumed a 25 °C/km geothermal gradient (Nathenson and Guffanti, 1988) and a 10.2 MPa/km pressure gradient. Storage efficiency was assumed to be 2%. With these parameters, we estimate that 21 Gt of supercritical CO₂ can be safely stored regionally in the Logan Canyon reservoirs and the correlative Potomac Formation onshore. This estimate is several orders of magnitude greater CO₂ storage volume than Monteverde et al. (2011) calculated since we include the whole region. Our estimates are also similar to Fukai et al. (2020) from a volume per square km perspective, for a given efficiency, but the spatial distribution of the storage volumes are better constrained by our seismic and well-log data for the NBCT.

Our studies show that the Logan Canyon sequence sandstones provide large reservoirs for carbon storage, but we do not resolve reservoirscale (km to tens of km) distributions of sands. We also do not evaluate the influence that supercritical injection will have on stimulation of earthquakes (Zoback and Gorelick, 2012) or concerns of migration of supercritical CO₂ away from the injection wells. An initial geocellular earth model and reservoir simulation model of carbon injection in the GSD region by Schlumberger Carbon Services (Brown et al., 2011) shows only modest increases in changes in geopressures caused by injection, minimizing concerns of rock fracturing raised by Zoback and Gorelick (2012). The initial Schlumberger model also shows that vertical and lateral migration pathways for the injected plume have only km-scale migration on a supradecadal time window. Our detailed sequence stratigraphic interpretations will allow more detailed studies of potential reservoirs. Future work will build on the broad regional-scale modeling of CO₂ storage in target reservoir sandstones, by making reservoirscale predictions of the flow of supercritical CO₂ through those formations.

We develop an updated sequence stratigraphic framework for the Cretaceous sequences of the BCT and identify well-log sequences and composite sequences in Cretaceous strata (Miss, LC3, LC2, LC1, and DCx) that correlate (±50 m) with seismic sequences (LK2, LK1, MK3, MK2, and MK1/UK2, respectively). The Aptian LC3/LK1 sequence contains thick, blocky sandstone on the GSD; it becomes finer grained downdip and to the southwestern part of our study area. This sequence is thickest south of the GSD, near the structural high of the DEH. The Albian LC2/MK3 sequence is thinner than the Aptian sequence and contains thick, blocky sandstones in the northern part of the study area that become finer grained to the southwest. The lower Cenomanian LC1/MK2 sequence is the thickest of the three sequences, is located on the inner continental shelf, and thins and loses sandstone content downdip.

The LC3/LK1 and LC2/MK3 sequences are ~10 m.y. in duration and reflect longerterm tectonic and global sea-level variations, while the LC1/MK2 sequence is ~2 m.y. in duration and has high sedimentation rates that resolve higher- frequency base-level changes, likely generated by global sea-level variations. We observe an overall long-term transgression from the Early to Late Cretaceous based on shifting fluvial to deltaic environments, deepening paleo water depths, changing clinoform geometries, and correlation to global events. This transgression culminated in maximum flooding corresponding with the deepest Cretaceous paleodepths of the BCT and global OAE2 at the Cenomanian/Turonian boundary. During this longterm Cretaceous transgression, we observe a shift from southern to northern sourced deposition. Sediment depocenters also shift from the outer to inner continental shelf. Following the OAE2 event, sediment depocenters shift seaward, and seismic facies exhibit flatter, shingled clinoform geometries, indicating a lowering of long-term sea level during the Late Cretaceous.

Although global sea level, thermal subsidence, and loading provide accommodation in the BCT, proximity to sediment input and tectonic controls determine regional and local deposition. Due to differential compaction around the GSD, we observe thinner sequences on the structural high and thicker sequences in the flanks. Depocenters are also found near major faults, such as the structural DEH. In addition, we find the thickest sandstone deposits on or near the GSD.

Our results show that the mid-Cretaceous strata on the GSD and nearshore Maryland within the BCT are a world-class target for carbon sequestration. Highstand sandstones within the Logan Canyon sequences are porous and permeable, capped by thick mudstones, and are constrained by the structure of the Barremian igneous intrusion and subsequent differential compaction that formed the GSD.

This work was supported by the U.S. Department of Energy under award numbers DE-FE0026087 (Mid-Atlantic U.S. Offshore Carbon Storage Resource Assessment Project) and DE-FC26-0NT42589 (Midwest Regional Carbon Sequestration Partnership Program) managed by Battelle. This study would not have been possible without the efforts of many to create a useable data set within the Baltimore Canyon Trough. We thank L. Cummings and N. Gupta of Battelle for their guidance over the course of this project. We acknowledge G. Lang (USGS and Haifa University) for first recognizing the Delaware High and thank him for his insight into regional seismic correlations. We thank W. Fortin for leading the effort to reprocess legacy seismic data. B. Dunst and K. Carter (Pennsylvania Geological Survey) were very helpful in providing insight into and digitizing well logs. We thank M. KunleDare and P. McLaughlin of the Delaware Geological Survey (DGS) for their work organizing and cataloging legacy well-log and core data used in this study. We would also like to acknowledge C. Lombardi (deceased), who began this work on the NBCT and first recognized a need for an updated sequence stratigraphic framework of the basin. We thank D. Hodgson and A. Findani and two anonymous reviewers for comments.