Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Issues getting synthetics from FK database #277

Open
ammcpherson opened this issue Oct 23, 2024 · 5 comments
Open

Issues getting synthetics from FK database #277

ammcpherson opened this issue Oct 23, 2024 · 5 comments

Comments

@ammcpherson
Copy link
Member

Good afternoon,

I am trying to save synthetics from MTUQ to visualize against real data. I have done this with SPECFEM-generated Green's functions, but when I attempt to do the same thing using an FK database, no synthetics are generated. I cannot share the FK database of course, but these are the basic commands of what I am doing:

fk_path = '/store/wf/FK_synthetics/'+model
db = open_db(fk_path,format='FK')

origin = Origin({
        'time': '2022-09-24T15:18:54.694000Z',
        'latitude': 61.4915,
        'longitude': -145.5887,
        'depth_in_m': 41000.0,
        'id': str(eid)
})

data = read(path_data, format='sac', 
            event_id=event_id,
            station_id_list=station_id_list,
            tags=['units:m', 'type:velocity']) 


data.sort_by_azimuth()
stations = data.get_stations()

greens = db.get_greens_tensors(stations, origin)

synthetics = greens.get_synthetics(MomentTensor([269000000000000,-15643000000000000,15374000000000000,5999000000000000,7576000000000000,4205000000000000]),components=['Z','R','T'])

os.makedirs(eid+'_synthetics',exist_ok=True)
synthetics.write(eid+'_synthetics/%s/' % model,format='sac') 

As I said, this works when I am using SPECFEM-generated GFs, and the only thing that changes is how I open the GF database:

SF_path = './GFs/%s/%s' % (eid,model)
db = open_db(SF_path,format='SPECFEM3D')

Is there something that needs to change about my approach to get synthetics from an FK database?

Thanks,
Amanda

@rmodrak
Copy link
Member

rmodrak commented Oct 24, 2024 via email

@rmodrak
Copy link
Member

rmodrak commented Oct 24, 2024

It looks like the network, station, and location codes are not getting copied to the synthetics. The result is the following output file names.
...R.sac
...T.sac
...Z.sac

@rmodrak
Copy link
Member

rmodrak commented Oct 24, 2024

One way to get the correct filenames is to manually pass a list of ObsPy stats objects via the stats keyword argument when calling get_synthetics().

The syntax would be something like

synthetics = greens.get_synthetics(mt, stats=stations, components=components, mode='map')

This gets complicated though, I would have to look closely in order to recall the exact correct syntax.

Of course, we want these station codes to get passed seamlessly... trying to figure out a larger fix

@ammcpherson
Copy link
Member Author

Thanks for your quick reply, @rmodrak!

I am having a little luck with a small code block along the lines of:

stats = {}
count = 0
for sta in stations:
    stats[count] = {}
    stats[count]['network'] = sta['network']
    stats[count]['station'] = sta['station']
    stats[count]['location'] = 'S3'
    for ii in range(3):
                if ii == 0:
                    stats[count]['channel'] = 'Z'
                elif ii == 1:
                    count += 1
                    stats[count] = stats[count-1]
                    stats[count]['channel'] = 'R'
                elif ii == 2:
                    count += 1
                    stats[count] = stats[count-1]
                    stats[count]['channel'] = 'T'     
    count += 1

In which I now get 1 file written out (there should be closer to 200 files) with the naming convention 'AK.PAX.S3.T.sac'. The 'S3' comes from the naming convention written out from the SPECFEM GF generated synthetics. This only works with mode='apply' in get_synthetics, though. When I try to use mode='map' I seem to encounter issues in _allocate_stream, where the function is making some stats object but also attempting to use the one I have passed to get_synthetics.

@rmodrak
Copy link
Member

rmodrak commented Oct 25, 2024

Very helpful troubleshooting, thanks.

If you needed a quick workaround, maybe try modifying the write() method:

def write(self, path, format='sac'):

It should be possible to pass the stations list as an argument to the write() method, then use that to build the filenames.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants