Active tectonics

The active tectonics of Iran is mainly affected by the northward convergence of the Arabian and Eurasian plates. The rate of this convergence relative to Eurasia is estimated at around ~25 mm/yr at longitude 58°E (^{Vernant et al., 2004}). Northward moving of Central Iran with respect to western Afghanistan at a rate of 6-13 mm/yr, obtained from GPS study (^{Walpersdorf et al., 2014}), produces N-S, right-lateral, shear across eastern Iran. This shear is accommodated on large, active N-NW trending, right-lateral strike-slip faults which slice eastern Iran (e.g., ^{Berberian, 1981}; ^{Walker and Jackson, 2004}; ^{Meyer and Le Dortz, 2007}; ^{Walker et al., 2009}; ^{Allen et al., 2011}). ^{Walpersdorf et al., (2014)} showed that these faults get additional NS shortening by rotating counterclockwise in the horizontal plane, at a rate of 1 ± 0.4°/Ma.

Main faults

Due to the arid and sparsely vegetated environment of southeastern Iran, the Quaternary activity of the faults and their related structures are well exposed. The faults in this region are either oblique-slip or strike-slip. The oblique-slip faults have both senses of motion i.e. dip-slip and strike-slip. The activity of the strike-slip fault systems has developed horsetail splay at the tips of the faults, releasing and restraining bends, which are accompanied by structures such as pull-apart basin and flower structure (Rashidi et al., 2017).

Major strike-slip fault systems in the study area (Fig. 1c) can be categorized based on their strike into N-S (Nayband, Gowk, and Sabzevaran), NW-SE (Jorjafk, Kuhbanan, Ravar, Mahan, and Rayen), and E-W (Lalehzar) directions. The NW-SE and E-W faults link the right-lateral N-S faults with other active faults in central Iran and distribute the right-lateral shear in this region (Rashidi et al., 2017).

Nayband fault system

The 300 km long Nayband right-lateral fault system with NS trending, extends along the western margin of the Lut block. Based on the main structural discontinuities, five right-stepping en-echelon pattern segments have been identified along this fault (^{Foroutan et al., 2014}). The 40 km long West Nayband fault is considered as a branch of the Nayband fault system that has right-lateral and revers movement (Rashidi et al., 2017). Despite showing a clear expression in the geomorphology along most of its length the Nayband fault is the only part of the western Lut border that has been seismically quiescent during the last few millennia (e.g. ^{Ambraseys and Melville, 1982}). The slip rate of 1.8 ± 0.7 mm/yr is estimated for Nayband fault.

Sabzevaran fault system

The N-trending Sabzevaran fault is a strike-slip fault with a small but significant reverse component that has an almost vertical fault plane (i.e. about 80°W) (^{Regard et al., 2004}). The Sabzevaran fault consists of several en-echelon segments with a right-step arrangement. Despite the evidence of high-level tectonic activity and the high slip rate, 5.7 ± 1.7 mm/yr (^{Regard et al., 2005}), there is no known destructive earthquake associated with the Sabzevaran fault. There are some minor faults at the northern tip of the Sabzevaran fault that are considered as its splays; such as the Sardueiyeh fault (around 100 km long fault that cuts Quaternary units in a reverse mechanism and still no earthquake-related deformation has found along it), Lalehzar fault (an active 85 km long right-lateral the strike-slip fault with some recorded earthquakes e.g. 1923/09/22 and 2010/07/31), Dalfard fault (a right-lateral strike-slip fault with reverse component), etc. It seems that shearing and thrusting of the recent rock units and alluvial fans along these faults are caused by the right-lateral strike-slip motion along the Sabzevaran fault (Rashidi et al., 2017).

Gowk fault system

The Gowk fault is an NNW right-lateral strike-slip fault that extends for 160-200 km along strike (^{Berberian et al., 2001}). It is a fault system that consists of segments in an en-echelon pattern and is one part of the Sabzevaran-Gowk-Nayband system of strike-slip faults (^{Walker et al., 2010}) which bounds the western edge of the Dasht-e-Lut desert (^{Berberian et al., 2001}). Pull apart basins with faulted boundaries (e.g. chaharfarsakh) have been generated along the Gowk fault system. By having a proven record of generating destructive earthquakes, (e.g. 1981.06.11 Mw 6.6; 1981.07.28 Mw 7.1; 1989.11.20 Mw 5.8; 1998.03.14 Mw 5.4 and 1998.11.18 Mw 6.6), this fault is one of the most active faults of Iran (e.g. ^{Berberian et al., 2001}). The current slip rate on the Gowk fault is estimated at 4.2 ± 0.7 mm/yr (^{Walpersdorf et al., 2014}).

Shahdad thrust and fold system

The Shahdad thrust and fold system is a curved NW trending along the East of the Gowk fault. Since the Gowk fault orientation is oblique to the NS right-lateral shear across eastern Iran, the overall motion is partitioned into strike-slip faulting along the Gowk fault and shortening on the Shahdad thrust and fold system (^{Berberian et al., 2001}). Unlike the Gowk fault that is associated with earthquakes since at least 1877, the Shahdad fault has no historical seismic record for the past millennia (^{Berberian et al., 2001}).

Kuhbanan fault system

The Kuhbanan fault is an active 300-km long, NW-SE-trending, right-lateral strike-slip fault, which consists of several en-echelon segments (^{Berberian, 2014}). The evidence shows that the sense of motion on the Kuhbanan fault may have changed from reverse to right-lateral slip through time (e.g. ^{Berberian, 2014}). The Kuhbanan fault currently slips laterally at 3.6 ± 1.3 mm/yr (^{Walpersdorf et al., 2014}). Although this fault has a long record of earthquakes (e.g. 1875, 1933 and 1977, ^{Ambraseys and Melville, 1982}), the largest earthquakes in the current cluster may not have yet occurred, due to the absence of a large-magnitude earthquake during the last 160 years of recorded seismic history along the Kuhbanan fault (^{Berberian, 2005}). Also, the GPS measurements imply that the fault is presently active and possibly slipping faster than over the Holocene and geological times. ^{Berberian (2005)} believes that all the faults in the southern tip of the kuhbanan fault (e.g. Dahuyiyeh, Chatrud, Bazargan, Dehbala, and East Kerman), are splay intermountain, cross-reverse faults which are developed at almost a right angle. Slip that is transferred from the Gowk fault to the Kuhbanan fault through the restraining bend between them, produces a wide zone of mountainous. Occurring of an NW-SE trend right-lateral strike-slip faulting within this zone (Dehu and Dehran faults) forms a flower structure in cross-section.

Rafsanjan fault

The Rafsanjan fault is an NW-SE-trending, right-lateral strike-slip fault ~200 km long. This fault is one of the major faults of southeast Iran, which shows a clear expression in the geomorphology along most of its length (e.g. ^{Walker, 2006}). Maybe it transfers an amount of the Sabzevaran fault right-lateral slip to the Anar fault. The overall NS right-lateral shearing across the region is accommodated by almost pure strike-slip on the Rafsanjan fault in the high, mountainous, regions and shortening component migrates to a series of parallel thrust faults in the plains north of the main Rafsanjan fault trace (e.g., ^{Walker, 2006} and ^{Walker et al., 2010}). An average slip rate of ~0.4 mm/yr is estimated for this fault (^{Fattahi et al., 2011}).

Jorjafk fault

The ~200-km-long Jorjafk fault with an NW-SE trending consists of five active main en-echelon segments (Rashidi et al., 2017). Although there is no record of historical seismicity, this fault shows clear signs of activity in the late Quaternary (^{Walker et al., 2010}). The overall motion on this fault seems to be oblique reverse and right-lateral. The dip-slip component migrates onto reverse faults in the adjacent lowlands whereas the strike-slip movement continues to occur along the main fault in the high land. (^{Walker et al., 2010}).

Ravar Fault

The Ravar fault is a right-lateral strike-slip with 150 km long. The Ravar fault consists of two main left stepping en- echelon pattern segments. Folds with different geometries can be seen throughout the Ravar fault, showing the presence of a reverse component. Flower structures were formed in the convergence zone between Lakarkuh and Ravar faults (at their southern tips). The river offsets are recognizable along it (Rashidi et al., 2017). Based on the horizontal slip rate of 0.54 mm/yr, the recurrence interval of an Mw 6.7 earthquake on the Ravar fault is estimated to be about 1400 years (^{Shafiei Bafti and Shahpasanzadeh, 2010}).

Rayen fault

Rayen fault is an NW-SE trending fault that shows a reverse component of motion in addition to its dominant strike-slip mechanism (Rashidi et al., 2017). There are a series of the NW-SE-trending folds along this fault that show signs of activity such as the truncation of the northern fold limb and the discrete northward-facing scarps developed in alluvial fans crossing the northern fold margin. This indicates the presence of underlying thrust faults that caused active folding in this region (^{Walker, 2006}).

Lakar Kuh fault

The N-S trending Lakar Kuh is a right-lateral strike-slip fault with a small reverse component. Its total length is reported to be about 130 km. Horsetail splays that formed at the northern and southern tips of this fault show a reverse mechanism with a small right-lateral strike-slip component. It is believed that the southern splay is the location of three events on 1 and 12 December 2017 with MN magnitudes of 6.0-6.1 (^{Nemati et al., 2020}).

Mahan fault-related folds

The NW-SE trending folds along the northern margin of the Kuh-e-Jupar which are situated close to the town of Mahan and extend northwards to beyond the village of Jupar, are considered as the Mahan fault-related folds. The geomorphic evidence such as the asymmetric pattern of river incision and height changes across the folds shows active folding in Mahan. In other words, these folds are underlain by active blind thrust faults that do not reach the Earth's surface (^{Walker, 2006} and ^{Walker et al., 2010}). A historical earthquake on 2nd January 1934 occurred in the district of Mahan which damaged villages along the Kuh-e-Sekonj range-front to the southeast of Mahan (e.g. ^{Berberian, 2005}). Mahan and Jorjafk faults are likely to be related due to the geometrical similarities and aligning along the same strike (^{Walker, 2006} and ^{Walker et al., 2010}).

Bafgh-Baghin fault system

The 170 km Bafgh-Baghin fault system is located along the Bafgh-Baghin mountain. This fault system with the NW-SE trend consists of eight segments. This is a right-lateral strike-slip fault with a small reverse component. There is no historical record of an earthquake in this area but the instrumental earthquakes with small to medium magnitude have been recorded at the depth of about 10 to 15 km (^{Sistanipour et al., 2012}).

Rafsanjan-Zarand fault system

It is a left lateral NE-SW trending fault that is extended from NW of Zarand city to SW of Rafsanjan city. This fault is created as a result of right-lateral strike-slip motion on the Bafgh-Baghin fault system in the SW and on the Kuhbanan fault system in the NE. ^{Sistanipour et al., (2012)} believe it is likely to have more intense earthquakes in this area especially at the intersection of Rafsanjan-Zarand and Bafgh-Baghin fault systems, in the future.

Zangi Abad-Kerman fault

It is a 30 km long fault that is situated close to the town of Zangi Abad and extends southwards to the city of Kerman which is considered as a probable fault on the geological map of Rafsanjan (^{GSI, 1992b}). ^{Hasanzadeh et al., (2011)}, who have identified this fault, have not found any evidence of Quaternary activity of this fault, i.e. deformation and incision of Quaternary sediments, the geometry of alluvial fans and rivers. Since this fault extends to the city of Kerman and is located along the active fault of Kuhbanan, it should be considered in assessing the seismic hazard of Kerman city (^{Hasanzadeh et al., 2011}).

Materials and Methods

The hydrological modeling process was done in the study area by Arc GIS 10.3 software. 30 m resolution Digital Elevation Model (DEM) and geological maps, four geological maps of Iran at 1: 250,000 Scale (i.e. Kerman, Bam, Rafsanjan, and Sirjan) and at least twenty four geological maps of Iran at 1: 100,000 Scale (e.g. Kerman, Gowk, Rayen, etc.) were used as basic data, to obtain the drainage network layer, basins area, and water dividing layer. The study area was divided into 51 basins and seven geomorphic indices were assessed including stream length-gradient index (SL), drainage basin asymmetry (Af), hypsometric integral (Hi), the ratio of valley-floor width to valley height (Vf), index of drainage basin shape (Bs), index of mountain front sinuosity (Smf), and Transverse Topographic Symmetry Factor (T). Then, the results were classified into three classes based on the ^{El Hamdouni et al., (2008)} classification (except for SL index and T factor (Mosavi and ^{Arian, 2015})) (Table 1). After that, for each basin, two kinds of relative tectonic activity index (Iat) was obtained by the average of the different classes of geomorphic indices (S/n) (El Hamdouni et al., 2008); Results from the analysis of all the 7 geomorphic indices (Af, Hi, SL, Smf, Vf, Bs, and T) are accumulated and expressed as Iat_T index, but the T factor was not considered in the Iat index calculation (Table 1). The study area was classified into three classes (including high (class 2), moderate (class 3), and low (class 4) tectonic activity) according to the Iat and Iat_T values. All the calculated results have been shown in table 2. Finally, to identify faults with seismic hazard potential to the city of Kerman, which is the aim of this paper, the plotted earthquake epicenters of the study area were being compared with the morphometric analysis derived from GIS.

Table 1 The geomorphic indices classifications were used in this study. * ^{El Hamdouni et al., (2008)}. ** Mosavi and ^{Arian (2015)}. 

Table 2 Values and classes of seven morphometric indices. 

Geomorphic indices of active tectonics

Asymmetric factor (Af)

The drainage basin asymmetry factor is a way to measure the tectonic tilting of catchments and is defined as follows (e.g. ^{Keller and Pinter, 2002}):

where Ar displays the basin area on the right side of the stream and At indicates the total area of the drainage basin. When a basin is developed under stable conditions with little or no tilting, Af approaches 50, whereas values of more or less than 50 would be considered as the asymmetry of a basin caused by active tectonic (^{Keller and Pinter, 2002}).

Hypsometric integral index (Hi)

This index reveals complex interactions between erosion and tectonics by indicating the volume of a basin which has not been eroded (^{El Hamdouni et al., 2008}). The distribution of elevation is defined by the hypsometric curve and the area below the curve is defined as the hypsometric integral (^{Strahler, 1952}). It is calculated as follows (e.g. ^{Keller and Pinter, 2002}):

Hi = (average elevation - min. elevation) / (max. elevation - min. elevation)) (2)

The maximum and minimum altitude is directly extracted from the topographic map and the average altitude is calculated from the average value of at least 50 scattered heights in the basin (^{Pike and Wilson, 1971}) or is obtained through the analysis of the Digital Elevation Model (DEM). The analysis of the Hi index includes two sections: the hypsometric curve and the integral. The shape of the hypsometric curves is related to the erosional stage. young landforms are characterized by high Hi values and convex curves. old landforms have low Hi values and are associated with S-shaped curves. Intermediate Hi values demonstrate mature landscapes (^{Keller and Pinter, 2002}; Pérez-Peña et al., 2009).

Valley floor width-valley height ratio (Vf)

The Vf index is defined as the ratio of the width of the valley floor to its average height (e.g. ^{Bull and McFadden, 1977}) and it is calculated as follows:

where Vfw represents the width of the valley, Eld and Erd demonstrate the height of the right and left sides, respectively and Esc shows the height of the valley floor. This index indicates the difference between V-shape valleys, that with the low Vf values represent tectonically uplifted areas, and U-shape valleys which have high Vf values represent the dominance of the erosional processes (e.g. ^{Bull and McFadden, 1977}). The Vf and Smf indices are so helpful to demonstrate the relationship between faulting and morphology in an area (^{Silva et al., 2003}).

Basin shape index (Bs)

The relatively young drainage basins in the active tectonic areas tend to be elongated perpendicular to the topographic slope of a mountain. Because of tectonic activity reduction, the elongated shape of the basin gradually transforms into a circular shape (^{Bull and McFadden, 1977}). The basin shape is converted into a quantitative index by using the elongation ratio which is computed by (e.g. ^{Ramirez-Herrera, 1998}):

Bs = Bl / Bw (4)

where Bl indicates the length of the basin from the headwaters to the mouth and Bw represents the basin width. basins with high tectonic activity have high values of the Bs indices and indicate elongated drainage, whereas low tectonic activity basins are characterized by low Bs values and a circular shape (^{El Hamdouni et al., 2008}).

Mountain front sinuosity index (Smf)

This index represents the balance between erosion forces, which tend to create irregular mountain front, and tectonic activities, which tend to create straight mountain fronts. The Mountain front sinuosity index is defined by equation 5 (e.g. ^{Bull and Fadden, 1977}):

Smf = Lmf / Ls (5)

where Lmf denotes the mountain front length along the connection of mountain and plain, and Ls is the straight-line length of the mountain front (^{Bull and Fadden, 1977}). Active mountain fronts have a straight front with low Smf values (Smf < 1.1); whereas inactive mountain fronts, due to adapting by erosive processes have irregular or sinuous fronts with high Smf values (Smf < 1.5) (^{Keller and Pinter, 2002}; ^{Silva et al., 2003}; Bull, 2007; ^{El Hamdouni et al., 2008}; Pérez-Peña et al., 2010).

The Stream length-gradient index (SL)

This index was first used to identify the recent tectonic activity by ^{Hack (1973)}. It evaluates the relative tectonic activity by concerning the relationships between rock resistance and topography along the drainage basin (^{Hack, 1973}). This index is calculated as follows:

SL = (ΔH / ΔL) L (6)

where ΔH/ΔL represents the gradient of the stream, and L displays the total lengthofthechannelfromupstreamtothemidpointofthechannelreachesforwhich the index is calculated (^{Hack, 1973}). The SL index is sensitive to sudden change along the stream profile, this anomalous change may indicate a tectonic uplift or sudden lithological change. The anomalously high Sl values for low resistance rocks show tectonic activity (^{Keller and Pinter, 2002}).

Transverse Topographic Symmetry Factor (T)

This index is appropriate for assessing drainage basin asymmetry. The changes in the tilting rate in different zones of an area demonstrate the difference in the uplifting rate. This index is obtained by the following equation:

T = Da / Dd (7)

where Da denotes the distance from the midline of the drainage basin to the midline of the active meander belt, and Dd indicates the distance from the midline to the basin divide (^{Cox, 1994}). The T factor is a vector that has direction and magnitude ranging from zero to one. T value equals zero for a complete symmetric basin and as the asymmetry of the basin increases, the T value rises and approaches 1. (^{Cox, 1994}; ^{Keller and Pinter, 2002}).

Figure 2 Distribution of the Af index in the Study area, |Af-50| > 15 (class 1), 7 < |Af-50| < 15 (class 2), |Af-50| < 7 (class 3) (^{El Hamdouni et al., 2008}) 

Results of Geomorphic indices analysis

Asymmetric factor (Af)

The asymmetry factor (Af) values were computed using DEM in a GIS environment. This index has been calculated as |Af-50|, which is the amount of difference between the observed value and the neutral value of 50 (^{El Hamdouni et al., 2008}). The calculated | Af-50| values in the study area vary between the range from 0.3 to 40.4. Figure 2 illustrates the distribution map of the Af index. According to ^{El Hamdouni et al., (2008)} classification, the calculated |Af-50| values along the study area have been categorized into 3 classes: Class 1 (|Af-50| >15); Class 2 (7 < |Af-50| < 15) and Class 3 (|Af-50|< 7) (Table 1). Class 1 which represents high tilting includes more basins than the other two classes, i.e., 24 out of 51 basins (Table 2). Af values of 11 basins belong to class 2 and 16 basins to class 3. Thus, most of the drainage basins are asymmetric in this study area. Basin No.33 has the lowest Af value (|Af-50|=0.3) and the highest Af value (|Af-50|=40.4) belongs to basin No.34 (Table 2). The three classes of Af index have scattered across the study area and most of the tilting is related to the active tectonic. For instance, evaluating the Af values indicates a dominant NNW tilting along drainage basins in the west of the study area which is the location of the major faults of eastern Iran, such as Gowk and Nayband faults, and the earthquake epicenters concentration (mostly destructive earthquakes with MW>5).

Hypsometric integral index (Hi)

By using DEM in a GIS environment, the Hi values were calculated. The Hi values which range from 0.06 to 0.48 are classified into two classes: Class 2 (0.4<Hi<0.5) and class3 (Hi<0.4) (Table 1). Figure 3b displays the distribution of the Hi index in the study area. Basins No. 2, 9, 16, 18, 20, 48, and 50 fall into the intermediate range of Hi value (Table 2) and display an S-shaped hypsometric curve (Fig.3a). So, these basins are considered as the mature stage of the Davis erosion cycle. The 44 other basins with a low range of Hi values (Table 2) and concave hypsometric curves (Fig.3a) demonstrate the old stage of the Davis erosion cycle (^{Davis, 1899}; ^{Keller and Pinter, 2002}; Pérez-Peña et al., 2009).

Figure 3a Hypsometry curves of 4 basins with their calculated Hi. Basins No. 16 and 20 with the moderate value of Hi (class 2) show an S-shaped hypsometric curve. Basins No. 3 and 5 with the low value of Hi (class 3) show concave hypsometric curves. (A) is the total surface of the basin. (a) is the surface area within the basin above a given elevation h. (H) is the highest elevation of the basin. 

Valley floor width-valley height ratio (Vf)

The parameters of Vf (Vfw, Eld, Erd, Esc) were measured for main valleys that cross mountain fronts using cross-sections drawn from a DEM-derived contour map in the GIS environment. The Vf values were calculated for as many as valleys that were possible, 300 m upstream from the mountain front. The mean of the Vf values is determined as the Vf index which ranges between 0.11 and 0.82 in the study area. According to this index calculation (Table 2), 22 basins are classified into class 1 (Vf < 0.5) which represent the deep V-shape valleys and are usually accompanied by tectonic uplift. 14 basins are categorized in class 2 (0.5 < Vf < 1). This index was not calculated for 15 basins due to the lack of appropriate conditions. These results in addition to the absence of class 3, imply that the values of Vf are relatively high for most of the study area. This activity is mainly caused by faults such as Gowk (in basin No. 16, 17, 22, 31, and 42), Nayband, Ravar, and Lakarkuh (in basin No. 2, 7, 16, and 17), Rafsanjan (in basin No. 10, 15, 18, 19, and 24), Mahan-Jupar (in basin No. 19, 21, and 22) the northern tip of Sabzevaran (in basin No.32, 33, 37, 41, and 49), and the southern tip of Kuhbanan (in basin No. 16, 17, and 22). The distribution of the Vf index and locations where calculations of the index are made are shown in figure 4.

Figure 3b Distribution of the Hi index in the study area, Hi > 0.5 (class 1), 0.4 < Hi < 0.5 (class 2), Hi < 0.4 (class 3) (^{El Hamdouni et al., 2008}) 

Basin shape index (Bs)

Bs index was computed using DEM in the GIS environment. According to the Bs index results, which vary between 0.33 and 5.25 (Table 2), most of the basins (46 out of 51) have low Bs values (class 3 Bs < 3) which is the characteristic of circular basins with low tectonic activity (^{El Hamdouni et al., 2008}). Basins No. 12, 15, and 28 with a medium value of Bs belong to class 2 (3 < Bs < 4). Basins No. 1 and 47 with a high value of Bs are placed in class 1 (Bs > 4), in other words, these are basins with more elongation and relatively higher tectonic activity (^{El Hamdouni et al., 2008}). The distribution of the Bs index for the study area is plotted in figure 5.

Mountain front sinuosity index (Smf)

In each basin, Lmf and Ls values are measured from a DEM-derived shaded-relief topographic map in a GIS environment. The Smf values were measured for proper prominent mountain fronts in each basin. The average of obtained Smf values is determined as the Smf index. Figure 6 displays the distribution of the calculated Smf values for the 51basins in the study area (Table 2). These values vary from 1.06 to 3.12. The basins No.7 and 17, which are in the NE of the study area, are categorized in class 1 (Smf <1.1). These basins indicate the straight state of mountain fronts which specify high active tectonic and uplifting (^{Bull, 2007}). Such high Smf values were expected for these two basins since important faults such as Nayband, Gowk, and Shahdad pass across them (Fig. 1b and Fig. 10a and b). Most of the basins (28 basins) are classified into class 2 (1.1<Smf<1.5). All ten basins of class 3 (Smf >1.5) are mainly located in the SW of the study area. Due to the lack of appropriate conditions, this index was not calculated for 11 basins.

Figure 4 Distribution of the Vf index in the study area, Vf < 0.5 (class 1), 0.5 < Vf < 1 (class 2), Vf >1 (class 3) (^{El Hamdouni et al., 2008}) 

Figure 5 Distribution of the Bs index in the study area, Bs > 4 (class 1), 3 < Bs < 4 (class 2), Bs < 3 (class 3) (^{El Hamdouni et al., 2008}) 

The Stream length-gradient index (SL)

The SL values were computed using DEM (30 m resolution) and a DEM-derived contour map, for every 100 m in the GIS environment. Elevation difference, distance along the main drainage of the study area, and total distance from the source were determined in the GIS environment. Using a methodology defined by ^{Hack (1973)}, the SL index was calculated with the obtained inputs. The SL values have been categorized according to Mosavi and ^{Arian's (2015)} classification: class 1 (SL>1000), class 2 (500<SL<1000), and class 3 (SL<500). The SL values range from 67.34 to 1419.42. The longitudinal river profile and SL index were analyzed to reveal that basins No. 7, 16, and 31 are placed in class 1, 7 basins in class 2, and 32 basins in class 3 (Table 2). This index was not calculated for 9 basins due to the lack of appropriate conditions.

SL values along the channel were plotted on the map of rock strength levels and shown in figure 7a. The Contour lines and distribution of the SL index are displayed in figures 7b and 7c respectively. Also, the longitudinal river profiles were drawn for each basin, some of which are shown in Figure 8. Most of the recognizable SL anomalies are related to active tectonic.

Figure 6 Distribution of the Smf index in the study area in addition to the mountain front and geological structures of the study area. Smf < 1.1 (class 1), 1.1 < Smf < 1.5 (class 2), Smf > 1.5 (class 3) (^{El Hamdouni et al., 2008}). 

Figure 7a SL index along the channels and rock strength levels of the study area. 

Figure 7b Contour lines were produced out of the SL index. 

Figure 7c Distribution of the stream-length gradient (SL) index in the study area, SL > 500 (class 1), 300 < SL< 500 (class 2), SL < 300 (class 3) (^{El Hamdouni et al., 2008}). 

Figure 8 Some longitudinal river profiles in the study area and the measured SL index. 

Transverse Topographic Symmetry Factor (T)

Using DEM (30 m resolution) in GIS environment values for Da and Dd were obtained. Then according to the ^{Cox (1994)} equation, the initial values of the T factor were calculated. After that, Tectonics FP software was used for statistical analysis of the T factor. Finally, the mean T values were computed for each basin and were shown as a vector in which its direction displays the tilting direction and its magnitude is T factor (Fig. 9a).

Figure 9 a) Steps to calculate the tilting index of basin No. 5 in the study area. Asymmetric vectors and polar representation of asymmetric vectors. the center of the net shows zero tiltings and its perimeter 1. b) Distribution of the tilting index (T) in the study area with the geological structures of the study area. 1 > T > 0.4 (class 1), 0.2 < T < 0.4 (class 2), 0 < T < 0.2 (class 3) (Mosavi and ^{Arian, 2015}). 

Based on the results (Table 2), the lowest value of T (0.06) belongs to basins No.29 and 32. Basins No.2, 26, 27, 31, 38, and 42 have the highest T value (0.95). The results have been grouped into 3 categories: class 1 (T>0.4), class 2 (0.2<T<0.4), and class 3 (T<0.2) (Mosavi and ^{Arian, 2015}). The T distribution map, shown in figure 9b, indicates that most basins (27 out of 51) have high tilting and are categorized in class 1 (Table 2). Class 2 and class 3 include 14 and 10 basins respectively (Table 2).