Indian Explosions of May 11, 1998 : An Analysis of Global Seismic Bodywave Magnitude Estimates

S. K. Sikka, Falguni Roy and G. J. Nair
High Pressure Physics Division, Bhabha Atomic Research Centre
Mumbai 400085, India
Seismic  waves  generated  by  the  Indian explosions of May 11, 1998, showed  large  variations  in  the  globally   estimated   body   wave magnitudes(mb). The estimated mb values were in general smaller at the stations  in  east  and  west directions with respect to these sources than in the north direction. Synthetic  seismograms  demonstrate  that such  pattern  is due to the cancellation and superposition of signals from these explosions separated in space by about 1  km.  In  view  of this,  the average mb estimates of the international data center (IDC) at Arlington, USA (mb=5.0) and the US geological Survey  (mb=5.2)  are smaller  than  the  true  mb  value.  After  taking  into  account the necessary corrections, a value of mb=5.39 is obtained  as  the  global average.  Gauribidanur  array  (GBA)  seismogram showed anomalous PP/P amplitude ratio. After correcting for this  anomaly,  an  estimate  of mb=5.4  is  obtained  at GBA. The revised mb estimate gives an average combined yield of 58 + 5 kt.
  At the Pokhran test site in Rajasthan, three nuclear explosions were detonated  by  India  on  11  May,  1998  at  1543  hours  IST . These explosions were triggered simultaneously and comprised a thermonuclear device, a fission device and a subkiloton device emplaced in spatially separated shafts(1). These explosions (hereafter  referred  as  POK-2) were  fully  contained  from radioactive point of view and the seismic body waves and  surface  waves  generated  by  these  explosions  were recorded  by  several  regional  and international seismic stations. A closer  examination  of  these  seismic  observations   revealed   the following  interesting  facts  : (1) A plot of body wave magnitude(mb) estimates of POK-2  from  51  global  stations,  as  reported  by  the International  Data  Center  (IDC),  Arlington,  USA, as a function of epicentral distance (Figure 1) together with that of the  Gauribidanur array,  India  (GBA)  shows large variations in the measured mb values (4.1 to 5.8) . It can be easily  deduced  from  the  figure  that  the maximum  to  minimum  P  wave  amplitude  ratio  for  POK-2 at a given distance is about 30 which is thrice as large  compared  to  that  for other  underground  explosions(2). (2) When the above global magnitude estimates are plotted as a function of event azimuth (figure 2)  ,  it is  seen  that  the  stations  between  340  deg.  and 20 deg. azimuth (distance ranging from 28.510 to 94.620) with  respect  to  the  POK-2 site  in  general showed higher mb values than the stations located in the  eastern  or  western  directions.  Adequate  data  from  southern direction  were  not  available as this region happened to be oceanic. This observed magnitude pattern is shown, below, to be consistent with the source geometry of POK-2. (3) The average mb value reported by the IDC is 5.0, while the US geological survey (USGS) reported the  mb  as 5.2  and  the surface wave magnitude (Ms) as 3.6. This Ms value is the same as that measured from the Bhabha Atomic Research  Centre?s  stand alone  seismic station at Jodhpur1. Figure 3 shows a plot of Ms versus mb for two  population  of  events  viz.  earthquakes  and  explosions together  with  those  of POK-2 in order to have a comparison of these estimates. It is seen that  the  USGS  estimate  falls  close  to  the explosion  population  on  this  plot whereas the IDC estimate lies in between the two populations indicating that the estimated  average  mb may  be  lower  than  the  true value. (4) On Gauribidanur array (GBA) seismogram (Figure 4), it is seen that the amplitude to period  ratio of the PP phase for POK-2 is about 1.5 times that of the P phase. This anomalous  behaviour, which was also observed for the Indian explosion of 1974( hereafter referred as POK-1) , indicates  that  the  P  phase from  Pokhran region undergoes much more attenuation as compared to PP phase while traveling to GBA.
The  first  three  of  these  observations  may  be  attributed to the modification of the body-waves due to the time  delays  introduced  by the  physical  separation  of the two large explosions of POK-2( their shafts were located 1km apart in  east-west  direction)  which  varied from  0.0s  to  about 0.125s. Thus, minimum time delay would be in the north and the south directions with respect  to  POK-2  while  maximum delay  would  be  in  the  east  and west directions. With the help of synthetic seismograms it will be demonstrated in the next section that the resultant amplitude of two explosions reduces from a maximum  true value  to  a minimum value as the delay is varied from 0.0s to 0.125s. However, time delays of this order do not  effect  the  Rayleigh  wave amplitudes  as  they  have  much larger periods in comparison to these delays(3). With the help of the observed  global  data  and  synthetic seismograms,  this  paper  aims to 1) highlight the effect of multiple explosions on global seismic  magnitude  estimates  and  2)  obtain  a realistic estimate of combined yield of these explosions.
Synthesis of teleseismic signals corresponding to  POK-2
The  synthesized  explosion waveform ,O(w), in the frequency domain is obtained as                    O(w)  =  S(w)D(w)M(w)R(w)I(w)                   (1) where  w  is  the  angular frequency , S(w) is the source function(4), D(w) is the  source  crust  function,  M(w)  is  the  mantle  transfer function(5)  ,  R(w)  is  the  receiver crust function and I(w) is the broad band seismograph response function. Fourier inverse transform of O(w)  gives  the  synthetic  explosion  seismogram   Y(t).   Composite seismogram, Z(t), corresponding to two explosions with amplitude ratio r is obtained as
   Z(t)  =   r H(t-tp) Y(t-tp) + H( t-tp-t1)  Y(t-tp-t1),          (2)
where  H(t)  is the Heaviside function, tp is the travel time of the P wave at a given epicentral distance and t1 is the time  delay  between the  two  P  arrivals  which  depends  primarily on the azimuth of the recording station and apparent phase velocity of the P wave.  For  the two  large  explosions  of  POK-2  with yield ratio 3:1 (determined by close-in ground shock measurements) , the value  of  r  is  2.33  (see equation 3).
     The  effect  of anelastic attenuation for a signal frequency f is given by exp(-pft*) where t* is the ratio of the travel  time  to  the average  Q  on the path from source to receiver. For paths from Nevada Test Site, USA, which lies in an orogenic belt to stations in  high  Q regions  t*  is  estimated(6)  to  be  0.35s - 0.45s. Further, beneath stable aseismic regions attenuation is expected to be lower  than  the orogenic  regions(7)  .  In  view of the above as well as to match the observed signal periods at various stations, we have used t* values of 0.4s  and  0.5s  for  synthesizing  the  broad  band  seismograms   at teleseismic  distances  (3000-10000  kms).  Figures 5 and 6 show the synthesized  waveforms  corresponding to t*=0.4s and 0.5s respectively for various time delays t1. For the minimum delay, i.e.  t1=0.0s,  the waveforms interfere constructively ( this is expected in the north and the  south  directions)  while  for  the  maximum  delay,  t1  =0.125s (obtained by taking an average value of 12  km/s  for  apparent  phase velocity  and  accounting  for  some  scatter  in  the  delays  due to anisotropy and  heterogeneities  in  the  source  region),  they  show cancellation  effect.  The  ratio  of  A/T  (where  A and T are signal amplitude and period respectively) for maximum and minimum time delays for t* value of 0.4s is obtained as 1:2 (see  Figure  5)  whereas  the same  for  t*=0.5s  is  1:1.76(see  Figure 6). These reductions in the values of A/T correspond to deviations of -0.30 and -0.25 respectively from the true mb value. Deviation in mb corresponding to a t* value of 0.3s is found as -0.46. Thus some of the worldwide mb  values  may  be underestimated up to these amounts.
Anomalous amplitude ratio of P and PP phases at Gauribidanur array
The  medium  aperture  seismic  array at Gauribidanur, GBA, has twenty short period seismometers arranged along two perpendicular  arms.  The digital  data  at  GBA are sampled at a rate of 20 per second for each sensor. More details and design characteristics of GBA  are  available elsewhere(8). Short period seismograms at GBA pertaining to POK-2 show PP/P  amplitude  ratio  close to 1.5. Similar feature was observed for POK-1 also. Figure 4 shows two array beams corresponding to POK-2 data at GBA, one tuned for the P phase with slowness 13.1 s/deg. (trace  b) and  the  other  one  for  the PP phase with a slowness of 14.0 s/deg. (trace a) for an azimuth of 332 deg. with respect to GBA. A comparison of P/PP amplitude ratio in trace a with that in trace b reveals that P amplitude has reduced in comparison to the PP  amplitude  in  trace  a thus  confirming  the  second  prominent  phase as PP which is also in agreement with the arrival time of PP at GBA from POK-2.
  A  comparison  of  amplitude  correction term B(D,h) (where D is the epicentral distance and h is the depth of  the  source)  developed  by Gutenberg(10)  for  P  waves  with  that  for  PP for various distance ranges(11) reveals that the reduction in amplitudes of P and PP phases at a given distance due to the geometrical spreading  are  comparable. However,  the  PP phase undergoes an additional attenuation due to the reflection at the free surface. The epicentral distance  of  GBA  from the  POK-2  site  is 14.4 deg. The PP phase should get attenuated by a factor of 2 as the reflection  coefficient(9)  for  this  distance  is around 0.5. Thus, after accounting for the reflection coefficient, the net  PP  to  P  ratio should be around 0.5. However, the observed PP/P amplitude ratio is about 1.5. From this observation,  it  is  inferred that  the  P  wave from POK-2 has undergone much more attenuation than what is accounted by B(14.4,0.0) value. In view of  this,  it  becomes essential  to  include  a  correction  term in the body wave magnitude relation to account for the additional attenuation of P wave amplitude by a factor of 3 ( obtained  by  comparing  observed  and  theoretical values  of  P/PP amplitude ratio) for this path. With the inclusion of this term , the mb estimate at GBA is obtained as 5.4. It may be noted that GBA is situated at an azimuth of 157 deg.  with  respect  to  the test  site  therefore the signals from POK-2 interferes constructively giving the true mb value of 5.4.
Revised magnitude estimate from global data
As pointed out in the introduction, the IDC through its reviewed event bulletin  announced  the  average  mb  for POK-2 as 5.0 whereas the mb estimate of USGS was 5.2. As shown above, with the help  of  synthetic seismograms we have demonstrated that the mb of POK-2 as seen globally will  vary  between a true value, when observed from the perpendicular direction with respect to the line joining the two  large  explosions, and  a minimum value when observations are made parallel to this line. However, there could be some effect on this value due  to  anisotropy. Signals from two explosive sources when superposed give the true mb as the  maximum  value whereas any time delay between them will result in partial cancellation of the amplitudes leading to lower mb values. Due to this a large number of global stations  could  have  underestimated the true mb.
   In  view  of the above, to obtain a true estimate of mb from global data it will be reasonable to  consider  only  those  data  which  are available  from  the perpendicular directions with respect to the line joining the sources  because  only  such  data  will  be  composed  of superposed  signals (t1 ~ 0.0s) . Using data from 12 such stations, 10 in the north between azimuth 340 deg. and 20 deg. and two in the south including GBA(Table 1), the average mb is obtained as  5.36.  However, it  may  be noted that MAW station showed a very small signal to noise ratio (SNR) for POK-2. At such low SNR values  ,  estimation  of  true signal  amplitudes  is  likely  to be less accurate. An estimate of mb using eleven stations data  after  excluding  that  of  MAW  gives  an average  value  of  mb  as 5.39 which is in good agreement with the mb estimate of GBA discussed in the previous section.  On  the  (Ms,  mb) plot (Figure 3), the point (3.6,5.4) , corresponding to our estimates, falls well within the explosion population.
  It  may  be  interesting  to  point  out  here that for the Pakistan nuclear explosion of May 30, 1998, the stations between 340  deg.  and 20  deg.  azimuth  with  respect  to POK-2 having similar azimuths and epicentral distances for Pakistan test site , gave an average mb value of 4.46 which is less by 0.14 compared to the  average  global  mb=4.6 reported by the IDC. This confirms that the higher average mb estimate obtained  from these stations for POK-2 was less dependent on the path effects and was mainly due to the source geometry of POK-2. Estimation of the yield using mb-yield regression relations
The  details  of  the  mb and medium properties of 17 peaceful nuclear explosions conducted at various locations in USA with precisely  known yields  from  radiochemical estimates are listed(12) in Table 2. Using the well known form of the relation between mb and yield viz.
                  mb   =   C1 + C2 log(Y) ,                        (3) the data points of table 2 are fitted to estimate the constants C1 and C2. Y is the yield in kilotons. The regression constants C1 and C2 are estimated  for average values of mb and the spreads in these constants are estimated by fitting maximum and minimum mb values  to  the  known yields.  The values of C1 and C2 obtained by least square fitting this data are 3.8+0.15 and 0.77+0.05 respectively. It is  well  known  that the  value  of  C1 depends on the site geology and its value is ~4 for hard rock formations like granite and is ~3.3 for alluvial basins(13). From literature it is seen that constant C2 varies(14) between 1.0 and 0.6 . For Nevada test site and Novaya Zemlya its values are  0.81  and 0.75 respectively(15). In order to estimate the site specific constant C1  for  Pokhran  , the known yield of 13 kt for POK-1 (estimated from the surface wave data at Quetta(16) , an yield  exponent  C2  =  0.77+ 0.05  obtained  from  the  above data fit and an average mb of 4.9 for POK-1 as reported by the International Seismological Centre,  UK,  are used in equation 3. The value of C1 for Pokhran region found from this analysis  is  4.04  + 0.04. The estimated yield of POK-2 for an mb=5.4 using the values of C1=4.04 and C2=0.77 turns out as 58kt (see  Figure 7)  with  a  spread  of  +  5kt. This is close to the value of 60 kt , estimated from a preliminary analysis1.  It  may  be  noted  that  the fitted  constants  for  US tests( see broken line in Figure 7) give an yield of 120kt for mb=5.4.
Conclusions
Large  variations  in  the estimated global magnitudes ranging between 4.1 to 5.8 were observed for POK-2. Further,  it  was  seen  that  the amplitudes  recorded  at  stations  having  similar distances from the explosion site varied up to a factor of about 30 which was  thrice  as larger than what is normally expected. Stations located in the eastern and  western direction from the POK-2 site, in general, recorded lower magnitudes compared to those situated in the north. With the  help  of synthetic  seismograms it is shown that the observed global mb pattern was consistent with the source geometry. In view of  such  a  pattern, averaging of the global mb estimates will result in an underestimation of  the true mb value. Therefore, in order to obtain a realistic value of mb, estimates from eleven stations comprising ten  from  the  north between azimuth 340 deg. and 20 deg. and GBA mb value of 5.4 have been used.  The  average  mb  obtained  by  this process is 5.39. Using the revised mb value of 5.4, the combined yield of POK-2  is  obtained  as 58+5 kt, which is very close to the preliminary estimate(1).
References
1.Sikka, S. K.  and  Kakodkar, Anil , BARC News Letter, BARC,1998, 172,1-4.
2.Dahlman, O. and  Israelson, H., Monitoring Underground Nuclear Explosions,    Elseviers, Amsterdam, 1977.
3.Kolar, O. C.,  and Pruvost, N. L., Nature, 1975, 253, 242-245.
4.Seggern, V. D. and Blandford, R., Geophys. J. R. astr. Soc, 1972, 31, 83-97.
5.Futterman, W. I.,  J. Geophys. Res., 1962,  67, 5279-5291.
6.Der. Z., Mcelfresh T.,  Wagner, R. and Burnetti, J., Bull. Seism. Soc. Am., 1985,
   75, 379-390. 7. Marshall, P. D., Springer, D. L. and Rodean, H. C., Geophys. J. R. astr. Soc.,
   1979, 57,609-638. 8.Varghese ,T. G., Roy , F.,  Rao, B. S. S., Suryavanshi, M. P. and Bharthur, R. N.,     Phys. Earth Planet. Inter., 1979, 18, 87-94.
9. Aki, K.  and Richards, P. G., Quantitative Seismology Theory and Methods,W. H.     Freeman and Co., San Francisco,1980.
10.Gutenberg, B, Bull. Seism. Soc. Am., 1945. 35, 57-69.
11.Richter, C. F., Elementary Seismology, W. H. Freeman and Co.,San     Francisco,1958.
12.Higgins, G. H., Proc. Panel  Vienna, IAEA, 1970, 111-121
13.Murphy, J. R., In:Identification of Seismic Sources- Earthquake or Explosion     (eds. E. S. Husebye and S. Mykkeltveit), 1981, 201-   205.
14. Springer, D. L. and Hannon,  W. J., Bull. Seismol. Soc. Am., 1973, 63, 477-500.
15.Chyba, C. F., van der Vink, G. E. and Hennet, C. B., Geophys. Res. Lett., 1998,      25 , 191-194.
16.Nair, G. J., Some Seismic Results of Rajasthan Explosion of 18 May, 1974,     AG224, Procurement Executive, MOD, UK.
17.Marshall, P. D. and Basham, P. W., Geophys. J. R. astr. Soc., 1972,28,431-458.
Acknowledgements.  We  are  grateful  to  Dr.  R.  Chidambaram  and Dr. Anil  Kakodkar for constant encouragement and support. Thanks are also due to  Dr.  A.  R.  Banghar  for  critically reading the manuscript. We thank Shri R. N.  Bharthur and other colleagues at  Gauribidanur  seismic  array  station  for  providing transcripts of digital data of GBA seismograms used in the present  study.