VALIDATION OF SPATIAL AUTOCORRELATIONMETHOD WITH L-SHAPE ARRAY IN JYOSO CITY,
JAPAN
Prithvi Lal SHRESTHAMEE09182

Supervisor:Toshiaki YOKOI

ABSTRACTThis study is aimed to validate the efficiency of an L-shape array for SPAC method using microtremorin estimating the shear wave velocity (Vs) structure. The experiment for validation was conducted in
the Toyota Community Baseball Ground, Jyoso City, Ibaraki Prefecture, Japan in March 2009 with anequilateral triangle array with side length of 40m and in June 2010 with an equilateral triangle arraywith side length of 50m, together with an L-shape array of the similar size. Multichannel Analysis ofSurface Waves (MASW) was also performed simultaneously in June 2010. In the same lot PS-logging
data are available from nearby IBRH10 station of the KIK-NET (NIED, Japan) that shows soft
sediment of about 20m thick with Vs of 110 m/s in the geological column of the site.
The comparison of the determined phase velocity and that calculated from PS logging datashows close matching of two sets of curves separately. One is between PS logging and the triangle
array (40m) and the other is between the triangle array (50m) and the L-shape array (50m). Former two
are of the almost same place whereas other two arrays are deployed about 200m away from the otherset. Some discrepancy between two sets is shown. This seems due to lateral variation of underground
velocity structure which is consistent with the result of MASW.
Based on the results of analysis we can say that an L-shape array can be applied to estimate
shear wave velocity for shallow depth so it can be layout in urban areas to determine phase velocity
information from microtremor. Therefore it may be feasible to apply it in the Kathmandu valley, Nepal
that is based upon the soft soil with high possibility of liquefaction or earthquake hazard.
Keywords: SPAC, PS logging data, Equilateral Triangle array, L-shape array.

1. INTRODUCTION
Nepal Himalaya lies in the active seismic belt. Seismicity in the Himalaya is the consequence ofunder-thrusting of the Indian plate towards the north lying Tibetan plate. Nepal has sufferedearthquake disasters through its history due to its high seismicity and highly vulnerable construction
practices, therefore one of the most earthquake disaster prone countries in the world.
1.1. Purpose of my study
Since we are expecting a great earthquake disaster in near future in our country Nepal, seismic hazard
assessment is required. Estimation of the amplification factor of the ground, i.e., microzonation on thebasis of the shallow Vs structure is a basic step of seismic hazard assessment. The microtremor arraymeasurement will be applied in Nepal in future, because this SPAC method seems to be reliable, easyNational Seismological Centre, Department of Mines & Geology, Kathmandu, Nepal.
Chief Research Scientist, International Institute of Seismology and Earthquake Engineering, BRI, Japan.

1to handle, comparatively affordable and do not cause any environmental problems, thus suitable as atool for seismic microzonation and earthquake disaster mitigation.
Conventional circular or equilateral triangle arrays, however, are sometimes difficult to
deploy in crowded urban areas where only L-shape or linear arrays layout can suit to road pattern. The
purpose of this study is to understand the effectiveness, limitations and advantages of L-shape arrayfor the SPAC method by applying and verifying the applicability and accuracy by comparison withstandard equilateral triangle array in Jyoso city, Ibaraki, Japan.

2. METHODOLOGY
2.1. SPAC method
This is a successful method to determine the phase velocity information from surface waves containedin microtremor (e.g., Aki 1957; Okada 2003). The SPAC coefficient for distance r between two
stations at the angular frequency ω provides the information about the phase velocity of the
propagating waves in the array. This is obtained by the azimuthal average of coherency between
microtremor records observed at two stations.
2 r 
Re[ S cx (; r ,  )] = 1
d ,

c
(

)
2



S cc (;0,0) S xx ( , r ,  )0ρ (  , r ) = Jo (1)

where, Jo( ) is the Bessel function of the first kind with zero order, c(ω) is the phase velocity (ω=2πf).
2.2. New interpretation
Shiraishi et al. (2006) proposed the following formula using a mathematical relationship betweencos(xcosθ) and the Bessel function of the first kind,n r 
 r  L  2  1 J 2 n  l cos 2n l ,Re  pq   J 0 n 1
 c    c   l 1(2)

where, Re (γpq) is equal to the integrand of the third member of Eq.(1), namely the coherency, L is thenumber of wave sources, λl is rate of the contribution of the l-th wave source to the power spectra at
the observation point,θl is azimuth of the l-th wave source, Jm( ) is m-th order Bessel function of first kind. Here J0 (ωr/c(ω))
is already included in the integrand. The application of L-shape array relies on this formula. L-shapearray has error of the 4th and higher order Bessel functions as follows.
L, r01  Re z 0,1  Re z 0,22

  J  r   J  r 2  cos4  OJ  r  .
010 c L c   l l 1
01

4

l




01

 8 c 

 

(3)

Equilateral triangle array averages out the terms of J2( ) and J4( ), whereas L-shape arraycancels out the term of J2( ) only. Therefore, L-shape array seems weaker against undesirable azimuth
dependent noise contained in microtremor. From Eq.(3), it is expected that the average over two pairsthat form the short sides can reduce the influence of J2( ) term drastically, whereas the oblique sides
follow Eq.(2). Therefore it is better to eliminate the oblique sides from the analysis.

22.3. Procedure of analysis
Triangle and L-shape arrays are set up
to get the microtremor data (Figure 1).
First SPAC coefficient is obtained
after filtering and resampling and thenscreening is done. Second, thedispersion curve of Rayleigh wave isdetermined
andfinally
velocity structure model is estimatedby the heuristic search method using
the very fast simulated annealing
method (VFSA) combined with thedownhill simplex method (DHSM)(VFSA-DHSM, Yokoi 2005).

Figure 1. Schematic figures of A) equilateraltriangle array and B) L-shape array.

3. DATA ACQUISITION
3.1. Observation site
The experiment was conducted in the
Toyota Community Baseball Ground,Jyoso City, Ibaraki Prefecture, Japan in
March 2009 and June 2010. IBRH10of KIK-NET (NIED, Japan) is locatedat the south-east corner of the same lot
where PS-logging data are available.
MASW (Hayashi and Suzuki 2004)
was also conducted this year. Softsediment as thick as about 20 m withvelocity of S-wave (Vs) 110 m/s isobserved in the geological column at
IBRH10. There are a national road
R294 about 300 m west and a
Figure 2. Location map of measurement site (Yokoi
prefectural road 24 about 100 m southand Hayashi 2009).
and both of them have heavy traffic all
day long (Figure 2).
The baseball ground itself, however, was almost quiet during the measurement (Yokoi and
Hayashi 2009).
3.2. Array deployment and instruments
Equilateral triangle array with side length of 40m with 7 sensors was deployed in March 2009 nearby
IBRH10 whereas 50 m of equilateral triangle with 10 sensors and L-shape array with 11 sensors were
deployed in June 2010 about 200m away from IBRH10. All sensors are vertical component. MASW
with total length of 105m was conducted in the same place.

34. ANALYSIS AND RESULTS
4.1. Preprocessing
Multiplexing and re-sampling are done applying digital anti-aliasing filter (Saito 1978). The screening
is conducted next with two steps. The data are divided into time blocks with 512 samples. The
consecutive time blocks are overlapped each other by 50% of their duration.
First step: If peak is greater than “ajudge” times of RMS amplitude then this time block is
not used in analysis. This is a countermeasure against impulsive noise due to traffic, i.e. vehicles
passing near by sensors.
Second step: If RMS amplitude in a time block deviates more than another given constant“a_sgm” times the standard deviation from the average, this time block is not used in analysis, where
the average and the standard deviation are calculated over the all time blocks that survived in theabove mentioned screening step 1. This is a countermeasure against outliers.
In this study ajudge=4 and a_sgm=2 are used for the screening of the obtained data.
4.2. SPAC coefficient calculation
Re-sampled and screened time block files are used to calculate SPAC coefficient that is an azimuthal
average of the coherency between microtremor records at two stations. The initial frequency range ofanalysis is set from 0.1Hz to 10Hz. Band width of Parzen window for smoothing power and crossspectra is set at 0.5Hz.
Figure 3 A) shows the SPAC coefficient of station pairs that compose the 50m triangle array
and B) shows SPAC coefficient with its standard deviation of the 50m triangle array (inter-station
distance 25m). It is clear that for the shorter distances the first positive lobe of SPAC coefficientshows higher values and has wider range of frequencies. It is also clearly seen that SPAC coefficient isdecreasing at low frequency side but theoretically it should reach 1 at the frequency 0.0Hz.
This variation is due to the available frequency range of seismometer. For Figure 3 thenatural frequency of the used seismometers is 2Hz.
1. 0

12. 5m

0. 8

14. 4m

0. 6

25. 0m

0. 4

28. 9m

0. 2

43. 3m

0. 0

50. 0m1.2
1.0
SPAC Co efficien t

SPAC Co efficien t

1. 2

-0. 2
-0. 4

25. 0m

0.2
0.0
-0.2
-0.4
-0.6

-0. 6

A 0

0.8
0.6
0.4

2

4

6

8

B

10

Fre q u e n cy(Hz)

0

2

4

6

8

10

F re q u e n cy(Hz)Figure 3. SPAC coefficients: A) Those of the station pairs that compose the
50m triangular, B) of the 50m triangular array (inter-station distance 25m).
4.3. Determination of dispersion curve
The dispersion curve of Rayleigh wave i.e., the phase velocity depending upon the frequency isdetermined from SPAC coefficient. First, SPAC coefficient ρ(ωr) is converted to the value of kr byapplying the following fifth order polynomial equation that approximates the inverse function of J0(kr)from kr = 0.0 to the first trough of J0(kr). If x = J0(y),y  6.0803x 5  9.2477 x 4  3.9322 x 3  0.1815x 2  1.7079 x  2.4121

4(4)

The first maximum of c(ω) from the low frequency side is recognized as a lower limit of the available
frequency range for the respective inter-station distance. The maximum from this lower limit to thehighest one is recognized as the high frequency limit of the available frequency range. Then thesevalues again are averaged and converted to c(ω) = rω/(kr). The weight coefficient used at averaging is
the reciprocal of the variance of SPAC coefficients at the respective inter-station distance and thefrequency.
4.4. Estimation of velocity structure by Heuristic search
A heuristic search method is conducted to obtain the optimum underground structure by fitting thetheoretical phase velocity of Rayleigh wave to the observed dispersion curve. Five layer models are
introduced with its search range. The used method is the downhill simplex method (DHSM, e.g., Press
et al. 2002) combined with the very fast simulated annealing (VFSA, Ingber 1989). Hereafter, thecombined methods is called the DHSM-VFSA. The optimum, namely, the fastest schedule for the
inversion of underground velocity structure from the dispersion curve of Rayleigh waves is used withthe parameters t0=1.0, a=0.6, and c=1.3 as given by Yokoi (2005).

5. DISCUSSION
5.1. Comparison of dispersion curves

Dispersion Curve
400

40mTriangle
50mTriangle

350

L-shapePS-logging

300
Phase Velocity(m/s)

Figure 4 shows a comparison of phasevelocity over frequency between 40m
triangle, 50m triangle and 50m L-shapearrays with the curve of PS logging data.
This shows close matching of two sets ofcurves separately, one is between PS
logging and 40m triangle and the other is
50m triangle and L-shape. PS logging and
40m triangle array are from the same place
where other two arrays are also from thesame place but about 200m away from theother set. A discrepancy is seen between
these two sets. This seems to be due tolateral variation of underground velocity
structure as MASW result implies.

250
200
150
100
50
0
0

1

2

4

5

Figure 4. Comparison of the dispersion curves withPS logging data.

5.2. Comparison of velocity structure

Vs StructureVs (m/s)
0
0

100

200

300

400

PS-logging
40m Triangle
50m L-shape
50m Triangle

10
Depth (m)Figure 5 shows the comparison of Vs
structure models determined by the data ofthe three arrays with the PS logging data. Vp
is fixed to 1500m/s for the 1st to 4th layerand 1956m/s for the underlying half space.
Density is calculated from Vp then fixed to
1.9g/cm3 and 2.13g/cm3 respectively(Ludwig et al. 1970). As seismometer’snatural frequency is 2.0Hz, the lower limitof the frequency range for analysis is set at
2.0Hz and the upper one at 5.0Hz byconsidering on the S/N of microtremor. The
value of misfit is as small as 3.0, 4.0 and 4.2

3
Frequency(Hz)

20

30

5

Figure 5. Comparison of velocity structureswith PS-logging datam/s for the triangular array (40m), triangular array (50m) and L-shape array respectively.
From Figure 5 it is clear that for the 20m depth Vs of all arrays are similar to the logging datai.e., around 110m/s to 150m/s that implies soft soil. Vs structures of the triangular array (50m) andL-shape array coincide to each other whereas 40m triangle shows a clear deviation from them andlower Vs than PS logging data. This shows that the lateral variation of underground velocity structure
implied by MASW is detected also using SPAC method.

6. CONCLUSION
In this study I have conducted the analysis of the data of microtremor array observation held in the
Toyota Community Baseball Ground, Jyoso city, Ibaraki prefecture, Japan for shallow depth using
L-shape and equilateral triangle array with the SPAC method for Rayleigh wave in comparison withPS logging data and MASW. The dispersion curve and Vs Structures determined by the data ofL-shape array coincide well to those of equilateral triangular array (50m). These deviate together fromVs structure of PS logging data and the dispersion curve calculated from PS logging data respectively.
This discrepancy, however, is due to lateral variation of Vs structure as the result of the triangular array(40m) shows and that of MASW implies. Based on the results of analysis I conclude that L-shapearray can be applied to estimate Vs structure for shallow depth so it can be layout in urban areas todetermine phase velocity information from surface wave in microtremor.

ACKNOWLEDGEMENT
I would like to express my sincere gratitude to Dr. Koichi Hayashi, OYO Corporation for his help
during the field work.

REFERENCES
Aki, K., 1957, Bulletin of Earthquake Research Institute, 35, 415-457.
Hayashi, K., and Suzuki, H., 2004, Exploration Geophysics, 35, 7-13.
Ingber, L., 1989, Mathematical and Computer Modeling, 12, 967-973.
Ludwig, W. J. et al., 1972, The Sea, Vol. 4, Wiley-Interscience, New York, 53-84
Okada, H., 2003, The Microtremor Survey Method. Society of exploration geophysicists.
Press, W. H., et al., 2002, Numerical recipes. Cambridge University Press.
Saito, M., 1978, BUTURI- TANSA, 31, 112-135 (in Japanese with English abstract).
Shiraishi, H., et al., 2006, Geophys. Res Let, Vol. 33, L18307.
Yamanaka, H., 2004, Proc. of 13th World Conference on Earthquake Engineering, Paper 1161.
Yokoi, T., 2005, Programme and Abstract, Seismological Society of Japan Fall Meeting, B049.
Yokoi, T., and Hayashi, K., 2009, Proc. of 9th Int. WS on Seismic Microzoning and Risk Reduction,
Cuernavaca, Mexico.

6