Morphotectonic analysis of the Upper Guajira, Colombia. A GIS and Remote Sensing approach

This study uses Landsat 8 and Digital Elevation Models (DEM) to show the dominant orientations of digital lineaments and morphotectonic features between measured trends and the tectonic evolution of the Upper Guajira, Colombia, in the northernmost region of the South American plate. Data from Landsat-8 and hillshaded images of three Digital Elevation Model (DEM) images with different resolutions (SRTM: 90m, ASTER-GDEM: 30m and Alos-Palsar: 12.5m), were used for the extraction and mapping of morpholineaments, drainage network and morphological features. Lineaments were analyzed by means of north azimuth frequency, length, density distributions, lithological distributions and geochronological periods. Tectonic control was supported by using the digitized geological map created by the Colombian Geological Service (SGC). Lineaments and faults were analyzed through the interpretation of a Riedel shear model as a result of the transtensional/transpressional tectonic arrangement of the Caribbean and South American plates. The directional trends of lineaments and faults indicate two dominant directions: NE-SW and E-W. The azimuth distribution analysis of measured structures and drainage channels show similar trends, except for some differences in the predominant directions of some drainage channels. The similarity in the orientation of lineaments, faults and drainage system highlights the degree of control exerted by underlying structures on the surface geomorphological features. Drainage channel classification illustrates the morphological and neo-tectonic complexity of the region. The extracted lineaments were divided into five geochronological groups based on the main ages of host rock formations according to the Colombian Geological Service (SGC) geological map. From the Cretaceous onward, the lineament azimuth frequency rotates from a NE-SW trend to a prominent E-W direction, which resembles the translation that Caribbean plate has been undergoing since the Cretaceous. Our results confirm that Remote Sensing techniques are reliable and useful to study the morphotectonic of an area and can be applied to zones of difficult access.


Introduction
Northern South America has a long history of orogeny and sedimentary basin formation driven by the influence of the Caribbean, Nazca and South American plates [1][2][3][4]. The southwestern margin of the Caribbean plate is a wide zone of approximately 600 km at the northernmost limit of the Andes cordillera [5][6][7][8][9][10][11][12], and its evolution has been described as a complex interaction between different tectonic terranes bounded by major strike-slip, and thrust systems [13][14][15]. The Upper Guajira (Fig.1) is a geologically complex area in the northernmost tip of South America, which developed as a result of several deformation phases since the Proterozoic, presenting several distinct geologic, tectonic and morphologic features [9,[16][17][18][19][20]. The most prominent feature of the study area is the Cuisa fault (CF), which is a high angle strike-slip fault with a right lateral displacement of approximately 15 km, extending from west to east [21,22], and that has been related to a Late Cretaceous-Paleogene convergence event [20]. We identify previously unknown morphotectonic features of  [70]. Units were divided by its chronological period.
Three main collisional orogenesis and a period of regional taphrogensis affected the NAB before the initiation of the Northern Andean Orogeny (e.g. The mid-Proterozoic Grenville Orogeny, the mid-Ordovician-Silurian Quetame Orogeny, a Permo-Triassic Orogeny and a taphrogensis event referred as "Bolivar Aulacogen" after Triassic Orogeny [13,[58][59][60]). The Andean orogeny started in the Cretaceous as a result of the subduction and convergence between the Caribbean and South American plates. The more recent theory asserts that the Caribbean Plate emerged from the Galapagos Hotspot, which sourced an oceanic plateau giving raised to a large coordinates, with a resolution of 1 arcsecond (30m). Finally, a DEM with a 12.5m resolution was used, which is one of the products of the ALOS mission carried out by the Japan Aerospace Exploration Agency (JAXA) from 2006 to 2011. PALSAR was one of the three instruments on the satellite ALOS which provided L-band synthetic aperture radar DEM files from ALOS mission were downloaded from the NASA website. Table 1 the specification of the downloaded satellite images. Preprocessing of DEM images is a required step to avoid artificial spikes and isolated pixels that will affect the generated surface [84,85]. These anomalies represent local artifacts, arising from different factors such as feature matching techniques, coarse spatial resolution or human activity. For this reason, it is important to consider how much the selected zone is affected by those conditions, in order to improve the image processing and interpretation. The study area presents few recent anthropic structures; in fact, most of the area is Uninhabited because of the arid conditions. An algorithm called Fill was applied to each DEM source using ArcGis10.5, removing the imperfections. The Colombian geological service (SGC) through its virtual portal (https://srvags.sgc.gov.co/arcgis/rest/services) provides freely available geological maps of the whole country, digitized in mxd format readable in ArcMap. The map of the Upper Guajira used in this study was created by Rodríguez and Londoño (2002) [70] at a 1:250000 scale. From this map, the lithological, structural and hydrological information was extracted in order to be compared with the results obtained from the automatic delineation of lineaments and drainage network.

Digital processing for delineation and identification of linear features
In order to identify geological lineaments, linear features from the DEM and Landsat-8 images were extracted. Hillshade techniques are widely used in tectonically active regions in order to delineate linear features with preferential trends [24,29,31,45,47,50,[86][87][88][89]. Hillshade images simulate a lighting effect from a base point with specific solar altitude and azimuth values, and the linear features showing the maximum enhancement and contrast are those that are perpendicular to the direction of the light source [50]. To extract lineaments with different orientations, eight hillshaded raster images were created for each DEM using the hillshade tool in ArcGIS 10.5 (Fig.2). The first shaded relief image was illuminated from the north (i.e. the solar azimuth was 0 • ) , with a solar elevation of 45 • (ArcGIS default); for the other seven relief images, an identical solar elevation was used along with seven other illumination directions, including: 45 • , 90 • , 135 • , 180 • , 225 • , 270 • and 315 • . Azimuth illumination and tilt values tend to hide or highlight morphological features in a hillshaded image, depending on their orientation with respect to the simulated light source. The next step consisted of applying the MDOW (multi-directional oblique weighted) shaded relief method [29,91], in which two multi-illumination hillshaded images were created from combinations of the individual shaded relief images. These combination images highlight lineaments and topographic features with a variety of directions, which would not be visible if only a single shaded image was used. The first four shaded relief images were superimposed and combined to produce one image with multi-illumination directions (0 • , 45 • , 90 • , 135 • -71°20'0"W 71°30'0"W 71°40'0"W 71°50'0"W 72°0'0"W 72°10'0"W 71°20'0"W 71°30'0"W 71°40'0"W 71°50'0"W 72°0'0"W 72°10'0"W 71°20'0"W 71°30'0"W 71°40'0"W 71°50'0"W 72°0'0"W 72°10'0"W Combination 1), and the four remaining images (180 • , 225 • , 270 • and 315 • -Combination 2) were processed in the same way (Fig.2). For a comparative and complementary analysis in the detection of geological lineaments, bands 7 and 8 of the Landsat-8 data were used. Band 7 has been useful for lithological discrimination and the extraction of morphotectonic features such as lineaments [29,47,48,93]. The panchromatic band (band 8) has a resolution of 15m allowing the extraction of a greater quantity of useful information, mapping of lineaments, and detection of physical characteristics of the terrain in greater detail.  [24,29,50] and the software manual. Each vector layer was upload in ArcGIS 10.5 in order to be processed and analyzed. Some lineaments were paired with anthropogenic structures (i.e. roads, buildings, etc.), coastal limits and/or possible errors related to the cloudiness in the images (i.e. ASTER gDEM and Landsat 8), so after a visual review, the lineaments corresponding to the mentioned pairings were deleted manually (Table.2).

Methodology used for the analysis of digital lineaments
The numerical values of lengths and directions were obtained for each set of lineaments using the ArcMap 10.5 COGO tool. For each group of lineaments rose diagrams were created to determine the dominant azimuth direction and its physiographic subdivisions. Determining the dominant direction based on the lineament frequency can include some errors given that the lineaments drawn automatically are usually short in length, and many individual samples could represent a long and unique line [29], so the rose diagrams allow the identification of the lineament lengths. Rose diagrams were done using dynamic tables tool in Microsoft Excel 2016 software. The parameters of the rose diagrams carried out in this research are [24,88,89]: A) The number of lineaments (N) determined the number of lines over a specific area. B) The total length of lineaments defined the sum of lengths of all lineaments in the study area. C) The morphotectonic network density D is defined by the following equation: [Density = N/A]], where A is the area of a particular region. D) The azimuth trend of lineaments determined the lineament orientation respect to a specific geographic coordinate system. The azimuth of lineaments is illustrated by rose diagrams, which are divided into 72 intervals of 5• (i.e. In total 360•) moving clockwise.
Different types of density maps of each lineament sample and their merging derivatives were created in order to understand their spatial distribution, the highlight faults, weakness zones and high topographic changes [30,44,51,94]. Furthermore, lineament samples were ordered and distributed with respect to the age of the geological units according to the SGC geological map (1C), thus evaluating the azimuth trend along the geological timeline [29]. Using the fault polylines of the SGC map, we created a Riedel Shear model which allowed us to understand the historical activity of faults due to strike-slip motion, such as the CF. This kind of model can only be applied in zones where strike-slip faults occur, so transtensional and transpressional regimes, such as the ones that have occurred along the boundary between the Caribbean and South American plates, are good areas to use this kind of model. Thus, some lineaments can be related to deformational elements (i.e. R and R' fractures) formed by the strike-slip motion due to the regional kinematics based in in the Riedel model classification, by assuming that regional linear features can be classified with this method within a zone that presents a transpressional tectonic arrangement.

Lithological discrimination from Landsat-8 multispectral data
In addition to the linear analysis of faults and lineaments, Landsat 8 images were used to show typical morphological features of the tectonic regime that have evolved over the area. For this research Landsat-8 multispectral images integrated with GIS techniques applied using ArcGIS 10.5, PCI Geomatica 2012 and ENVI 5.3 to map and discriminate the lithology and structural features of selected zones of the Upper Guajira. Multispectral images consists of the three visual primary color bands (red, green, blue) in addition to the longer-wavelength bands. These three bands (432) may be combined to produce a TCC (True Color image), which is one of the techniques used to select bands for a visual interpretation. The Principal Component Analysis (PCA) transformation is a multivariate statistical technique that selects uncorrelated linear combinations (eigenvector loadings) of variables, in such a way that each component successively extract linear combinations of all the digital numbers of the other band images and has a smaller variance. PCA were computed for the Landsat-8 data to discriminate physical features and lithological units in the study area. These PCA images were assigned to the display colors (Red, Green and Blue) respectively to form PC-color composite image in order to differentiate between the different rocks assemblage of the area as these are the three components that best reflect lithological variations [48]. The PC-5, PC-3 and PC-2 composite proved to be significant to discriminate between geological features.

Delineation of drainage networks
In tectonically active regions, drainage networks are often strongly influenced by the orientation of geological structures [29,32,88], reflecting the interaction between surface processes with the faults and folds that have been formed over a tectonically active area [96][97][98][99]. Drainage networks may indicate local deviations from the predominant patterns, which help identify neotectonic activity in different tectonic settings [29]; thus, the automatic identification of drainage networks and drainage basins from DEMs has become an important tool to include in a morphotectonic analysis [100]. Modern GIS software packages make it possible to extract an automatic drainage network and drainage basins from a DEM [26,[28][29][30]101]. In this work, the SRTM DEM was selected for the extraction of drainage network and basins, because its relatively low spatial resolution (90 m) allows for rapid processing [30]. This automated delineation of the drainage network was done by using the GIS Flow accumulation method [42] based on the hydrological manual made by ESRI (2011) [102]. The tool Hydrology of ArcMap 10.5 was used for the stream network delineation. The extraction of drainage patterns is a multi-step process and includes: (1) Sink filling (Pre-processing), (2) flow direction identification, (3) flow accumulation calculation, (4) stream definition and segmentation, and (5) Strahler Stream order classification. The threshold raster defines the cell values that are determined to have enough accumulation to be classified as a stream, and the threshold raster is generated from a raster created by the Flow Accumulation tool. Flow accumulation raster is in simplest form the number of upslope cells that flow into each cell. By applying a threshold value to the results of the Flow Accumulation, a stream network can be delineated. Low value thresholds lead to an increase in the level of detail of the stream networks [29], so for this research the value was 500, value that provides to be workable for drainage network deliniation. The network streams were split into several segments based on their orientations using the "Split Line at Vertices" tool in order to segment the vector stream polylines, thus avoiding errors in the analysis due to the sinuosity and topography of the area. The COGO toolbar of ArcMap was also used to extract the orientation and length for each drainage channel segment. Knickpoints (sometimes referred to as knick zones), are discontinuities or steep segments in a river profile which are usually a response to stream power, climate changes, lithological instability, or tectonic deformities [100,103,104]. In this study we also attempt to identify knickpoints on longitudinal stream profiles in order to assess whether they have lithological and/or tectonical disturbances [107,108,140].

Morphometric analysis of lineaments
Here we evaluate the automatically extracted lineaments in order to analyze their tectonic implications. Lineament density maps of frequency, intersection and length are used in order to visualize and characterize the spatial patterns of the lineaments [44,109].

Length Analysis
Morphometric results from the three DEM images are shown in Figure 3, in which it can be seen that the Alos Palsar scenes detect more than twice the number of lineaments detected by the other two DEM images, with lengths mainly between 0-1 km. The length of the polylines from ASTER and SRTM images are similar, but SRTM lineaments show a greater number of samples due the cloud errors over the Macuira zone in the ASTER image. Lineaments from Landsat 8 have the shortest lengths, but band 8 samples are greater in number due to the higher resolution of band 8 compared to the other bands. Table 2 summarizes the results of the lineament samples detected from the satellite images and the lenght of the lines that represent the faults of the Upper Guajira, according to the SGC digital geological map.   All density maps were created using a radius of 1 km. The zones of maximum density anomalies probably indicate a higher intensity of deformation. The legend below the maps indicates the lengths of the digital lineaments, and that the percentage in the pie diagrams indicate the percentage the lineaments for each combination)

Lineament Frequency Density Analysis
The purpose of the lineament frequency density analysis is to calculate the distribution of number of lineaments per unit area [43]. This analysis is also known as lineament frequency analysis [94]. We produced No. Lin/km 2 ) and maximum (> 8. 9 No. Lin/km 2 ). In in a large part of the study area of the study area (37%), the density of lineaments at a maximum, corresponding to the zones with notable topographical changes (i.e. physiographic divisions of the area), while the areas with low (9%) and moderate (17%) densities correspond to flatter areas. Areas with high and very high densities comprise 19% and 18% of the region respectively. Studies have shown that areas with a high density of lineaments tend to have a high intensity of deformation [110], high degree of rock fracturing [111], shearing [112], permeability [113], higher ground water yield [114], mineral occurrences associated with hydrothermal alteration zones [115], higher soil erodibility [116], slope failures [115], steeper slopes [31] and/or regional faults [29].

Lineament Length Density Analysis
A lineament length density map is a useful guide to understanding the zonation of the lineament samples ( Fig.4B1) [44,117]. This analysis is also known as lineament-length density [94] and consists of calculating the lengths of the lines completely contained within a cell limit [31,43,118]. Thus, the total length per unit area is calculated using all the lengths of the lineaments that are completely contained within a specific cell, and is expressed in km/km 2 [44]. For each group of lineaments from each DEM image, a lineament length density map was produced, so the length density maps were extracted by merging the two lineament combinations of each DEM due to the similar lengths observed for both combinations in the lineaments extracted from each DEM (Tab. 2). A general length density map was created by merging all the lineament combinations from all DEM samples (Fig.4B1). The most distinguished feature observed in the four maps is the high concentration along the CF, which is geometrically enhanced on the maps, meaning that the longest fractures segments in the area are observed along the CF. Morphometrically, the lineament-length density maps show some variations depending on the type of image used: (1) The density measures were lowest for the Alos Palsar DEM, which is the one that presents the biggest population of lineaments (Fig.4B2). (2) In the ASTER gDEM, over the Macuira serrania there is a blank zone due to the clouds, which is reflected in the density measures ( lineaments extracted from each DEM, as well as two general intersection density maps were generated from the union of the all lineament lineament intersections of combinations 1 and 2 respectively. Based on a visual comparison between all density maps, a high correlation was found between areas with high density levels. The purpose of using intersection density maps is to identify areas of greatest fracture [44,119]. When the lineament segments do not intersect in the cell, the resultant map will be planar, with almost no density contours. The zones of high lineament intersection over the study area are feasible zones for the evaluation of fractures or the Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 23 October 2020 doi:10.20944/preprints202010.0476.v1   evaluation of potential water reservoirs because the areas with the greatest intersection lineaments tend to be areas in where water has a greater possibility of flowing [44]. The Alos Palsar DEM shows a greater quantity of intersections, enhancing more specifically the zones with highest densities compared to other DEMs intersections. ASTER gDEM lineament intersections show a lower quantity of intersections, but the density spatial distribution is more zonated than the SRTM patterns. The intersection density measures of the SRTM and ASTER individuals are very similar, unlike the Alos Palsar density measure which shows higher measure values due the quantity of the specimens. Combining the individuals of the combination 1 and 2 individually allowed us to refine more zones with high concentrations of intersections.  40  45  50  55  60  65  70  75  80  85  90  95  100  105  110  115  120  125  130  135  140  145  150  155  160  165  170  175  180  185  190  195  200  205  210  215  220  225  230 35  40  45  50  55  60  65  70  75  80  85  90  95  100  105  110  115  120  125  130  135  140  145  150  155  160  165  170  175  180  185  190  195  200  205  210  215  220  225  230 35  40  45  50  55  60  65  70  75  80  85  90  95  100  105  110  115  120  125  130  135  140  145  150  155  160  165  170  175  180  185  190  195  200  205  210  215  220 35  40  45  50  55  60  65  70  75  80  85  90  95  100  105  110  115  120  125  130  135  140  145  150  155  160  165  170  175  180  185  190  195  200  205  210  215 0 5 10 15 20 25 30 35   40  45  50  55  60  65  70  75  80  85  90  95  100  105  110  115  120  125  130  135  140  145  150  155  160  165  170  175  180  185  190  195  200  205  210  215 35  40  45  50  55  60  65  70  75  80  85  90  95  100  105  110  115  120  125  130  135  140  145  150  155  160  165  170  175  180  185  190  195  200  205  210  Orientation is one of the most revealing characteristics of linear features over the surface. Thus, linear features along the study area were interpreted using rose diagrams, whereby length and frequency of these linear features are also represented according to the azimuths tendencies. A detailed analysis of azimuth distributions of extracted lineaments and faults was performed for the whole area (Fig.6A) and individual subregions were delimited by physiographic regions. The Upper Guajira was divided into four physiographic units (Table 3), in order to compare them with the diagrams of the faults and lineaments. The units are as follows from west to east: a) Cabo de la Vela and Carpintero subregion (Fig.6B): Located at the northwestern zone of the Upper Guajira. This region presents the most unequal azimuth directions between the lineaments (NE-SW) and the reported faults (NW-SE). The low number of faults in the zone allows a visual descriptive analysis of this area, so our hypotheses for dissimilar azimuth directions are: (1) Faults are not very abundant and are only detected in metasediments and metaconglomerates rocks from late Cretaceous; with similar ages as the ones observed in Cabo de la Vela ultra-mafic complex (i.e. Ophiolite zone) from late Cretaceous (Campanian), in contact with poorly consolidated tertiary deposits of gravel and sand, so those faults could be controlled due the ophiolite accretion [57,120,121]. (2) The dextral movement of the Caribbean plate between Late Cretaceous-Paleogene may be related to the lineament azimuths in the zone, so the position of the faults could be associated to the rotation that the Guajira peninsula and SNSM underwent during these periods [3,76]. The faults over these zones could be geometrically positioned due to this rotation, so their original azimuth could be NE-SW but rotated (3) Lineaments in this subregion are not directly related to faults, being evident in the different linear orientations as they are probably related to the convergence of the Caribbean plate rather than by the motion of these faults, which do not seem to be identified through the analysis of lineaments. So lineaments do not highlight the faults reported by the SGC map. A visual interpretation of the lineament maps shows that the DEM with higher resolution better delineates the lineaments over this region. b) Cosinas and Simarua subregion (Fig.6C): Is the largest zone of all, merging the Cosinas and Simarua serranias where the CF deformation is more evident. Lineaments and faults show a strong relationship between their azimuth trendings. The biggest peak in the fault rose diagram of the faults corresponds to the CF system, the second one corresponds to the Cosinas and Simarua thrust faults system. c) Jarara subregion (Fig.6D): This zone also reflects a good relationship between the azimuth trending of faults and lineaments. The only difference between them lies in the particular distribution faults in the NW-SE direction, whose origin will be explained in next sections. d) Macuira subregion (Fig.6E): In this zone the NE-SW azimuth trends have predominance in faults and lineaments, but in the fault rose diagram a NW-SE trend appears. This differences in trends are a result of the normal faults family, which will be discussed later. Landsat 8 lineaments could not be identified due to the cloudiness over the area.
The preferential direction of all the extracted lineaments in all the physiographic regions and for each digital image is a NE-SW azimuth for the whole area (Fig.6A) Table 3. Statistics of the trend direction of lineaments in each studied area.

Chronological implications of lineaments
One of the most distinctive characteristics of the Upper Guajira peninsula is the series of isolated crystalline massifs between flat Cenozoic basins [57,122,123]. Geological units with ages ranging from Proterozoic to Neogene are exposed and can provide relative constraints on the lineaments that cross them. The method consists of detecting changes in lineament orientations through geologic periods and correlating these pattern changes with the geological evolution of the Caribbean plate. So the reliability of lineaments ages depends on the precision of the geological map [70]. All lineament samples of all DEM data were merged, and for each geological period units, lineaments were analyzed using the same morphometrical parameters that were used in the previous section. Only the lineaments that were totally within a polygon with an established period (according to the geological map) were used; the lineaments that were found cutting more than one geological unit were not taken into account for this analysis. The classification of all lineaments was performed and six populations were obtained: the ones observed in a Precambrian basement (Proterozoic eon), two Mesozoic families from Jurassic and Cretaceous and three from all Cenozoic periods (i.e. Paleogene, Neogene and Quaternary). Temporal evolution of the lineament trends in various geologic units is shown in Figure 7A. Its important to take into account that the age of the rock determines only a maximum age of the lineament formation. Real age can be much younger, but younger rocks should be high related to younger lineaments. One important consideration is that the lineaments of the Quaternary and Neogene periods were interpreted in the same population group. The lineament samples show a clockwise rotation of nearly 5 o per period since the Cretaceous (Figs.7A1 & A2). The evolution of this trend resembles the clockwise ENE rotation of the Caribbean plate [124,125] (Figs.7A3).The quantity of lineaments within each polygon that contain a group of lineament samples, allowing us to produce quantitative distribution per geological periods (Fig.7B). This classification consists of enhancing the rock units that show higher density and greatest number of lineaments [29]. The metamorphic units have the highest population densities of measured lineaments in lithological units and chronological units, however, as is shown in the maps, not all cases in which there is a greater quantity of lineaments will present greater density measures (Fig.7B). It could be observed that in spite of the fact that the sedimentary rocks have a large territorial area and a greater population of lineaments, the population density in this type of rocks is lowest, this could be related because the unconsolidated units present in the Upper Guajira ( [3,70]). Cretaceous units present higher densities and larger number of lineaments, suggesting that these units could present significant amount of deformation over the area, this could be also related to the collision of the Caribbean plate with the North and South American plate margins during the uppermost Cretaceous and Tertiary, process that may have been accompanied by a change in polarity of subduction, from eastward-dipping subduction of the Pacific below the Proto-Caribbean to westward-dipping subduction of the Proto-Caribbean below the Pacific (Caribbean) [3,4,61]. So, the highest densities (No. Lin/Km2) of lineaments (Fig.7B1) and the highest amount of lineaments (Fig.7B2) in the Upper Guajira are presented in cretaceous metamorphic units, such as phyllites, quartzites with serpentines and chert, and gabbros [70] (i.e. K2p, K2j, K2e, K2i and K2c geological units) .

Geometric analysis of the tectonic complex, a Riedel shear system approach
Faults are as a rule, indicated by lineaments [33]. Moreover, lineaments are generally associated with faults and linear zones of fracturing, indicating deformation patterns which may be expressions of a certain tectonic regime. Furthermore, geometric characteristics of faults and lineaments reflect the tectonic motions that have been occurring recently over an area. A semi-qualitative approach called Riedel Shear System, was applied in order to understand the geometric relations between faults and lineaments in the Upper Guajira [29,[126][127][128]. Riedel structures were obtained for each physiographic subregions and for the whole Upper Guajira region based on the geometry of the rose diagrams. Riedel shears are subsidiary shear fractures that propagate a relatively short distance out of the major fault, being coeval with it. The geometrical arrangement of Riedel shears indicates the sense of movement within the main wrench zone, widely used to interpret the kinematic evolution. The basic geometry of the Riedel structures consists of conjugate shear bands denoted by R and R'. R-bands are synthetic to the sense of the major fault and generally build the most prominent set. They develop at an acute angle (15 o ) clockwise to a sinistral main fault or anticlockwise to a dextral main fault. R'-bands are antithetic faults oriented at high angle (75 o ) conjugate with R-bands. The R-and R'-bands create an angle of about θ/2 and 90 oθ/2 to the general shear-zone direction, respectively, and intersect in an acute angle of β = 90 oθ , where θ is the angle of internal friction [127,129,130]. The Caribbean plate shows a dextral oblique movement with respect to the South American plate in northern Colombia, producing the strike-slip faults of OF and CF (Fig 1B). The deformation features are organized in a set of structures over or near the wrench zone that may form simultaneously, and which are geometrically oriented by the horizontal strain-ellipse, which are supposed to occur at 45 o in regard to the strike-slip major fault. Thus, the displacement of the CF produced a set of structures such as folds and thrust (in a direction perpendicular to the SHmax), normal faults and tension fractures and Riedel structures (in a direction parallel to the SHmax), all of them interacting in a geometrical arrangement. The present-day stress field shows how pre-existing faults can be reactivated, controlling the deformation of Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 23 October 2020 doi:10.20944/preprints202010.0476.v1 younger units. As explained above, the study area is characterized by a compressive maximum horizontal stress (SHmax) with a NW-SE direction, in which the master fault (CF) can therefore be interpreted as a dextral fault. The spatial pattern of lineaments can be interpreted as a dextral simple shear tectonic model. Thus, considering the SHmax with a NW-SE direction and considering the CF strike of N85 o E, the dominant system of NW-SE azimuth lineaments can be interpreted as synthetic dextral R-bands fractures (N85 o E). The NNW-SSE azimuth can be interpreted as antithetic sinistral R'-bands fractures (N20 o W). Based on the above facts, strain ellipsoids were defined for each area to understand the bulk geometry of deformation to create a model of the general kinematic system of the area (Fig.8). It is important to denote that thrust and normal faults can be also discriminated due to their geometries as well as fracture distributions by using this method, showing the reliability of the Riedel analysis for a big area, such as the Upper Guajira region.

Lithological and structural discrimination from multispectral data
Digital processing of Landsat-8 OLI/TIRS images using a variety of composites and enhancement techniques have been tested on selected structural or lithological elements in order to understand the landform and its structures due to the transpressive tectonic system. The Principal Component Analysis composite (PC5, PC3, PC2) has been shown to be accurate and useful for the identification and discrimination of lithology and structural features [48]. We used the geological map to better constrain our spectral result. Distinct structures such as faults, fold and particularly the Cuisa strike-slip master fault, in where a notorious displacement between geological units can be observed and inferred (Figs.9A, 9a1 & 9C). It was found that the CF is highly enhanced and shows a good lithological discrimination in the PCA composite image. Folds near the CF were also identified, showing at the north of the fault the axial surface in NW direction over Neogene limestones and clays (E3s), and at the south of the fault an axial surface in a NE trend over Jurassic arcillolites and shales units (J3sc and J1ch) (Fig.9a1). Some structural features such as faults and lithological contacts can be observed in the Jarara subregion, enhancing features related to deformation history of the study area. It was observed, by comparing the geological map of the Jarara area with the PCA combination of this zone, that faults are highly enhanced, especially over the Cretaceous unit (K2j), which is the one that contains this type of fault. Normal faults are also enhanced,  (Fig.9B &  9b1). A sedimentary basin in the southeast of the study area was also identified using the Landsat band composite (432) (Fig.9C). According to the geological map [70] (Fig.1C) the area exposed consist mainly by arcillolites and alluvial deposits with ages ranging from Neogene to Quaternary; located at the north of the western zone of the CF, being highly folded as it was showed over the Cosinas and Simarua zones (Figs.9a1), and denoting that the axial traces of the folds to the south of the CF (Figs.9a1) show similar NE orientation similar to the preferent lineament azimuthal trends, whereas the axial trace of folds to the north of the CF resembles a NNE trend more like R'-band Riedel shears (Fig.9C), suggesting that CF is still active by highlighting deformation structures over Jurassic to Neogene units from multispectral composite images. This pattern coincides with the dextral movement of the master CF. The drainage network consists of 44568 streams organized into 5 stream orders, according to the Strahler order [131,132]. Table 4 shows that the Stream order seems to increase at the same rate as the average length increase; except for the III order stream, which shows a slight deviation from what is expected. The general drainage frecuency density and length density over the Upper Guajira region is 5, 83No/km 2 and respectively 0.85Km/km 2 [133]. Drainage network results from variations in the porosity, topography, geology, structure and chemical composition of soils and rocks [134], which means that the drainage pattern records evidence of tectonic deformation [28,29,134,135]. In that way, drainage patterns are normally controlled by the tectonic activity, creating anomalies along the surface [136]; therefore, drainage anomalies analysis provide valuable information on the actual tectonic deformation. Drainage networks have been classified in a variety of studies, as it was done in this research (Fig.10) [28,29,134,136]. The dendritic pattern is the most common pattern along in a network, and is located in the flattest areas. Dentritic patterns have a treelike branching system, which consists of the junction of channels that flow in the same direction, but the branches are oriented in several directions [134]. Generally, this type of pattern occurs in flat-lying rocks, meaning that they do not have a significant structural or slope influence [29,137]. Parallel drainage patterns were also found, particularly in the SE hillside of the Macuira serrania and formed in areas where the water flow is highly controlled, indicating a zone with flat lowlands. Elongated forms are also found in areas of parallel landforms [136] such as the valley between the Macuira and Jarara serranias, where the Macuira fault separates the Jarara and Macuira serrania. Angular patterns are also found in the NE slope of the Macuira serrania, where the slopes are steep. This kind of pattern shows a transition from angular to dendritic patterns to the north. Annular patterns are less common than patterns mentioned above and reveal the influence of primary streams developed in a semi-concentric circular domelike hills on sedimentary rocks [134,136], such as in the Carpintero serrania. These patterns are probably associated with the accretion of the Ophiolites on the Late Cretaceous, as discussed above. To validate the drainage network extracted from the SRTM DEM, the rivers 3)

Val ue
High : 846m Low : ~1m SGC Rivers from the geological map were also used. The comparison shows that the morphology of the patterns between the extracted streams and the ones from the map are very similar, corroborating our results.

Drainage network density analysis
A density frequency map (No.streams/km 2 ) was done for the extracted drainage network (Fig.11A). A visual interpretation of the map shows that the zones with higher concentrations of drainages correspond to the plane zones that surround the nearby serranias. Moreover the combination of density and geological map that higher concentrations of drainage patterns are located over Paleogene(K), Quaternary(Q) and Neogene(N) lithotypes. With the use of the Buffer Tools on ArcGIS we classified the length of streams under the influence of lineaments in order to estimate the streams controlled by lineaments. Both lineament combinations from the SRTM data were buffered with the streams by using a 150m threshold. All the streams within the buffer zone were extracted and density and analysis of the segments was performed (Fig.11B). Results showed that about 1652.49m (25%) of the streams are influenced by lineaments. From the streams influenced by SRTM lineaments, 49% correspond to the first-order samples, and about 30% , 15%, 4.48% and 1.52% influence was found in second, third, fourth and fifth orders respectively (Fig.11B) .
These results indicate that most of the influence occurred within the first order streams, which suggests that the impact of lineament is strong in the areas that have steepest slopes and rough terrain (Fig.11B). The dominance of lineament over the streams is diminished over the sedimentary basins in where streams do not seem strongly dominated by tectonic activity.

Domain analysis of drainage network and aligned morphological features
A detailed analysis of azimuth stream channels was performed with GIS tools on individual domains. As for the lineaments and faults, every stream channel was segmented in straight lines (i.e. between two vertexes) in order to use rose diagrams to calculate its preferential azimuth direction. The line segments corresponding to streams were weighted by their azimuth frequency. 10 grid domains covering the area were analyzed (Figs.10A  & 12). The grids were chosen selecting the ones with a significant quantity of linear features within their area. Hence, Domains 8, 9, 10 are in the SE of the zone, the 14,15,16 domains are within the middle part and the remaining four (22,23, 24 and 29) are in the NE sector. For each picked grid, faults and SRTM lineaments frequency diagrams were constructed as well. A structural and tectonic control of the drainage network along the study area is observed.
The drainage networks tend to have azimuthal direction in all directions, especially in intervals of 45 o . Nevertheless, the interest of decipher the principal azimuthal trends of the streams is to observe if there is a relation with the morphological trends, such as faults and lineaments. Table 5 show the principal azimuth trending that linear samples have over each domain. The NE-SW and E-W trends are the most common directions of linear features. Domains 8 and 24 are the ones that present identical azimuth trends (Table 5) Examining only the higher order channels (i.e. 4th and 5th order), the zone is drained by five main streams: The Masayumachama in the NE sector, the Itapara draining in the north and the Silimahana draining at south; both in the middle sector of the study area, and the Urash and Manash in the SW sector (Fig.13). The Northern stream flows from SW-NE, starting its path in the Jarara serrania, bordering the Macuira serrania until the Caribbean Sea. The slope degree map (Fig.13) shows that the principal flow paths of the selected tributaries Streams influenced by SRTM lineaments. Violet colored segments of streams are controlled by available lineaments in the study area. A density map of the controlled streams showing control of lineaments over the streams is exposed here.
are above the areas with lowest slopes, and from this raster were extracted the different slope measures of the tributary paths. Table 6 show the morphometrical results of the selected tributaries. Gradient variations along the biggest streams can be identified by analyzing their path profiles, highlighting drainage anomalies [100,104,138,139]. These anomalies are commonly known as knickpoints, which generally correspond to the terrane response of climate changes, lithological variability, changes in the stream power and/or tectonic deformations [29,40,139,140]. Therefore, the study of longitudinal profiles gives information on the part over the streams that undergo higher influence due to the lithology and/or the tectonic environment. Is important to consider that these measures disregard factors such as stream power or climate changes. Another aspect to bear in mind is that the origin of knickpoints may not always be possible without a detailed investigation (Figs. 13 & 14). Figure 14 shows longitudinal profiles of the five picked tributaries, labeled form A to E as seen in table 6. As expected, these profiles have different and variable curves and gradients. The slope changes vary from 0 o to 20 o Stream C starts its pathway in the SW foothill of the Jarara serrania, and is about 4 km east of stream B in the middle of its path. Knickpoints are located in the first 5km where the topography starts to decrease. Other

Surface geomorphological features Vs. Subsurface structures, Tectonic Implications
Extracted lineaments from remote sensing datasets allowed us to identify several populations of lineaments with various lengths and orientations. Two main orientations were found: I) NE-SW and II) E-W. The lineament samples were mapped and analyzed in several ways such as density maps or rose diagrams, which enhanced the recognition of structural features such as the CF and structures related to a transpressive system in the zone. High lineament density zones over the density maps are generally related to zones of high topographic relief (i.e. physiographic subregions) and high fracture zones, moreover due to most of the high density areas are not covered by Quaternary and Neogene surficial deposits, the linear structural features of older geological units can be easily detected, allowing to observe older deformations over those units. The Upper Guajira exhibits different geological units and several prominent structures. Lineaments can provide relative constraints on the evolution of the zone by organizing them according to the geological units of different ages that are crossed by lineaments. All the extracted lineaments from remote images were grouped according to the age of the host rock. Our results show that the NE-SW trend continued from the Precambrian to Quaternary unit; in the Cretaceous age, the NE-SW trend becomes more prominent because the evolution converges from the late Cretaceous onward due to the drifting by the North America and South America plates since the late Jurassic [3,4,7]. The E-W trend is better documented from the Paleogene formations. It was during this period that the CF began to displace the Upper Guajira block from the SNSM-Lower Guajira block [74]; then the Guajira Peninsula started its northeastward translation, which may explain its current geographical position, despite having been formed tens of kilometers away [141]. The E-W azimuth trend also affected all exposed rocks from Precambrian to Quaternary, but it is very conspicuous in the Paleogene. Neogene to Quaternary units expose a similar azimuth trend, but a noteworthy point is that during these periods probably represents the lineaments associated with the continued oblique collision of the margin between the Caribbean and South American Plates. Occurrence and relative abundance of prominent trends (NE-SW and E-W) through geological time is proposed to be generated by repeated activation of pre-existing crustal structures during different tectonic periods as it was previously addressed. We also observe that in spite of the fact that the sedimentary rocks have a large territorial area and a greater population of lineaments, the population density in this type of rocks is the lowest. Cretaceous units are the ones that present higher densities and greater number of lineaments, suggesting that these units present important features of deformation over the area. Cenozoic tectonic activity in the study zone has been dominated by mainly compressional/strike-slip stress regimes, which produce the reactivation of crustal structures and recent deformation. In order to understand the present-day behavior of structures over the zone, the extracted faults from the geological map and lineaments were organized in a Riedel shear model, using the CF as the major fault. Results show that the current stress field of the area is characterized by a compressive maximum horizontal stress in NW-SE direction, so the spatial pattern distribution of lineaments can be interpreted as fractures related to the dextral movement of the CF. The E-W trending of lineament samples can be interpreted as the most recent ones, being parallel to the major strike fault; The NW-SE trending lineaments can be interpreted as synthetic dextral R-shear fractures; and the NNW-SSE trendings are interpreted as antithetic sinistral R'-band fractures, which are the less discernible family of fractures. Distinct geomorphological features typical of a transpressional system were portrayed by processing DEM and Landsat-8 images. Faults, folds and recent folded basins over the study area help to better understand the results generated by the lineament analysis. The current position of these tectonic blocks, together with their tectonic structures and geochronological associations of lithotypes, support a dextral transpressional tectonic regime over the study area that has been suggested for various authors [3,20,56]. It is important to denote that we can not assure all digital lineaments have a tectonic nature, however we suggest must of them have a tectonic significance because: 1) the azimuth of lineaments have a strong correlation with the azimuth direction of faults, so the major geological contacts have a NE trend similar to the major trends of lineaments. 2) Another fact is that the detected lineaments are in the areas in where the main faults occur, contrasting the lack of lineaments in the flat areas as it can be seen in the density maps (Figs. 3, 4 and 5).
3) The Riedel shear scheme shows that the principal azimuth of lineaments and faults are in the same orientation as the θ 3 tensor (Fig.8). We recommend more future work on the tectonic activity and physical conditions of the Upper Guajira region, by using geophysical data such as the EGM-2008 gravity data and SAR sensor data, with a seismic analysis in order to understand better the actual activity of the faults, specifically for the strike-slip CF.

Drainage network implications
Within the 10 selected domains where a clear drainage orientation is found, the azimuth distributions of streams are similar to those of the lineaments and faults, with few differences in the prevailing directions. In some cases, the lineaments and faults have a secondary NW-SE trend, which is a main trend in the stream distributions; nevertheless, the NE-SW is a conspicuous azimuth trend of various stream paths. This similarity in the azimuthal trends indicates that faults and lineaments are preferred zones for the water movement, because the water runs preferably over tectonic fractures. There are differences in the spatial distribution of lineament and stream densities over the area, showing that the higher density areas of streams are in the planar zones different than the lineament distribution. Streams were organized by a Strahler order showing that most of the drainage influence over the study area is controlled by first order streams, indicating that the tectonic impact over the streams occurred in the physiographic subregions. Distinct drainage patterns may be related to the movement or deformation produced by morphotectonic features and landforms. It was shown that the best correlation between lineament and drainage streams are presented in the zones with higher slopes, but its important to mention that there is necessary to improve this results within the specific basins (Fig.13) and a particular stream path, by using the normalized steepness index (ksn), SL Index and/or hypsometric profiles. Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 23 October 2020 doi:10.20944/preprints202010.0476.v1

Conclusions
This study shows the results of a combined use of remote sensing scenes and GIS techniques in order to interpret the complex structures of the Upper Guajira landscape, which has a geological history since the Precambrian age. Samples of lineaments were used to analyze the tectonic complexity of the study area. The analysis of quantity and length of extracted lineaments from eight different azimuth angles shows that results are more precise and accurate than just using one azimuth angle. The present paper reveals that physiographic zones over the study area, which contains the serrania ranges, are characterized by a high lineament density. Consequently, these areas have greater slope changes, groundwater infiltration potential and influence by recent tectonic activity, which implies that these zones should be avoided for engineering constructions, unless a detailed study of the terrane is performed. Major faults are located over the physiographic zones as well. The azimuth trends of all the groups of lineaments are predominantly in NE-SW and E-W directions, indicating an homogeneous structural deformation of geological units. The analysis of lineaments with regard to the geochronological units highlights of the movement of the Caribbean plate since the Cretaceous. The response of streams to the lineaments and active faults is observed in the form of clear changes in longitudinal profiles, modifying the gradient and stream patterns. Greater influence of lineaments on drainage networks is found on first order streams, characterized by a decrease of the dominance of lineaments as the stream order increases. A general similarity between lineament trends and in the direction patterns of stream segments and surface structural features suggest that the landform and the structural features of the study zone are highly controlled by a tectonic transpressive system. Evolution of lineaments is related to the reactivation of preexisting faults; probably due to the dextral movement of the CF, in response to the stress field. Therefore, NE-SW trends highlight oblique movement of the Caribbean plate under South American plate, whereas E-W trends highlights the rotation of the stress field that allows the CF to be notorious structural feature of the southern boundary of the Caribbean plate. Processing of Landsat-8 images allowed us to identify and enhance lithological units and landform structures which are related to the tectonic system that has been described. In addition this study provides valuable information that might be significant for the management and development of the region. Thus, our study represents a contribution to the fields of tectonics, structural geology, engineering geology and hydrogeology. Therefore, this study contribute in the knowledge of the morphological arrangement of the Upper Guajira, by exposing a practical guide that analyse different types of remote sensing data with GIS techniques.

Conflicts of Interest
The authors declare no conflict of interest.