[7378] | 1 | """ Created by: Jonathan Griffin |
---|
| 2 | Date: 29 July 2009 |
---|
| 3 | Purpose: This script searches output from Ross Wilson's event selector to find events that are within a geographic area. |
---|
| 4 | Maximum and Minimum longitude and latitude values are given, and a search distance within which the max/min lon/lat of |
---|
| 5 | the quake must fall within to be selected. |
---|
| 6 | """ |
---|
| 7 | |
---|
| 8 | |
---|
| 9 | search_dist = 0.75 # radius in degrees to search |
---|
| 10 | |
---|
| 11 | # List coordinates to search for |
---|
| 12 | MaxLon = 120.2605 |
---|
| 13 | MinLon = 113.7402 |
---|
| 14 | MaxLat = -9.4653 |
---|
| 15 | MinLat = -11.2732 |
---|
| 16 | |
---|
| 17 | # Specify input Folder |
---|
| 18 | inDir = r'/nas/gemd/georisk_models/inundation/data/western_australia/broome_tsunami_scenario_2009/anuga/boundaries/Results_Australia_1629_0.01_15.00/' |
---|
| 19 | inFile = 'fault.xy' |
---|
| 20 | #inFile = 'fault.txt' |
---|
| 21 | inFilePath = file(inDir + inFile) |
---|
| 22 | |
---|
| 23 | # Output file to list I.D. of possible quakes |
---|
| 24 | outFile = 'possible_quakes.txt' |
---|
| 25 | outFilePath = file(inDir + outFile, 'w') |
---|
| 26 | |
---|
| 27 | #Establish variable and lists |
---|
| 28 | quake = 0 |
---|
| 29 | possible_quakes = [] |
---|
| 30 | quake_list = [] |
---|
| 31 | count = 0 |
---|
| 32 | # Read data from fault.xy |
---|
| 33 | for line in inFilePath.readlines(): |
---|
| 34 | #Skip header row |
---|
| 35 | if "Lon" in line: |
---|
| 36 | continue |
---|
| 37 | #Split data and read into variables |
---|
| 38 | data = line.split(',') |
---|
| 39 | Lon = float(data[0]) |
---|
| 40 | Lat = float(data[1]) |
---|
| 41 | Quake_ID = int(data[2]) |
---|
| 42 | subfault_ID = data[3].strip('\r\n') |
---|
| 43 | subfault_ID = int(subfault_ID) |
---|
| 44 | |
---|
| 45 | # if new quake_ID, find latitude/longitude range and determine if |
---|
| 46 | # within search distance |
---|
| 47 | if Quake_ID != quake: |
---|
| 48 | #print 'quake list', quake_list |
---|
| 49 | if len(quake_list) != 0: |
---|
| 50 | lonlist = [] |
---|
| 51 | latlist = [] |
---|
| 52 | for i in range(len(quake_list)): |
---|
| 53 | lonlist.append(quake_list[i][0]) |
---|
| 54 | latlist.append(quake_list[i][1]) |
---|
| 55 | #print 'lonlist', lonlist |
---|
| 56 | maxlon = max(lonlist) # find max/min lat/lon here |
---|
| 57 | minlon = min(lonlist) |
---|
| 58 | maxlat = max(latlist) |
---|
| 59 | minlat = min(latlist) |
---|
| 60 | if ((abs(maxlon-MaxLon) < search_dist) & (abs(minlon-MinLon) < search_dist) & (abs(maxlat-MaxLat) < search_dist) & (abs(minlat-MinLat) < search_dist)): |
---|
| 61 | possible_quakes.append(Quake_ID) |
---|
| 62 | quake = Quake_ID |
---|
| 63 | quake_list = [] #reset quake_list for next earthquake |
---|
| 64 | quake_list.append([Lon,Lat,Quake_ID,subfault_ID]) |
---|
| 65 | count += 1 |
---|
| 66 | else: #append data to list of subfault coordinates for that quake |
---|
| 67 | quake_list.append([Lon,Lat,Quake_ID,subfault_ID]) |
---|
| 68 | |
---|
| 69 | # Print out all earthqaukes that fall within search area. |
---|
| 70 | print 'possible quakes', possible_quakes, 'out of ',count, ' earthquakes' |
---|
| 71 | |
---|
| 72 | #write out possible earthquakes to text file |
---|
| 73 | outFilePath.write('Possible earthquakes\n') |
---|
| 74 | for quakes in possible_quakes: |
---|
| 75 | print quakes |
---|
| 76 | outFilePath.write(str(quakes) + '\n') |
---|
| 77 | |
---|
| 78 | |
---|
| 79 | inFilePath.close() |
---|
| 80 | outFilePath.close() |
---|