Locations of the 1982 Miramichi (Canada) Aftershocks: Implication of Two Rupture Regions Activated

Abstract

On 9 January 1982, in the Miramichi region of New Brunswick, Canada, an earthquake with mb 5.7 occurred. It was followed by extensive aftershocks and felt throughout eastern Canada and north-eastern USA. Since this earthquake occurred in an uninhabited region, the damage was minor. Due to an mb 5.7 event is rare in north-eastern America, investigating it and its aftershocks is important for understanding intraplate seismicity. Digital seismic stations were not yet common by 1982. Fortunately, four seismic phases at three stations could be used to locate larger aftershocks. A simplified master-event location method combined with regional depth-phase modeling was used to locate aftershocks. For each aftershock its focal depth was first determined using a depth phase; then, with the depth fixed, the epicenter was determined using the four arrival time readings measured at the same three stations. The located aftershocks were divided into three groups. In each group the earthquake numbers are similar, but the majority of the energy was released in one group. The epicenters formed two trends in the NE–SW direction, implying that the Miramichi earthquake sequence activated two rupture regions.

1. Introduction

The 9 January 1982 Miramichi, New Brunswick, magnitude (mb) 5.7 earthquake was a rare case in North America (Figure 1). It was felt throughout eastern Canada and northeastern USA and intrigued scientists and the public as it was the largest one in eastern Canadian and eastern US in recent 100 years. The mainshock (mb 5.7; M_W 5.6) occurred at 12:53 UT on the 9th January, followed 3.5 hours later by a large mb 5.1 (M_W 4.9) aftershock. On the January 11th the largest aftershock (mb 5.4, M_W 5.0) followed, then on March 31st another large aftershock (mb 5.0, M_W 4.9) occurred. The above three large aftershocks were called principal aftershocks. The mb magnitude was used for the mainshock and the three principal aftershocks in the majority of the publications and the media for many years. The moment magnitudes can be found in the report by Bent [1].

Figure 1.

The location of the 1982 Miramichi earthquake sequence, the seismicity occurred from 1980 to March 2022, and the geological background in its surrounding regions. The triangles show the locations of seismic stations, KLN, EBN, GGN, LMN, LPQ , and GSQ. The diamonds show the locations of cities or towns. The symbols CSZ and LSZ are abbreviations for the Charlevoix seismic zone, and lower St. Lawrence seismic zone. A solid circle color coded with focal depth shows an earthquake epicenter, which were retrieved from the incorporated research institutions for seismology (IRIS). A focal depth value is indicated by the depth scale on the right. The color scale at the bottom shows the height of the topographic locations above the sea level. The St. Lawrence River runs through Quebec City, Charlevoix seismic zone, and lower St, Lawrence seismic zone. The St. Lawrence faults system also runs along this trend (e.g., [2]); for clarity it is not plotted. At the north of CSZ is the Saguenay Graben. The 1988 M_W 5.9 earthquake occurred along this Graben. Figure 1 was plotted using the GMT program [3].

Three field surveys were conducted in 1982 by the Geological Survey of Canada (GSC) to investigate the aftershock sequence. In the January survey (S1), the most detailed coverage of the aftershock activity was from 19 to 22 January (the temperatures were below -25C°) when aftershocks were recorded by analog MEQ-800 seismographs at four sites within 10 km of the active zone. The hypocenter of the mainshock was estimated using the hypocenters of the detected small aftershocks. The April survey (S2) was conducted in response to the 31 March mb 5.0 aftershock, whose hypocenter was also estimated using the hypocenters of the detected small aftershocks. The survey in June (S3) followed the 16 June mb 4.7 earthquake (Figure 1). As this event was located about 30 km west of the Miramichi mainshock (e. g. Wetmiller et al., [4]), it is not discussed in this article.

Responding to a request from Canada, the U.S. Geological Survey (USGS) installed a portable digital network. This network located about 40 aftershocks between the 15 and 22 January 1982 [5]. Among the 40 aftershocks, 4 larger ones were relocated and their focal mechanisms were studied by Saikia and Herrmann [6].

The focal mechanism of the mainshock was a thrust type (e.g., [7]). The rupture was inferred to be updip on a west dipping NNE striking fault plane (Choy et al., [8]). The “beach-ball” is plotted using gCMT data (the global CMT project; globalcmt.org), and the inferred rupture plane is labeled p2 (seeFigure 2).

Figure 2.

The shift-corrected epicentral distribution of the located 68 aftershocks. The size of a solid circle is proportional to the magnitude, while the color matches the focal depth (see the depth scale on the right). Star S1 marks epicenter of the mainshock and star S2 marks the epicenter of the mb 5.0 aftershock, determined by Wetmiller et al. [4]. The coordinate point (0, 0) is at (47.0°N; 66.6°W). The aftershocks in the upper part of the figure were separated into two groups (the left and right groups) by a gap region indicated by a dashed-line with an arrow at Az 38°. The aftershocks in the lower part were included in the bottom group. The diamond shows the epicenter of the mb 5.4 aftershock, located using phase readings at stations EBN, GGN, LPQ , and GSQ. The triangles indicated with 5.7, 5.1, 5.4, and 5.0 show the epicenters of the mainshock and its three principal aftershocks, located by Choy et al. [8]. The epicenter of the mainshock was calibrated to that obtained by Wetmiller et al. [4]; accordingly, the epicenters of the 3 principal aftershocks were moved with the same amount of the distance and direction as those of the mainshock. The “beachball” shows focal mechanism of the mainshock calculated from the moment tensor solution (the global CMT project; globalcmt.org). The nodal plane indicated with p2 is the inferred the rupture plane. The two arrows pointed to the “beachball” show the compressive force direction in the source region.

As this earthquake sequence occurred in an almost completely uninhabited region, the damage was minor. However, investigating this mainshock and its aftershocks is important for understanding intraplate seismic activity, and assessing the seismic hazard in the source region and its vicinity for the future.

Since there was no close Canadian digital seismic station, a new station (KLN) was installed by GSC on 23 January 1982, 14 days after the mainshock, to better monitor the sequence. The station belonged to the Eastern Canadian Telemetered Network (ECTN). KLN recorded hundreds of aftershocks, and the waveform record quality was excellent. Two existing ECTN stations, EBN and GGN, also had clear records for the larger aftershocks (m_N ≥ 2.8; the magnitude m_N was defined by Nuttli, [9]). Figure 1 shows the locations of these three stations, as well as those of stations LPQ and GSQ.

Between latitudes 46.88° N – 47.16° N and longitudes 66.35° W–66.80° W, there were about 700 aftershocks (the smaller aftershocks detected in the field surveys are not included) in the Natural Resources Canada (NRCan) catalog database. Ma and Motazedian [10] determined the focal depths for more than 100 aftershocks with m_N ≥ 2.8, using depth phase sPmP recorded at EBN, but left the epicenters unchanged. Most of the aftershocks in the database were assigned the same epicenter (47.00° N, 66.60° W), which is the epicenter of the mainshock, determined by Wetmiller et al. [4].

For an earthquake with Pg and Sg arrival readings at KLN, Pg at EBN, and Pn at GGN, a conventional location method can in principle be used to determine its hypocenter. However, the value of an earthquake’s focal depth is usually much smaller than those of the station distances; for a seismic phase recorded at a station, the depth typically has a much smaller contribution to the travel time than that from the station distance. As such, the error (uncertainty) in the depth is much larger than the uncertainties in the latitude and longitude of an epicenter. To reduce the error in an epicenter, a focal depth can be first determined using a depth phase; then, the epicenter is located at the depth determined.

Since the station coverage for the sequence was not good (Figure 1), and regional velocity models are not good either, it was not possible to determine an epicenter for an aftershock with small absolute errors. However, errors in the relative locations in a small aftershock group should be smaller and can be obtained using a master-event method (e.g., [11]).

Since the energy released by the mb 5.4 aftershock is of the same order as that of the mainshock, the Miramichi earthquake is also called a double-earthquake. The mainshock and its 3 principal aftershocks were relocated by Choy et al. [8]. The 4 focal depths were determined using a waveform modeling method [10]. The depth of mb 5.7 was 6.8 km, mb 5.1, 5.5 km, mb 5.4, 5.2 km, and mb 5.0, 2.0 km. The 4 depths were progressively shallower with occurrence times.

The durations of the surveys conducted by GSC or USGS were shorter than one week; the epicenters of the aftershocks available from IRIS database have large uncertainties. This implies that no clear patterns for the Miramichi earthquake sequence are available yet.

Our motivation was to obtain a reliable pattern of the hypocentral distribution by locating larger aftershocks. In the following sections we briefly introduce the seismicity and geological background in the vicinity of the Miramichi earthquakes; analyze the waveform data; briefly introduce the methods for locating aftershocks; present the hypocentral distribution features of the located aftershocks; display the time series and strength of the three earthquake groups obtained; analyze the errors in the relative locations between two adjacent aftershocks; and discuss some related issues.

2. Seismicity and geological background in the vicinity of the Miramichi earthquakes

The Miramichi earthquake sequence locates in the eastern part of the Appalachian Mountain range (Figure 1). The mountain range is mostly in the United States (US). It forms a zone from 160 to 480 km wide, running from the island of Newfoundland southwestward through New Brunswick, Canada to Central Alabama, US. In the north side, past the mountain range, is located the St. Lawrence River. Along this river, the St. Lawrence Faults system developed (e.g., [2]). In the northeast side of the Quebec City, about 100 km away is the Charlevoix seismic Zone (CSZ), which is the most active earthquake zone in eastern Canada. In the history, some destructive earthquakes occurred in the CSZ (seismescanada.rncan.gc.ca/historic-historique/events/18701020-en.php). Adjacent to the CSZ at the northern side is the Saguenay Graben. The 25 November 1988 M_W 5.9 earthquake occurred along this Graben (e.g., Ma et al., [12]). Downstream from the CSZ, is the Lower St. Lawrence Seismic Zone (LSZ). The largest earthquake detected is the 16 March 1999 m_N 5.1 earthquake (e.g., Lamontagne et al., [13]).

3. Waveform data analysis

Figure 1 shows the stations we used to locate the Miramichi aftershocks. The closest one to the sequence (about 25 km) was KLN. At this station clear Pg- and Sg-phase were recorded for almost all aftershocks. Figure 3 shows 10 seismograms generated by 10 aftershocks, recorded at KLN. The top 5 traces were generated by the aftershocks in the right group (see Figure 2), while the bottom 5 traces by aftershocks in the left group. The T_Sg– T_Pg times along the top 5 traces are shorter than those along the bottom 5 traces; due to those epicenters to KLN are shorter (refer to Figures 1 and 2).

Figure 3.

Vertical component short period displacement seismograms recorded at KLN (along a displacement trace, sPg phase is easier to identify than along a velocity trace), generated by 10 aftershocks in Table 1. All records are aligned at the Pg phase. The time on the left side of each record is the raw record start time. For each aftershock, the hypocentral distance to KLN is mainly constrained by the time difference T_Sg – T_Pg, while the focal depth is mainly constrained by T_sPg– T_Pg. The top 5 records were generated by the aftershocks in the right group, while the bottom 5 records by aftershocks in the left group (see Figure 2). The 10 aftershocks were selected along a line which is orthogonal to the located epicenter distribution trend (about Az 38°, see Figure 2). The T_Sg – T_Pg times along the top 5 records are shorter than those along the bottom 5 records; due to those epicenters to KLN are shorter. The time difference δt corresponds to the spatial gap between the two groups.

The second closest station was EBN (about 135 km from the sequence). At EBN, the onset of the Pg phase was usually clear on the seismograms generated by aftershocks with m_N ≥ 2.8, so the arrival time of the Pg phase could be measured. Figure 4 shows vertical displacement seismograms recorded at EBN. The trace indicated with mb 5.4 was generated by the mb 5.4 aftershock. The trace indicated with mb 5.0 was generated by the mb 5.0 aftershock. The other 4 traces in the top panel were generated by 4 aftershocks in the left group; the 5 traces in the bottom panel were generated by 5 aftershocks in the right group (see Figure 2). The T_Sg– T_Pg times (about 16.54 s) along the traces 2, 3, 4, and 5 in the top panel, are shorter than those (about 17.04 s) along the bottom 5 traces; due to those aftershocks in the left group are closer to EBN. The T_Sg– T_Pg times along the traces generated by the two principal aftershocks in the bottom panel are approximately equal; implying that they are in the same group.

Figure 4.

Vertical component short period displacement seismograms recorded at EBN (station distance ~ 135 km) generated by 9 aftershocks in the Miramichi sequence. The trace indicated with mb 5.4 was generated by the mb 5.4 aftershock. The trace indicated with mb 5.0 by the mb 5.0 aftershock. The other 4 traces in the top panel were generated by 4 aftershocks in the left group; the other 3 traces in the bottom panel were generated by 3 aftershocks in the right group (see Figure 2). For each aftershock, the hypocentral distance to EBN is mainly constrained by the time difference T_Sg – T_Pg. The T_Sg – T_Pg times (about 16.54 s) along the traces 2, 3, 4, and 5 in the top panel are shorter than those (about 17.04 s) along the bottom 5 traces; due to those epicenters to EBN are shorter.

The third closest station was GGN (station distance ~200 km). At this station, the Pn phase generated by aftershocks with m_N ≥ 2.8, was clear (see Figure 5).

Figure 5.

Travel time comparisons. The upper three traces were generated by the master event (ME, No. 36 in Table 1); the bottom three traces were generated by a secondary event (SE, No. 44 in Table 1). The time marked along the bottom axis is relative for conveniently aligning the traces. The symbol t0 means origin time; ∆Pg, ∆Sg, and ∆Pn are travel time differences. The origin times of the ME and SE were aligned for comparison (the determined values of the two origin times are in Table 1; here they were calculated using Sg - Pg times). The travel times for the ME at KLN: Pg phase is 4.00 s and Sg phase is 6.97 s; at EBN: Pg is 22.15 s; at GGN: Pn is 31.07 s. the travel times for the SE at KLN: Pg phase is 3.58 s and Sg phase is 6.22 s; at EBN: Pg is 22.65 s; at GGN: Pn is 30.71 s. At KLN: ∆Pg is + 0.42 s (4.00–3.58) and ∆Sg is + 0.75 s (6.97–6.22); at EBN: ∆Pg is − 0.50 s (22.15–22.65); at GGN: ∆Pn is + 0.36 s (31.07–30.71). These differences determine the relative positions of the two aftershocks. These traces are vertically enlarged, horizontally expanded, and interpolated to 100 points/s, the arrival times could be measured with a precision of 0.01 s.

The waveforms at station LMN (station distance ~200 km; Figure 1) were also analyzed. The Pn phase, generated by aftershocks with m_N ≥ 2.8, was clear when the seismograph at this station operated normally. Unfortunately, in 1982, the seismograms at LMN were often rectangular pulses with different amplitudes (probably due to instrument malfunction), so the Pn arrival times could only be measured for some aftershocks. As such, the waveform records were not used for the locations.

The Pn phase at stations LPQ and GSQ are also clear for aftershocks with magnitude m_N ≥ 3.5. The Pn phase at these two stations and at GGN, as well as Pg and Sg at EBN were used to locate the mb 5.4 aftershock.

4. Methods

4.1 A simplified master-event relocation method

Master-event relocation methods have been studied and used by many scientists (e.g., [11, 14, 15, 16]). The steps in the master-event method described by Havskov and Ottemöller [11] are as follows: (1) Locate the master-event (ME) using a conventional location method; (2) Select stations and phases which are common to the ME and the slave-events (SEs); (3) Calculate the residuals at the selected stations for the ME; (4) Add the residual of a certain phase to the arrival time readings of the same phase for the ME and SEs; and (5) Relocate all events (the ME and the SEs) using a conventional event locating program.

The step (2) is necessary. We tested that once this necessary step was kept the uncertainty between two adjacent epicenters can be small, and the epicenter distribution pattern can be reliable, so we only kept this step in our work. We used the program in SEISAN (e.g., [17] to locate the aftershocks.

4.2 A depth phase modeling procedure to determine a focal depth

When the Pg- and Sg-phase arrival times at KLN, Pg at EBN, and Pn at GGN are available for an aftershock, theoretically the four source parameters (origin time, latitude, longitude, and focal depth) can be determined using the four time readings. However, the uncertainty in the four parameters could be large, especially that in focal depth. In practice, if only a few arrival time readings are available or station coverage is poor, the focal depth is assigned a nominal value to stabilize the location when locating an event.

To reduce uncertainties in epicenters due to uncertainty in focal depth, the time difference along a trace between depth phase sPg and its reference phase Pg (Figure 3) is used to first retrieve a reliable focal depth for the aftershock that generated the trace; then, the epicenter of the aftershock is located at the focal depth retrieved. In this way, the trade-off between the epicenter and the focal depth is removed, so the uncertainty in the epicenter can be dramatically reduced.

The crucial step in the procedure to retrieve a focal depth using a depth phase is the generation of the synthetic traces along which the depth phase appears. In generating synthetic waveforms, a crustal model, station distance, focal mechanism, and focal depth are needed input parameters. Since the crustal structures through which the waves travel are related to travel times, the crustal model is a key input parameter. The reflectivity method [18], the centroid moment tensor solution for the Miramichi mainshock from the gCMT Catalog (see data and resources section), and the crustal model introduced in next section were used to generate the synthetic traces. The details of the depth phase studies can be found in e.g., [19, 20, 21, 22].

5. Crustal model

There are several studies on crustal models for eastern Canada (e.g., [23, 24]). Rayleigh-wave dispersion data from the 23 June 2010 M_W 5.2 earthquake about 60 km northeast of Ottawa, Ontario [25], were used to obtain 14 crustal velocity models around the epicenter [24]. The Rayleigh wave travel paths for model No. 8 ran approximately through the Miramichi region.

In recent years, some shallow, small earthquakes occurred in the Miramichi region. Those small earthquakes generated Rg–wave (Rayleigh waves traveling in the crust) records. Using the Rg–wave dispersion data, models for the shallow part of the crust (0–10 km) were obtained [26]. A velocity model was formed for the Miramichi region by replacing the shallow part of the model by Motazedian et al. [24] with the model reported by Ma [26]. The Vp/Vs ratio is assumed to be 1.74.

6. Aftershock locations

6.1 Location of the ME and the adjustment to the crustal model

We selected a m_N 3.5 aftershock (No. 36 in Table 1) as an ME. Its waveforms were presented in the upper three traces in Figure 5. This aftershock had very clear onsets of phases Pg, Sg, and Pn; the arrival time readings can be accurate. A focal depth of 4.5 km was previously estimated using the depth phase sPmP, recorded at EBN, and the epicenter (47.0°N, 66.6°W) of the mainshock [10]. However, the sPmP – PmP time could not be accurately measured at EBN (see Figures 5 and 6 in [10]), resulting in an uncertainty of about 1.0 km in focal depth. In this article, we increased the depth accuracy by using the depth phase sPg at KLN. Figure 6 demonstrates the depth phase sPg modeling for the ME. The top trace U_D/5.5 was generated using a depth of 5.5 km; other synthetic traces were generated with a depth increment of 0.1 km. The sPg - Pg time along trace U_D/5.9 and the time difference along the observed trace were approximately equal, so the modeled focal depth was 5.9 km. As the arrival times of Pg and sPg could be precisely compared, the uncertainty in the focal depth obtained using sPg was reduced.

No	Date			Time			Lat.	Long.	H	mN	er_n	er_e	Dis	Modu	Perc
							(°)	(°)	km		km	km	km	km
1	1982	02	24	04	43	01.0	46.977	−66.623	3.1	2.9	3.3	3.0	0.135	0.100	74
2	1982	02	27	17	34	57.6	47.033	−66.566	5.6	3.4	4.2	5.8	1.732	0.608	35
3	1982	03	01	09	33	57.1	46.955	−66.607	4.8	3.4	3.0	2.7	0.451	0.000	00
4	1982	03	03	00	28	32.7	46.996	−66.602	5.7	2.8	3.7	3.8	0.648	0.283	44
5	1982	03	04	06	06	31.1	47.019	−66.576	5.5	2.9	4.3	5.2	1.482	1.000	67
6	1982	03	13	11	38	13.0	46.995	−66.616	5.0	2.9	3.5	3.5	0.470	0.100	21
7	1982	03	13	23	27	51.6	46.974	−66.545	2.8	2.9	3.8	4.3	1.408	0.224	16
8	1982	03	16	11	14	01.8	46.919	−66.620	3.8	3.5	2.6	1.9	0.811	0.141	17
9	1982	03	26	05	36	39.6	46.988	−66.553	5.4	2.8	3.7	4.5	1.212	0.400	33
10	1982	03	26	13	38	07.3	47.006	−66.590	5.3	2.9	3.9	4.3	0.254	0.200	79
11	1982	03	31	21	02	21.2	46.979	−66.562	3.1	5.0	3.7	4.1	0.254	0.100	39
12	1982	03	31	21	29	19.3	46.971	−66.575	2.8	2.9	3.5	3.5	0.254	0.141	56
13	1982	04	02	13	50	12.3	46.977	−66.567	3.5	4.3	3.6	3.9	0.367	0.141	39
14	1982	04	02	19	49	35.4	46.912	−66.623	3.1	3.1	2.5	1.8	2.272	0.283	12
15	1982	04	08	04	54	33.9	46.980	−66.565	3.7	3.4	3.7	4.0	0.367	0.141	39
16	1982	04	11	18	00	53.4	46.980	−66.615	5.7	4.0	3.3	3.1	0.076	0.000	00
17	1982	04	11	18	27	19.3	46.984	−66.609	5.7	3.2	3.4	3.3	0.377	0.141	38
18	1982	04	11	20	07	00.1	46.898	−66.673	5.5	2.9	2.2	1.5	0.397	0.100	25
19	1982	04	18	22	47	21.3	46.958	−66.592	3.7	4.1	3.1	2.9	0.189	0.141	75
20	1982	05	02	23	31	36.9	46.961	−66.616	4.3	3.3	3.1	2.8	0.793	0.141	18
21	1982	05	06	16	28	07.7	46.942	−66.628	4.7	4.0	2.7	2.2	0.720	0.100	14
22	1982	06	16	11	43	30.5	47.001	−66.933	8.7	4.7	3.1	2.2	20.01	1.655	08
23	1982	06	18	11	24	36.0	46.899	−66.668	4.5	3.0	2.3	1.5	1.244	0.000	00
24	1982	06	25	06	47	10.3	46.931	−66.634	6.0	2.9	2.7	2.0	1.518	0.224	15
25	1982	07	18	15	01	04.8	46.940	−66.619	4.6	2.9	2.8	2.2	1.518	0.224	15
26	1982	08	12	20	43	18.2	46.982	−66.613	4.6	3.3	3.3	3.2	0.235	0.100	43
27	1982	09	19	01	37	17.5	46.957	−66.594	6.3	3.1	3.0	2.8	0.319	0.000	00
28	1982	10	18	04	37	48.9	46.967	−66.602	5.7	3.0	3.1	3.0	0.188	0.100	53
29	1982	10	21	18	12	47.6	46.982	−66.617	3.9	2.8	3.3	3.2	0.222	0.000	00
30	1982	10	26	15	31	33.0	46.985	−66.618	4.8	3.5	3.3	3.2	0.135	0.000	00
31	1982	10	28	06	35	10.9	46.842	−66.677	5.9	2.8	2.1	1.4	1.674	0.141	08
32	1982	10	31	12	44	41.4	46.999	−66.561	5.9	2.9	3.8	4.7	2.087	0.412	20
33	1982	12	22	12	53	26.4	46.994	−66.610	4.8	3.0	3.5	3.6	1.114	0.316	28
34	1983	02	12	18	00	25.8	46.886	−66.661	4.5	2.8	2.2	1.4	0.377	0.141	37
35	1983	05	12	20	42	25.4	46.853	−66.656	3.0	3.0	2.0	1.3	1.564	0.141	09
36	1983	05	13	17	26	02.1	46.977	−66.612	5.9	3.5	3.2	3.1	0.269	0.100	37
37	1983	05	13	23	40	57.4	46.982	−66.613	5.8	3.9	3.3	3.2	0.235	0.100	43
38	1983	06	10	04	22	39.2	46.964	−66.626	5.9	3.3	3.0	2.6	1.202	0.361	30
39	1983	06	11	13	47	58.6	46.979	−66.592	6.0	3.4	3.3	3.5	0.601	0.200	33
40	1983	06	28	08	05	49.3	47.033	−66.674	2.7	3.3	3.8	3.7	6.042	0.510	08
41	1983	11	02	06	02	00.4	46.857	−66.679	4.2	2.8	2.0	1.3	1.552	0.000	00
42	1983	11	17	15	32	18.2	46.975	−66.614	5.6	3.7	3.2	3.0	0.470	0.100	21
43	1983	11	18	10	28	39.5	46.862	−66.698	2.2	3.0	2.0	1.3	2.246	0.141	06
44	1984	02	24	03	17	13.8	46.967	−66.586	5.9	3.7	3.2	3.2	0.539	0.283	52
45	1984	03	27	22	56	24.6	46.899	−66.502	3.9	3.0	2.8	2.4	9.261	0.539	06
46	1984	04	13	15	35	51.3	46.974	−66.595	4.7	3.1	3.3	3.3	0.956	0.316	33
47	1984	07	02	05	24	54.0	46.888	−66.665	5.9	3.0	2.3	1.5	1.918	0.224	12
48	1984	08	04	05	11	13.0	46.953	−66.638	2.0	2.9	2.9	2.4	1.526	0.224	15
49	1984	10	13	01	45	15.8	46.823	−66.758	2.0	3.0	1.8	1.2	3.931	0.000	00
50	1984	11	07	19	44	31.5	46.959	−66.606	4.7	2.8	3.0	2.7	0.685	0.100	15
51	1984	11	30	05	54	22.4	46.978	−66.622	5.0	3.8	3.2	3.0	0.470	0.100	21
52	1985	05	13	18	46	19.3	46.854	−66.631	2.0	2.8	2.0	1.3	2.420	0.141	06
53	1985	10	05	05	34	13.6	47.008	−66.587	6.0	3.9	3.7	4.4	0.334	0.100	30
54	1985	10	05	06	17	33.3	47.005	−66.587	6.0	2.8	3.7	4.3	0.334	0.100	30
55	1985	12	21	06	03	10.9	47.010	−66.602	2.3	3.1	3.9	4.2	1.162	0.283	24
56	1986	01	21	02	32	26.2	46.959	−66.597	5.2	3.4	3.0	2.8	0.685	0.100	15
57	1986	03	16	05	01	46.7	46.980	−66.614	4.1	2.8	3.3	3.1	0.235	0.100	43
58	1986	06	01	14	53	14.3	46.968	−66.604	4.3	3.4	3.2	3.0	0.956	0.316	33
59	1986	10	17	14	48	00.1	46.970	−66.578	4.5	4.1	3.4	3.4	0.324	0.000	00
60	1986	10	18	12	24	30.2	46.888	−66.638	4.6	2.8	2.2	1.5	2.058	0.100	05
61	1986	10	23	12	58	04.3	46.839	−66.654	4.4	3.4	2.1	1.4	2.420	0.141	06
62	1987	04	22	14	32	53.0	46.812	−66.689	3.3	2.8	1.8	1.3	2.247	0.100	04
63	1988	03	06	18	13	18.1	46.828	−66.707	2.0	3.2	1.8	1.2	3.841	0.224	06
64	1988	05	09	01	23	05.2	46.974	−66.620	4.4	3.5	3.2	2.9	0.470	0.100	21
65	1988	06	12	18	10	15.9	46.917	−66.671	5.4	2.8	2.4	1.7	3.219	0.424	13
66	1988	08	26	05	59	10.7	46.984	−66.617	4.6	3.8	3.3	3.2	0.377	0.000	00
67	1989	01	16	02	33	56.4	46.971	−66.582	4.2	3.0	3.4	3.4	0.539	0.283	52
68	1989	06	10	10	39	49.5	46.881	−66.688	4.6	2.8	2.1	1.4	2.246	0.141	06
average											3.0	2.8	1.233	0.183	24
69	1982	01	20	10	00		47.003	−66.619	5.4	2.8
S1	1982	01	09	12	53	52.0	47.000	−66.600	7.0	5.7
S2	1982	03	31	21	02	20.0	47.010	−66.570	2.5	5.0
S3	1982	06	16	11	43	00.0	47.010	−66.970	7.0	4.7
	1982	01	17	13	33	55.7	46.98	−66.61	4.9
	1982	01	17	13	32	59.9	46.99	−66.63	2.4
	1982	01	18	19	34	49.2	47.01	−66.61	3.5
	1982	01	21	00	39	55.7	46.98	−66.60	4.8

Table 1.

Catalog of the located 68 aftershocks.

Note: the epicenters and origin times are those outputs from computer program (no shift-correction yet); the magnitudes m_N came from the NRCan database. The magnitude type for No. 11 and No. 22 is mb. Column H lists the focal depths determined using the depth–phase method. Columns er_n and er_e are the errors in km for the N–S (latitude) and E–W (longitude), respectively. These error values are in the output files. Dis stands for the distance between an aftershock and its closest neighbor event; modu stands for the modular of a vector formed using the errors of one aftershock and those of its closest neighbor event; perc stands for percentage (100 × modu/dis). No. 11 (mb 5.0) is indicated in Figure 2; No. 22 (mb 4.7) is indicated in Figure 1; No. 36 and No. 44 are used in the Section 6. The earthquake No. 69 is for a small aftershock, S1 for the mainshock, S2 and S3 for the two principal aftershocks [4]; and the last four are small aftershocks located by Saikia and Herrmann [6].

Figure 6.

Regional depth phase sPg modeling at KLN (distance 23.6 km) for the ME. The top synthetic vertical trace U_D/5.5 was generated using a depth of 5.5 km. Other synthetic traces were generated with a depth increment of 0.1 km. Trace KLN/SHZ is the observed vertical short period seismogram at KLN. The synthetic and the observed Pg are aligned. The time difference between sPg and Pg along trace U_D/5.9 and the time difference along the observed trace is approximately equal. Therefore, the modeled depth for the ME is 5.9 km.

After the focal depth for No. 36 was obtained, the SEISAN computer program was run at the newly obtained focal depth value. During the first trials of the epicentral location, the residuals between the arrival times of the observed and the calculated Pn phases at GGN were not small, so the P wave velocity value in the crustal model beneath the Moho was adjusted to reduce the residuals.

6.2 Location of the 68 aftershocks

The aftershocks with m_N ≥ 2.8 that occurred after KLN installation usually had clear onsets of Pg and Sg phases at KLN, Pg at EBN, and Pn at GGN. Of the 113 aftershocks for which the focal depths were determined [10], 68 satisfied the requirements for using the ME relocation method. Therefore, 68 aftershocks were located at the focal depths determined using depth phase sPg. The epicenters of the 68 aftershocks are plotted in Figure 2 and listed in Table 1.

6.3 Epicenter corrections

Since the station coverage was poor, the available arrival time readings were limited, and a 1-D crustal velocity model used, an epicentral shift relative to its true location was unavoidable. To obtain an epicentral distribution with absolute errors as small as possible, we performed an epicentral shift correction.

The star labeled S2 in Figure 2 marks the epicenter of the mb 5.0 aftershock, determined by Wetmiller et al. [4] using the centre of the small aftershocks they detected. It was assumed that the epicenter had smaller absolute errors compared to the epicenter for the same aftershock we obtained. The reason is that Wetmiller et al. [4] used the arrival times at portable stations which were less than 10 km from the earthquake sequence the absolute errors in the small aftershocks they detected were small. Accordingly, the absolute errors in the epicenter of the mb 5.0 they obtained were smaller. The differences between the two epicenters for the same mb 5.0 earthquake were subtracted from the epicenters of all the 68 aftershocks inTable 1. Figure 2 shows the corrected epicentral distribution.

6.4 Location of the mb 5.4 aftershock which did not have KLN record

Since the mb 5.4 aftershock was strong, it had clear Pg and Sg at EBN (Figure 4), clear Pn at GGN, LPQ , and GSQ , and it already had an accurate focal depth solution [10], so it can be located without a record at KLN. The aftershock mb 5.0 also had common phase arrival time reading as those of the mb 5.4, and had an accurate focal depth solution, so the mb 5.0 can also be located with the same precision as that of the mb 5.4. After the two epicenters were obtained using the arrival time readings at the above 4 stations, the epicenter of the mb 5.0 was corrected to that obtained by Wetmiller et al. [4]; accordingly, the epicenter of the mb 5.4 was moved the same amounts in latitude and longitude as did for the mb 5.0. The diamond symbol indicated with 5.4 in Figure 2 shows the corrected epicenter for the mb 5.4 aftershock.

6.5 Distribution features of the located hypocenters

After the epicenter shift correction, the mainshock is located within the southern part of the located aftershock cluster (Figure 2). Most aftershocks occurred within a 5 × 5 km² area, with the remaining ones scattering to the southwest. Overall, the aftershocks trend in a northeasterly direction (about Az 38°). This trend is close to the strike of one of the nodal planes for the mainshock obtained by the CMT group (202°) and by Choy et al. ([8]; 195°). Within the overall trend, the located aftershocks appear to form a pair of northeast trends separated by a gap region indicated by a dashed-line at Az 38°.

To observe more distribution features, the epicenters were divided into three groups: the left group, right group, and the bottom group (Figure 2). The hypocenters in the left and right groups were projected onto a vertical plane at Az 128° (Figure 7). The gap region indicated by a vertical dashed line in the figure separates the aftershocks into clearly two groups. The hypocenters on the left side were clustered together, and most were in a depth range of 4 to 6 km. In this left group only one aftershock with magnitude ≥5.0, of which the epicenter (the triangle) was determined by Choy et al. [8]. This epicenter and its focal depth have been corrected for comparison. The hypocenters on the right side were distributed from about depth 7 km to about 2.5 km. In this right group the mainshock and aftershocks mb 5.4 and mb 5.0 located. The triangles show the projections of hypocenters determined by Choy et al. [8].

Figure 7.

Projections of the shift-corrected hypocenters in the left and right groups onto a vertical plane, striking at Az 128° (NW–SE) indicated in Figure 2. The diamond shows the projection of the hypocenter for the mb 5.4 aftershock, determined in this article. Most of the aftershocks occurred at depths between 3 km and 6 km. Stars S1 and S2 represent the hypocenter projections of the mainshock and the mb 5.0 aftershock, respectively, obtained by Wetmiller et al. [4]. The triangles indicated with 5.7, 5.1, 5.4, and 5.0 show the projections of the mainshock and its three principal aftershocks, located by Choy et al. [8]. The epicenter of the mainshock was calibrated to that obtained by Wetmiller et al. [4]; the epicenters of the 3 principal aftershocks were moved with the same amount of the distance and direction as those of the mainshock, accordingly. The 4 focal depths were calibrated to those obtained by Ma and Motazedian [10].

We also projected all the corrected hypocenters onto a vertical plane at Az 38° (Figure 8). The aftershocks were separated into two groups by a gap region indicated by a vertical dashed line. The aftershocks at the left side of the vertical line are those in the bottom group in Figure 2. The aftershocks at the right side are those in the left and right group. They occupied a region of around 5 × 5 km². Comparing to the region the left and right groups occupied in Figure 2, it can be inferred that the source volume formed by the left and right group is about 5 × 5 × 5 km³.

Figure 8.

Projection of the located 68 hypocenters (the epicenters were shift-corrected) onto a vertical plane, striking at Az 38° (NE–SW), indicated by a dashed line with an arrow in Figure 2. Most of the aftershocks occurred at depths between 2 km and 6 km. Stars S1 and S2 represent the hypocenter projections of the mainshock and the mb 5.0 aftershock, obtained by Wetmiller et al., [4]. The triangles indicated with 5.7, 5.1, 5.4, and 5.0, show the hypocenter projections of the mainshock and its 3 principal aftershocks, relocated by Choy et al., [8]. The aftershocks were separated into two groups by a gap region indicated by a vertical dashed line. Most of the aftershocks occupy a region of around 5 × 5 km².

7. Time series and strength of the three earthquake groups

In Section 6 the spatial features of the located aftershocks were analyzed. In this section the time series and the strength of the three earthquake groups are very briefly analyzed. After the events were sorted and plotted in Figure 9, we found in the upper panel (a), the right group (see Figure 2) there are 23 events, 6 of them with magnitude ≥4.0; in (b), the left group, 26 events, 2 of them with magnitude ≥4.0; and in (c), the bottom group, 22 events, only 1 of them with magnitude ≥4.0. Rough calculations show that the energy released in the right group is more 10 times than that in the left group; more 200 times than that in the bottom group. A common feature is that after 7 years since the mainshock, aftershocks still occurred in all groups.

Figure 9.

Time series and strength of the 3 earthquake groups. Each vertical line segment stands for an earthquake; the length shows the magnitude; the position along horizontal axis shows the origin time (the origin time of the mainshock was set to 0). (a) In the right group (see Figure 2) there are 23 events, 6 of them with magnitude ≥ 4.0. (b) In the left group, 26 events, 2 of them with magnitude ≥ 4.0. (c) In the bottom group, 22 events, 1 of them with magnitude ≥ 4.0. All the origin times of the earthquakes were converted into Unix epoch times, and the origin time of the mainshock was set to 0. The time duration is from 1982-0109 12:53 52 to 1989–0109 12:53 52.

8. Error estimates in the relative locations of epicenters

When the azimuthal station coverage of an event is not good, the calculated location can move away from the true location because of the inaccuracy in the used velocity model [27]. Figure 1 shows that the station coverage was not good, and the shifts were unavoidable in the epicenters we obtained. In this section, we examine the uncertainties (errors) in the relative locations between one epicenter and another. If the errors are small, the patterns of the epicentral distribution are reliable.

Errors in locations could be caused by inaccuracies in the used crustal model. To examine if the errors in the relative locations between two epicenters in a group were small, a location test for the ME and one SE was performed. We chose two aftershocks (the ME, No. 36, and an SE, No. 44, in Table 1; both have clear onsets, Figure 5). Figure 10 shows the obtained epicenters of the two aftershocks using the same two data sets of arrival time readings but different crustal models and crustal thicknesses. The solid circle with mb 5.7 shows the location of the mainshock, determined by Wetmiller et al., [4]; the two solid squares show the locations of the ME and SE, obtained using our crustal model and a crustal thickness of 32 km. The two solid circles labeled 32, 34, 36, and 38 mark the epicenters of the same two aftershocks obtained using the GSC crustal model (Vp = 6.2 km/s) and crustal thicknesses of 32, 34, 36, and 38 km, respectively.

Figure 10.

The epicenters located using different crustal models for the same 2 aftershocks. The solid circle with mb 5.7 represents the epicenter of the mainshock, determined by Wetmiller et al., [4]; the two solid squares with ME and SE represent the located epicenters using our crustal model with thickness of 32 km. The 2 solid circles side by side show the same two aftershocks, located with the GSC crustal model (Vp = 6.2 km) and crustal thicknesses of 32, 34, 36, and 38 km, respectively. When the crustal thickness was changed, the 2 epicenters moved, but their relative positions were visually unchanged.

As shown in Figure 10, the epicenters of the two aftershocks mainly shifted southwards, and the relative locations of the two aftershocks were visually unchanged.

When the two epicenters were over plotted, only subtle changes in relative distances and relative azimuths between the two epicenters could be found. The absolute locations of two epicenters were determined using the calculated travel times, the same crustal model, and the same two sets of observed arrival times. When the parameters in the crustal model change, an increase or decrease in the calculated travel times causes the epicenter to move accordingly. Since the relative locations are mainly determined by the differences in the two sets of observed travel times (Figure 5), when the observed travel time readings were not changed, the relative distances and relative azimuths between the two epicenters could only have subtle changes, due to the change in used crustal model. Therefore, this test shows that errors in a crustal model cause systematic errors in the epicenters of an earthquake group, the errors in the relative locations in a group are very small.

Qualitatively speaking, two major types of errors -- errors in the crustal model and in the phase arrival time readings, cause the errors in the epicenters. The errors in the epicenters caused by arrival time reading errors may be roughly estimated. Along the top trace in Figure 5, the Sg – Pg time is $\delta$ t = 6.97–4.00 = 2.97 s. If the P-wave traveled to station KLN with Vp = 6.2 km/s, and Vp/Vs = 1.74, the distance between the station and the epicenter is $\Delta =\delta$ t $\times Vp/\left(1.74-1\right)=$ 24.88 km. Since the precision of arrival time readings is to 2 decimal places, the reading error in $\delta$t is ±0.02 s, and the error in the station distance is δΔ = ±0.02 $\times \mathit{Vp}/\left(1.74-1\right)$ = ±0.17 km. If the error in the crustal model causes a +3.5 km error in the latitude of an aftershock, the total error (caused by the error in the crustal model and the error in arrival time readings) is +3.5 ±0.17 = +3.67 or +3.33 km.

Since the same crustal model and arrival time readings of the common phases at the common stations were used to locate the aftershocks, the signs for the absolute errors in the output files should be same. For example, the epicenter of aftershock a is (46.977° ±3.2 km, −66.612 ±3.1 km; #36 in Table 1), and the epicenter of aftershock b is (46.982° ±3.3 km, −66.613 ±3.2 km; #37). For aftershock a, if we take (46.977° +3.2 km, −66.612-3.1 km), then b is (46.982° +3.3 km, −66.613-3.2 km). In other words, for all aftershocks, the same sign of the errors (+ or -) needs to be assigned because the same crustal model was used. The error caused by the crustal model dominates the total error in an epicenter in the output files of the location program.

The errors in the relative locations of two adjacent aftershocks may be mathematically estimated using the absolute errors in their epicenters. Assume that the epicenters of any two adjacent earthquakes a and b are (latitude_a, longitude_a) with errors (err_na, err_ea) and (latitude_b, longitude_b) with errors (err_nb, err_eb). If vectors A = A₀ + ∆A = (latitude_a, longitude_a) + (err_na, err_ea) and B = B₀ + ∆B = (latitude_b, longitude_b) + (err_nb, err_eb), then the vector difference is C = B - A = (B₀ – A₀) + (∆B – ∆A) = C₀ + ∆C. Then we obtain C₀ = (latitude_b - latitude_a, longitude_b - longitude_a) and ∆C = (err_nb – err_na, err_eb - err_ea). Vector C₀ shows the location of b relative to that of a_, while vector ∆C shows the errors in the location of b relative to a. Assume aftershock No. 8 in Table 1 is a and one of its adjacent aftershocks, No. 14, is b, the errors in the location of No. 14 relative to that of No. 8 are −0.1 km (2.5–2.6) in latitude and − 0.1 km (1.8–1.9) in longitude. No. 13 and No. 15 are adjacent, and the errors in their relative location are (−0.1 km, −0.1 km).

For any given aftershock in Table 1, its neighbor events can be found by comparing its distances to other aftershocks. For a given aftershock the distance from its closest neighbor event, the module of its ∆C with its closest neighbor, and the ratio of this module over this distance are listed in Table 1. The average of the modules is 0.183 km, which is much smaller than the gap indicated by the line with an arrow at Az 38° in Figure 2 or that indicated by the vertical line in Figure 7; therefore, the pattern of the obtained hypocenters is reliable.

9. Discussion and conclusion

More than 40 years ago, on 9 January 1982 in the Miramichi region of north-central New Brunswick, an earthquake with magnitude mb 5.7 occurred. Since the digital seismographs had not been widely deployed at that time, the source parameters of the mainshock and its aftershocks were not well determined. We analyzed the seismograms and found that at station KLN, there were very clear onsets of Pg- and Sg-waves, at EBN, clear onsets of Pg-waves, and at GGN, clear onsets of Pn phase for the larger aftershocks. We also unexpectedly found that the depth phase sPg was well developed and recorded at KLN. Once the velocity records were converted into displacement records, the onsets of the depth phase could be read correctly and accurately. The depth phase information can be used to determine focal depth accurately, and the Pg, Sg, and Pn arrival time readings at the three stations can be used to stably determine epicenters and the origin times at fixed focal depth using a conventional location method. To obtain a reliable epicentral distribution pattern, the uncertainties (errors) in the relative locations of the epicenters need to be small. To reduce errors in the relative locations, we used a simplified master-event method; specifically, we used the elite part in the master-event method. By using arrival time readings of the common four phases at the same three stations, the errors in the relative locations of adjacent aftershocks can be small.

In our study on the locations of larger aftershocks in the Miramichi sequence, we made a great effort to reduce the errors in the relative hypocenters. The errors in the relative locations between an epicenter and its closest neighbor event are listed in Table 1. Most of the values in the modu column are less than or equal to 0.3 km, which is smaller than the narrowest part of the gap indicated by the dashed vertical line in Figure 7. In another word if any epicenter in the left group is moved by 0.3 km in any direction, it cannot go into the right group. This implies that the pattern of two groups is reliable.

The left and right group phenomenon could be observed by arranging the waveform records at stations KLN and EBN. Figure 3 shows the vertical component seismograms recorded at KLN, generated by aftershocks along a line approximately running through KLN at about Az 128°. The top 5 records were generated by the aftershocks in the right group, while the bottom 5 records by aftershocks in the left group (see Figure 2). The T_Sg – T_Pg times along the top 5 records are shorter than those along the bottom 5 records; due to those epicenters to KLN are shorter. The time δt indicated along the bottom trace in Figure 3 corresponds to the spatial gap between the left and right group. Similarly Figure 4 shows vertical component seismograms recorded at EBN generated by 9 aftershocks, occurred along a line running through EBN at about Az 128°. The T_Sg – T_Pg times (about 16.54 s) along the traces 2, 3, 4, and 5 in the top panel (left group) are shorter than those (about 17.04 s) along the bottom 5 traces (right group). The time difference 0.5 s (17.04–16.54) also corresponds to the spatial gap between the left and right group in Figures 2 or 7.

In Figure 2 the located aftershocks were divided into three groups. The mainshock, its largest (mb 5.4), and the mb 5.0 aftershocks are in the right group. In the left group only one principal aftershock (the mb 5.1) occurred there. The total energy released in the right group is more 10 times than that released in the left group. The fault in the right group is larger than that in the left group. However, the located aftershock number is 26 in the left group, more than that (No. 23) in the right group (Figure 9). The left group looks like an earthquake swarm.

The aftershocks we located were in about 8 years when KLN station was operated. Since the mainshock occurrence more than 40 years has past, the aftershock activity in the source region still continues. The mystery related to the Miramichi earthquake sequence, such as why there are so many aftershocks, and several principal aftershocks followed the mainshock, is still waiting for being explored.

Based on the analyses of the located aftershocks in previous sections, the following can be concluded: (1) the major source volume was about 5 × 5 × 5 km³; (2) the focal depths ranged from about 2 km to 7 km; (3) two separate fault systems (the left and right group in Figure 2) were activated; the right one, activated by the mainshock and its two principal aftershocks, was large, and most energy was released there; the left one was small; (4) the trend in the aftershock epicenters was close to the northeast strike of the nodal plane for the mainshock; (5) the epicenter distribution trend was parallel to the trend of the Appalachian Mountain range (NE–SW); and (6) for more than 40 years the aftershocks have been occurring in the mainshock source region.

The procedure used to locate the Miramichi aftershocks has been successfully used for other earthquakes (e.g., [28, 29]). It is applicable for any earthquake sequence that has depth phase records.

The reliable epicentral distribution trends we obtained are crucial information for the seismic hazard assessment in the source region and its vicinity.

Acknowledgments

This research was supported by the Natural Sciences and Engineering Research Council of Canada under the Discovery Grant program. We gratefully acknowledge the book editor Dr. Walter Salazar for INTECHOPEN LIMITED, and Dr. John Adams at Natural Resources Canada for their constructive comments, suggestions, and text revisions. The waveform records were processed using SAC, redseed, and geotool programs. The location program HYPOCENTER in the SEISAN program package was used; Figure 1 was prepared using the software GMT [3].

Data and resources

The seismograms used in this study were collected from the Natural Resources Canada (NRCan) at earthquakescanada.nrcan.gc.ca (last accessed in November 2016). The SEISAN program tutorial by Havskov, Ottemöller, and Voss (2014) was downloaded from seisan.info/seisan-tutorial.pdf. The Global CMT Catalog is available at globalcmt.org/CMTsearch.html.