Abstract
One of the most important challenges in seismic evaluation of structures comes from wide variability of possible motions that a specific site may experience. This variability arises from the intrinsic randomness of seismic motions (aleatory uncertainty) and lack of knowledge about the physical process behind an earthquake (epistemic uncertainty). The epistemic uncertainty in estimating ground motion can be properly characterized by providing an understanding of the physical aspects of earthquake source models. In this paper, a combination of the theoreticalbased Green’s function method and hybrid slip source model is used to reduce the source uncertainty. The procedure is implemented on the 1994 Northridge earthquake as a case study. For this purpose, the method is validated against recordings stations through demonstrating a good agreement between the elastic response spectra corresponding to the simulated and recorded data, confirming the reliability of the procedure. A simple additive weighting algorithm is used to select the best fit amongst simulated waveforms and optimal source model.
1. Introduction
Due to the increasing number of applications of nonlinear analysis techniques aimed at assessing the structural response and estimating the damage from future seismic events, a credible and complete characterization of ground motion is required for earthquake engineeringanalysis. Engineers usually use a suite of recorded accelerograms from past earthquakes for timeseries analysis. They select these accelerograms according to certain criteria. However the use of recorded data for prediction of future ground motions produce some problems. For example in some regions that seismicity is low and therefore the number of recorded motions is not sufficient to consider the variability of earthquake ground motion, engineers have to use records from other regions and scale them to match with the elastic design spectrum prescribed by the adopted seismic code. When the structure is located on a particular soil condition, scaling the spectral amplitude of a signal recorded on a rock site does not provide information on the long duration actually generated [1]. Similarly, earthquake ground motions recorded in regions near to faults are different from those observed in the farfault regions. In nearfault regions, faulting mechanism, direction of rupture propagation (e.g., forward directivity), and the static deformation of the ground surface (e.g., flingstep) significantly affect ground motions. These specifications also may have serious implications for flexible structures [2]. Scaling the records and matching them with the target elastic design spectrum does not consider the specifications of nearfield earthquake ground motion.
As an alternative approach, engineers can simulate realistic synthetic seismograms for engineering purposes rather than using records from past earthquakes [35]. Applying simulation techniques, engineers can directly use the synthetic time series, so scaling problems is eliminated. However, these simulation techniques need the definition of a large number of input physical parameters. The accurate values of input parameters rarely can be fixed a priori, and this fact limits the use of simulation proceduresin engineering practices. Because of the uncertainty in input parameters, we have to generate a large number of strong motion scenarios to manage the input values in possible ranges. So, a number of ground motion parameters for each site can be generated and statistically analyzed to provide some insights on the probability distributions of groundmotion values and related statistical estimators. Therefore expanding knowledge about each input parameter can reduce overall uncertainty in a seismic scenario.
The objective of this article is to present the procedure to combine the hybrid slip source model of earthquake scenario source with theoreticalbased Green’s function simulation method.We then include the effect of uncertainties of different earthquake source parameters in simulation of nearfield events.
2. Theoretical method for simulating nearfield strong motions
During the last decade, researchers have proposed various numerical approaches for simulating ground shaking. The three wellknown procedures for simulation of ground motion are stochastic procedure based on $\omega $Square source model, empirical Green’s function, and theoreticalbased Green’s function methods [6, 7]. Point source [8] and finite fault [9] techniques are $\omega $Squared based procedures in which the Fourier’s spectrum of the motion at a site is broken into contributions from earthquake source, path, site and type of motion. The difference between these two procedures is that in the former the causative source of earthquake is modeled as a point source, while in the latter the fault plane is discretized into several equal rectangular elements, each of which is treated as point source. Empiricalbased Green’s function method [10, 11] uses the recordings of small events to reproduce large ground motions. Theoreticalbased Green’s function method [12] computes the displacement components of ground motion based on Green’s function analytically. We can also find use of above methods in the form of hybrid to reproduce broadband seismogram in engineering problems [13, 14].
However, there are limited number of consistent procedures to model the specifications of near field earthquakes such as fling step and forward directivity. There is a modified form of Green’s function solution proposed by Hisada and Bielak [15] that considers fling step and forward directivity effects at nearfault sites. In this approach the singularities of the original Green’s functions are eliminated by subtracting and adding the static Green's function in the following form:
where, ${U}_{k}$ is the $k$th component of displacement in the Cartesian coordinate system, $Y$ and $X$ are an observation point (i.e. the point location where the simulated ground motion characteristics are monitored) and a source point on the fault plane, respectively, $S$ is the fault plane, ${D}_{i}$ is the $i$th component of the fault slip, and ${T}_{ik}^{D}$ and ${T}_{ik}^{S}$ are the dynamic and static traction Green’s functions of layered halfspaces. The simulation of nearfault ground motions through Eq. (1) is much more efficiently than the original representation theorem [16] which is represented by the total traction Green’s function alone. The original representation theorem requires a lot of CPU time to numerically evaluate the nearsingularities of the dynamic Green's function, which show sharp peaks centered at the source point. The first integral of Eq. (1), in which the static Green's function is subtracted from the dynamic function, models the directivity effect. This integral can also eliminate the nearsingularities, because the static function includes all the sources of the singularities of the dynamic function. The second integral involves the static Green’s function and models the fling step effect. The near singularities appear only in this integral. The evaluation of this integral necessitates distribution of a set of integration points in the neighborhood of the observation point. However, this integration needs to be performed only once, as the values of the static functions remain invariant for all frequencies. Also in the near field ground motion simulation, in which the source is close to the observation point, it’s needed to compute a lot of Green’s functions. So, the wavenumber integrations of Green’s functions require a lot of CPU time, because the integrands do not converge to zero when the wavenumbers increased [17]. To cope with this problem, Hisada and Bielak [15] used the following wavenumber integrations to evaluate the integral of the first integrand of Eq. (1), since the dynamic Green’s functions are already subtracted by the static Green’s functions:
where ${t}_{ik}^{D}$ and ${t}_{ik}^{S}$ are the dynamic and static terms of the integrands of the traction Green’s functions, respectively. Since the dynamic integrand converges exactly to the static integrand with increasing wavenumber, the integrand of Eq. (2) attenuates rapidly. A schematic diagram illustrating the simulation process is shown in Fig. 1. References [15, 17] are recommended for researchers who interested in more details. This extended kinematic source method is implemented in the GRFLT12S code by Hisada and Bielak [15]. A number of researchers [18, 19] has recently used this approach and has demonstrated the efficiency of this procedure in simulation of nearfield strong ground motion.
Fig. 1Simulation of near field ground motion
3. Hybrid slip model considering uncertainty in source parameters
Advances in understanding the source of past earthquakes led to efforts to estimate future earthquake source properties, and a large amount of work has been done in recent years to estimate the distribution of slip on the fault plane during earthquakes. However, accurate determination of heterogeneous slip on the fault plane before an earthquake is not possibleat the present time. The asperity model and the $k$square slip modelare two wellknown probabilistic approaches for prediction of heterogeneous slip on fault plane. The asperity model first quantified by Somerville et al. [6] using slip distributions of 15 events, and developed by Wang [20] and Wang and Tao [21] by using slip distributions of 29 events. Different relationships between moment magnitude of earthquake and source/asperity parameters is developed based on this model. Herrero and Bernard [22] proposed the $k$square slip model, in which the slip spectrum decays as ${k}^{2}$slip in the wave numbers larger than the corner wave number. Later, Gallovic and Brokesove [23] developed twodimensional $k$square slip model.
Recently, the hybrid slip model is presented using combination the asperity and $k$square models [24, 25]. An example of the hybrid slip model containing two asperities is shown in Fig. 2. To generate the fault model, the moment magnitude of earthquake scenario is at first determined from seismic hazard analysis [26]. Based on this moment magnitude, the fault model parameters (such as the fault dimensions, the number of asperities, the dimensions of asperities, the location of asperities on the fault plane, the average slip on the fault plane, and the average slip on the asperities) are determined by statistical empiricalrelations. The fault plane is divided into a number of subfaults and the deterministic slips are assigned to each subfault (Fig 2(a)). The slip on the fault plane is transferred from spatial domain into wave number domain using twodimensional Fourier transform. Then ksquare slip model is used to generate random slip for wave numbers higher than corner wave number (Eq. (3)):
where ${k}_{cx}$ and ${k}_{cy}$ are corner wave numbers along strike and down dip, $\mathrm{\Phi}({k}_{x},{k}_{y})$ is phase spectrum, $\stackrel{}{D}$ is the average slip on the fault, $L$ and $W$ are the length and width of fault. Finally, deterministic and random part are superposed in wave number domain and transferred back into spatial domainwith adding noises to deterministicpart of slip distribution (Fig. 2(b)).
Previous studies show that in the case of using the hybrid slip source model and stochastic finite fault simulation technique to generate synthetic accelerograms, the estimated acceleration response spectra is consistent with the recorded response spectra at periods lower than 1.0 sec, but for intermediate and high periods the results are overestimated [27]. This discrepancy may result from the essence of the stochastic finite fault models that tend to overestimate the structural response in long period [28]. To overcome this problem, we use the theoreticalbased Green’s function method that models the low frequency motions more precisely. So, we evaluate the efficacy of hybrid slip model for the prediction of motions in longer periods.
Fig. 2An example of the hybrid slip model with two asperities: a) the deterministic part of the slip distribution, b) the resulting hybrid slip
a)
b)
4. Simulation of the 1994 Northridge earthquake
In this section we assume the 1994 Northridge earthquake (${M}_{W}=$6.7) as a future seismic event, considering that no previous studies (as waveform inversion analysis) are available to constrain the source model. We generate the plausible sources by hybrid slip procedure and simulate several possible shaking scenarios. Table 1 shows the relations and their deviations of fault model parameters and location of hypocenter within fault plane.We use the empirical relations proposed by Wang [20] and Wang and Tao [21] in which the uncertainties in both global and local source parameters are considered in fault model simulation.Assuming one up to four different asperities on fault plane [6] and varying hypocenter location within fault plane, one hundred fault models are generated using hybrid slip source model.
Table 1Relation and standard deviation of fault model parameters and location of hypocenter [20, 21]
Fault model parameters  Relation  $\sigma $ 
Rupture area $S$ / km^{2}  $\mathrm{l}\mathrm{o}\mathrm{g}S={M}_{W}$ 4.05  0.29 
Fault length $L$ / km  $\mathrm{l}\mathrm{o}\mathrm{g}L=0.5{M}_{W}$1.9  0.18 
Fault width $W$ / km  $W=S/L$  – 
Average slip on the fault $\stackrel{}{D}$ / cm  $\mathrm{l}\mathrm{o}\mathrm{g}\stackrel{}{D}=0.5{M}_{W}$1.35  –0.29 
Corner wave number along strike ${k}_{cx}$  $\mathrm{l}\mathrm{o}\mathrm{g}{k}_{cx}=$1.89$$0.5${M}_{W}$  0.18 
Corner wave number down dip ${k}_{cy}$  $\mathrm{l}\mathrm{o}\mathrm{g}{k}_{cy}=$2.18$$0.5${M}_{W}$  0.16 
The whole area of asperities ${S}_{all}$ / km^{2}  $\mathrm{l}\mathrm{o}\mathrm{g}{S}_{all}=\mathrm{l}\mathrm{o}\mathrm{g}S$ 0.67  0.16 
Area of the max. asperity ${S}_{m}$ / km^{2}  $\mathrm{l}\mathrm{o}\mathrm{g}{S}_{m}=\mathrm{l}\mathrm{o}\mathrm{g}S$ 0.83  0.18 
Length of the max. asperity ${L}_{m}$ / km  $\mathrm{l}\mathrm{o}\mathrm{g}{L}_{m}=\mathrm{l}\mathrm{o}\mathrm{g}L$0.44  0.15 
Coordinate of the max. asperity center along strike ${X}_{m}$ / km  $\mathrm{l}\mathrm{o}\mathrm{g}{X}_{m}=\mathrm{l}\mathrm{o}\mathrm{g}L$0.53  0.18 
Coordinate of the max. asperity center down dip ${Y}_{m}$ / km  $\mathrm{l}\mathrm{o}\mathrm{g}{Y}_{m}=\mathrm{l}\mathrm{o}\mathrm{g}W$0.30  0.2 
Average slip on the max. asperity $\stackrel{}{D}$_{m} / cm  $\mathrm{l}\mathrm{o}\mathrm{g}{\stackrel{}{D}}_{m}=\mathrm{l}\mathrm{o}\mathrm{g}\stackrel{}{D}+$0.34  0.07 
Average slip on the other asperities ${\stackrel{}{D}}_{0}$ / cm  $\mathrm{l}\mathrm{o}\mathrm{g}{\stackrel{}{D}}_{0}=\mathrm{l}\mathrm{o}\mathrm{g}\stackrel{}{D}+$0.28  0.06 
Location of hypocenter within fault plane  Relation  $\sigma $ 
Relative coordinate of rupture start point along strike ${X}_{s}$ / km  $\mathrm{l}\mathrm{o}\mathrm{g}{X}_{s}=\mathrm{l}\mathrm{o}\mathrm{g}L$0.37  0.24 
Relative coordinate of rupture start point down dip ${Y}_{s}$ / km  $\mathrm{l}\mathrm{o}\mathrm{g}{Y}_{s}=\mathrm{l}\mathrm{o}\mathrm{g}W$0.09  0.17 
Rake angle is assumed variable over the fault plane and is allowed to vary within the range of 55 to 145 degrees [29]. Slipvelocity function is the functional form of slipvelocity with time and its details influence ground motions. In this paper, three different source time functions namely rectangular, triangle and exponential slipvelocity functions [30] are used with multiple time window for each subfault on the simulated ground motion. Rise time is the time period between rupture initiation time and the time of peak moment release rate. This parameteris highly variable on the fault plane. We use the empirical relationship between the value of rise time and the seismic moment of target earthquake, ${M}_{0}$, proposed by Somerville et al. [6]:
This equation yields a rise time of 1.01 second for Northridge earthquake with ${M}_{w}=$ 6.7. So, we selected rise time of 1.01±1 second with uniform distribution to consider the uncertainty of this parameter. The rupture velocity is defined as the speed at which rupture front propagates on the fault plane and controls the signal duration and contributes to the directivity effect. We assumed the rupture velocity in the narrow range of 0.65–0.85 of the shearwave velocity (${V}_{S}$) which is inferred with the recorded strong motion data in earthquakes around the world [6]. Two crustal velocity models [29] were employed for rock and soil stations (Table 2). For ${V}_{s}<$1500 m/s and ${V}_{s}>$1500 m/s, the Swave qualityfactor is assumed as ${Q}_{s}=$0.02${V}_{s}$ and ${Q}_{s}=$0.1${V}_{s}$, respectively. As well as it is assumed that the Pwave qualityfactor is ${Q}_{p}=$1.5${Q}_{s}$ [31].
Table 2Velocity and density model for rock and soil stations
Depth (km)  ${V}_{p}$ (km/sec)  ${V}_{s}$ (km/sec)  Density (g/cm^{3})  
Soil  Rock  Soil  Rock  Soil  Rock  
0.0  0.8  1.9  0.3  1.0  1.7  2.1 
0.1  1.2  1.9  0.5  1.0  1.8  2.1 
0.3  1.9  1.9  1  1.0  2.1  2.1 
0.5  4  4  2  2  2.4  2.4 
1.5  5.5  5.5  3.2  3.2  2.7  2.7 
4  6.3  6.3  3.6  3.6  2.8  2.8 
27  6.8  6.8  3.9  3.9  2.9  2.9 
40  7.8  7.8  4.5  4.5  3.3  3.3 
Acceleration time histories during the 1994 Northridge earthquake are numerically simulated by combining the Hisada and Bielak procedure [15] with the hybrid slip source model. We generated the synthesized wave forms up to 2.0 Hz, since this limit of frequency is sufficient to reflect the directivity and fling step effects. Two different rock stations of Santa Susana (SSU) and Simi Valley (U55), and a soil station of Newhall Los Angeles county fire department (NHL) in nearfield were used to compare the simulated ground motions with the observed ones.
Fig. 3Comparison the simulated with recorded acceleration response spectra at three nearfield stations. The gray lines show the acceleration response spectra of 100 simulated ground motions, the black solid lines show the median of simulated acceleration response spectra and the dashed lines indicate the acceleration response spectra of recorded motions in the 1994 Northridge earthquake
One hundred ground motions were synthesized by using the 100 fault models mentioned above. We discretized the fault planes to 1 km×1 km patches. In each simulation run, the source parameters, namely, the location of hypocenter, rupture velocity, rise time and rake angle were randomly generated from specified distributions and ranges. So, the synthetic database contains 100 earthquake scenarios, all of which are probable. It should be noted that in each simulation run, two horizontal and one vertical component of ground vibration were calculated at each observation point. At each of the three stations, the corresponding acceleration response spectra of the simulated ground motions were worked out with a damping constant of 5 percent. Fig. 3 demonstrates the comparison of the acceleration response spectra corresponding to the simulated and recorded data at three different stations. In this figure, a large scatter compared to the median exists at each station, so the median spectrum actually represents a general level of estimated ground motion.The median acceleration response spectrum is generally consistent with the recorded response spectra for the three stations at periods longer than 1.0 sec. Because of the generation of the synthesized waveforms in narrow band frequencies (0 $\le f\le $ 2 Hz), the mismatch between response spectra extracted by the observed and the synthetic accelerograms is visible. As mentioned before, the results of previous studies show good agreement between recorded ground motions and the simulated ones in frequency domain for periods less than 1.0 sec. Here, in fact we show the efficiency of hybrid slip source model for estimation of nearfault ground motions whose often contain strong coherent dynamic long period pulses and permanent ground displacements and the dynamic motions are dominated by a large long period pulse of motion.
Different source scenarios are examined to find the best fitting between simulated and observed motions recorded at SSU, U55 and NHL stations for different orthogonal directions. Since the bestfitted motions at these stations are from different scenarios of the 100 source scenarios, a systematic search is performed to find the optimum source parameters that provide the best fitting between observed and simulated lowfrequency ground motions. We estimate the optimal values of source parameters based on simple additive weighting (SAW) method [32]. Because of its simplicity, SAW is the most popular method inmultiple attribute decision making. In SAW technique, the best alternative, ${A}^{*}$, can be derived by the following equations:
where ${u}_{i}\left(x\right)$ denotes the utility of $i$th alternative, ${w}_{j}$ denotes the weight of the $j$th criterion, ${r}_{ij}=(1/{x}_{ij})/(1/{x}_{j}^{*})$ are normalized values of decision matrix elements. Using this method, we can derive the best alternative of the problem response surface.
In order to find the optimum source parameters using SAW method, we select peak ground acceleration and acceleration response spectrum as attributes and assign equal weights to the attributes. The optimal fault model obtained via SAW algorithm is shown in Fig. 4. Also the optimal rise time and rupture velocity are 1.43 sec and 2960 m/s respectively which agree with the waveform inversion results [29]. The simulated acceleration time histories from the optimal fault model at the three stations, and their corresponding response spectra are shown in Fig. 5. This figure also provides a comparison in time and frequency domain at NHL, SSU, U55 stations between the simulated motions from the optimal source scenario and the recorded motions in the 1994 Northridge earthquake. Since the simulated data are depicted within the frequency range of 0.0 to 2.0 Hz, recorded acceleration time histories is filtered in this frequency range in order to compare simulated and recorded data in the time domain. Good approximate match of the observed and the simulated waveform, in time and frequency domain at three stations, confirms that the hybrid slip model is an applicable procedure to model seismogenic sources.
Fig. 4Optimal finite fault model for the 1994 Northridge earthquake
Fig. 5Comparison between the simulated motions from the optimal fault model and the observed motions in the 1994 Northridge earthquake at: a) NHL, b) SSU, c) U55 stations
a)
b)
c)
5. Sensitivity of ground motion to source parameters
In this section we perform a sensitivity study to determine which source parameters mostly affect the ground motion. Sensitivity analysis is the study of how the uncertainty in the output of a model can be apportioned to different sources of uncertainty in its inputs [33]. Sensitivity analysis is used to increase understanding of the relationships between input variables and output in a simulation model and to know the relative contribution of each input variable to uncertainty in output of the model. So, the variables that do not contribute much to uncertainty in model output can be fixed at their best estimates, rather than treated as random, in order to simplifying the analysis problem. Also, the variables that cause significant uncertainty in the output can then be the focus of attention to understand them better and reduce uncertainty. There are a large number of approaches to performing a sensitivity analysis, such as scatter plots, regressionbased methods, variancebased methods, etc. We used the local sensitivity analysis approach [34] in which the sensitivity of the model output to the model inputs has been examined by varying selected input oneatatime around its expected values and observing its effect on the desired model output. In this approach, one parameter is varied a time while the other parameters are kept constant at their nominal values, and this is repeated for each of the other inputs in the same way.
The sensitivity of ground motion to the source parameters has been studied by several authors [1, 23]. However, only few studies attempt to assess how the uncertainties in the different source parameters (especially the fault model and location of hypocenter within fault plane) propagate to nearfield ground motion.
Earthquake intensity measure is an attribute of the strength of an earthquake ground motion and is used to predict earthquake damage and the seismic responses of structures. Studies have shown that the common intensity measures (such as peak ground acceleration) cannot provide an adequate description of the damaging effects of nearfault ground motions which arestrongly influenced by forward directivity and fling step phenomena. Boore [35] proposed the maximum direction (i.e., 100th percentile) value of spectral accelerations computed over all rotation angles for two orthogonal horizontal ground motion components, Sarot100, as intensity measure for nearfield earthquake ground motions. The effect of uncertainties in source parameters on the variability of Sarot100 is investigated in this section.
Fig. 6Variaition of maximum value of spectral response for individually randomized source parameters at: a) 3 km and b) 100 km distances
a)
b)
The impact of source parameters on ground response variability is assessed using Hisada and Bielak [15] model, examining the effect on standard deviation of the maximum direction value of spectral accelerations caused by alternatively varying and fixing different parameters. Fig. 6(a) and 6(b) present fluctuation plots of standard deviation of Sarot100 for distances of 3 and 100 km. For generating each spectrum, one of the source parameters is varied while holding the remaining parameters fixed to base case values. These figures show that fault model mostly contribute to variability of ground motion at near and farfield ground motions. The location of hypocenter within fault plane has strong impact on ground motion variability near fault, but its impact on response variabilitydeclinesas faulttosite distance increases. Other source parameters (i.e. rake angle, rise time, rupture velocity and slip velocity function) do not contribute significantly on the ground response compared to fault model and hypocenter location, and the sensitivity of them depends on the structural periods of interest. The results show that the increasing knowledge about the source model can significantly reduce the uncertainty in prediction of future earthquake ground motion.
6. Conclusions
In this paper, the theoreticalbased Green’s function method of synthesized ground motion is combined with the hybrid slip source model and simple additive weighting method to predict nearfield ground motion. Also, we introduce uncertainties of global and local source parameters in the proposed procedure. The proposed procedure is applied to simulate the strong motion of three different nearfield stations in the 1994 Northridge earthquake. Good agreement between the synthesized and recordings, in terms of response spectra and timeseries, conﬁrms the reliability and ability of these models.
In the hybrid broadband ground motion simulation method, low and high frequency ground motions are simulated by different seismological models, and are then added in the time domain to produce a broadband ground motion. So, utilization of hybrid broadband ground motion simulation method and itscombination with hybrid slip source model can be useful in estimation of near fault ground motion for engineering practices and conﬁrms the reliability of over models in extrapolating procedure. Local sensitivity analysis is performed to investigate the importance of each source parameter. Results show that fault model and location of hypocenter within fault plane are the main important parameters on variability of ground motion.
References

Ameri G., Gallovič F., Pacor F., Emolo A. Uncertainties in strong groundmotion prediction with finitefault synthetic seismograms: an application to the 1984 M 5.7 Gubbio, central Italy, earthquake. Bulletin of the Seismological Society of America, Vol. 99, Issue 2A, 2009, p. 647663.

Kalkan E., Kunnath S. K. Effects of fling step and forward directivity on seismic response of buildings. Earthquake Spectra, Vol. 22, Issue 2, 2006, p. 367390.

Nicknam A., Abbasnia R., Bozorgnasab M., Eslamian Y. Synthesizing strong motion using empirical Green's function and genetic algorithm approach. Journal of Earthquake Engineering, Vol. 14, Issue 4, 2010, p. 512526.

Rezaeian S., Der Kiureghian A. A stochastic ground motion model with separable temporal and spectral nonstationarities. Earthquake engineering and structural dynamics, Vol. 37, Issue 13, 2008, p. 15651584.

Yamamoto Y., Baker J. W. Stochastic model for earthquake ground motion using wavelet packets. Bulletin of the Seismological Society of America, Vol. 103, Issue 6, 2013, p. 30443056.

Somerville P., Irikura K., Graves R., Sawada S., Wald D., Abrahamson N., et al. Characterizing crustal earthquake slip models for the prediction of strong ground motion. Seismological Research Letters, Vol. 70, Issue 1, 1999, p. 5980.

Kurahashi S., Irikura K. Source model for generating strong ground motions during the 2011 off the Pacific coast of Tohoku Earthquake. Earth, Planets and Space, Vol. 63, Issue 7, 2011, p. 571576.

Boore D. M. Stochastic simulation of highfrequency ground motions based on seismological models of the radiated spectra. Bulletin of the Seismological Society of America, Vol. 73, Issue 6A, 1983, p. 18651894.

Beresnev I. A., Atkinson G. M. Source parameters of earthquakes in eastern and western North America based on finitefault modeling. Bulletin of the Seismological Society of America, Vol. 92, Issue 2, 2002, p. 695710.

Hutchings L. J. Modeling nearsource earthquake ground motion with empirical Green’s functions. Doctoral dissertation, State University of New York at Binghamton, 1987.

Wennerberg L. Stochastic summation of empirical Green’s functions. Bulletin of the Seismological Society of America, Vol. 80, Issue 6A, 1990, p. 14181432.

Bouchon M. A simple method to calculate Green's functions for elastic layered media. Bulletin of the Seismological Society of America, Vol. 71, Issue 4, 1981, p. 959971.

Nicknam A., Eslamian Y. A hybrid method for simulating nearsource, broadband seismograms: Application to the 2003 Bam earthquake (Mw 6.5). Tectonophysics, Vol. 487, Issue 1, 2010, p. 4658.

Martin Mai P., Beroza G. C. A hybrid method for calculating nearsource, broadband seismograms: Application to strong motion prediction. Physics of the Earth and Planetary Interiors, Vol. 137, Issue 1, 2003, p. 183199.

Hisada Y., Bielak J. A theoretical method for computing nearfault ground motions in layered halfspaces considering static offset due to surface faulting, with a physical interpretation of fling step and rupture directivity. Bulletin of the Seismological Society of America, Vol. 90, Issue 2, 2003, p. 387400.

Aki K., Richards P. G. Quantitative seismology: Theory and methods. WH Freeman and Co, 1980.

Hisada Y. An efficient method for computing Green's functions for a layered halfspace with sources and receivers at close depths. Bulletin of the Seismological Society of America, Vol. 84, 1993, p. 14561472.

Hisada Y., Bielak J. Effects of Sedimentary Layers on Directivity Pulse and Fling Step. In Proceeding of the 13th World Conference on Earthquake Engineering, 2004.

Ghayamghamian M. R., Hisada Y. Nearfault strong motion complexity of the 2003 Bam earthquake (Iran) and lowfrequency ground motion simulation. Geophysical Journal International, Vol. 170, Issue 2, 2007, p. 679686.

Wang H. Finite fault source model for predicting nearfield strong ground motion. Doctoral dissertation, Institute of engineering mechanics, china earthquake administration, 2004, (in Chinese).

Wang H., Tao X. Relationships between moment magnitude and fault parameters: theoretical and semiempirical relationships. Earthquake Engineering and Engineering Vibration, Vol. 2, Issue 2, 2003, p. 201211.

Herrero A., Bernard P. A kinematic selfsimilar rupture process for earthquakes. Bulletin of the Seismological Society of America, Vol. 84, Issue 4, 1994, p. 12161228.

Gallovič F., Brokešová J. On strong ground motion synthesis with ksquare slip distributions. Journal of Seismology, Vol. 8, Issue 2, 2004, p. 211224.

Wang H., Tao X. Hybrid slip model for predicting nearfield strong ground motion. Annual Meeting: Networking of Young Earthquake Engineering Researchers and Professionals, Honolulu, Hawaii, 2004.

Gallovič F., Brokešová J. Hybrid ksquared source model for strong ground motion simulations: Introduction. Physics of the Earth and Planetary Interiors, Vol. 160, Issue 1, 2007, p. 3450.

Yazdani A., Nicknam A., Khanzadi M., Motaghed S. An artificial statistical method to estimate seismicity parameter from incomplete earthquake catalogs, a case study in metropolitan Tehran, Iran. Scientia Iranica, In Press, 2014.

Sun X., Tao X., Tang A., Lu J. Hybrid slip model for nearfield ground motion estimation based on uncertainty of source parameters. Transactions of Tianjin University, Vol. 16, 2010, p. 6167.

Motazedian D., Atkinson G. M. Stochastic finitefault modeling based on a dynamic corner frequency. Bulletin of the Seismological Society of America, Vol. 95, Issue 3, 2005, p. 9951010.

Wald D. J., Heaton T. H., Hudnut K. W. The slip history of the 1994 Northridge, California, earthquake determined from strongmotion, teleseismic, GPS, and leveling data. Bulletin of the Seismological Society of America, Vol. 86, Issue 1B, 1996, p. 4970.

Guatteri M., Mai P. M., Beroza G. C., Boatwright J. Strong groundmotionprediction from stochasticdynamic source models. Bulletin of the Seismological Society of America, Vol. 93, 2003, p. 301313.

Olsen K. B., Day S. M., Breadley C. R. Estimation of Q far longperiod (>2s) waves in the Los Angles. Bulletin of the Seismological Society of America, Vol. 38, 2003, p. 93627.

Tzeng G. H., Huang, J. J. Multiple Attribute Decision Making: Methods and Applications. CRC Press, 2011.

Saltelli A., Ratto M., Andres T., Campolongo F., Cariboni J., Gatelli D., Saisana M., Tarantola S. Global Sensitivity Analysis. The Primer, John Wiley & Sons, 2008.

Saltelli A., Chan K., Scott E. M. Sensitivity analysis. Vol. 134, New York, Wiley, 2000, p. 475.

Boore D. M. Orientationindependent, nongeometricmean measures of seismic intensity from two horizontal components of motion. Bulletin of the Seismological Society of America, Vol. 100, 2010, p. 18301835.
About this article
The authors would like to thank anonymous reviewers for comments which helped to improve the manuscript.