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

Qbo wavelet #836

Merged
merged 14 commits into from
Sep 17, 2024
Merged

Qbo wavelet #836

merged 14 commits into from
Sep 17, 2024

Conversation

justin-richling
Copy link
Collaborator

@justin-richling justin-richling commented Aug 22, 2024

Description

Add additional QBO diagnostics for power spectral density wavelet calculations and plot, see comment from #816 (comment)

Checklist

  • My code follows the style guidelines of this project
  • I have performed a self-review of my own code
  • My changes generate no new warnings
  • Any dependent changes have been merged and published in downstream modules

If applicable:

  • New and existing unit tests pass with my changes (locally and CI/CD build)
  • I have added tests that prove my fix is effective or that my feature works
  • I have commented my code, particularly in hard-to-understand areas
  • I have made corresponding changes to the documentation
  • I have noted that this is a breaking change for a major release (fix or feature that would cause existing functionality to not work as expected)

Comment on lines 161 to 169
def get_psd_from_wavelet(data):
deg = 6
period = np.arange(1,55+1)
freq = 1/period
widths = deg / (2*np.pi*freq)
cwtmatr = scipy.signal.cwt( data, scipy.signal.morlet2, widths=widths, w=deg )
psd = np.mean(np.square(np.abs(cwtmatr)),axis=1)
return ( period, psd )

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not sure what the syntax guidance is for the diags - but I think it would be good to add a comment section for these routines. This one especially seems like it would benefit from a description. Something like this:

def get_psd_from_wavelet(data):
   '''
   Return power spectral density using a complex Morlet wavelet spectrum of degree 6
   '''

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@whannah1 thanks for the review. I agree that inline documentation about the new function will be very helpful!

@chengzhuzhang
Copy link
Contributor

Thank you for this PR @justin-richling and @whannah1
I will review and test soon!

@chengzhuzhang
Copy link
Contributor

It looks like I need to update the json.dump call to bypass the TypeError: Object of type int64 is not JSON serializable:

  •            json.dump(json_dict, outfile)
    
  •            json.dump(json_dict, outfile, default=str)
    

Also the 4th plot I got seems problematic:
image
@justin-richling could you share the ERA data you are testing with? This looks like a data problem.. Thank you!

@justin-richling
Copy link
Collaborator Author

@chengzhuzhang Thanks for pointing this out, the code was getting the power spectra level based on the test case only, not the reference case so I updated this to the qbo_driver.py script.

* Add vertical line to indicate wavelet spectra peak and x-axis label for corresponding periods
* Use square root values for clarity
@justin-richling
Copy link
Collaborator Author

justin-richling commented Sep 10, 2024

@chengzhuzhang @tomvothecoder After discussion with @whannah1 it seems for now we will not be including the external wavelet analysis package until a further date. However, I added some clean up and small enhancements to the plot. Here is a snippet of the QBO diags with updated wavet plot:

qbo_diags (2)

With that, this PR should be good to review again, but I seem to be failing the auto checks. I can't seem to narrow down why these small changes are failing but I'm happy to fix them if there easy solutions.

@chengzhuzhang
Copy link
Contributor

chengzhuzhang commented Sep 12, 2024

@justin-richling and @whannah1 thank you for the nice enhancement to identify the value of the peak period. I think we can do something similar to the original method. One thing I would suggest is that the x-axis is also updated and, now the wavelet method has different ticks compare to the original method, as shown below:
Screenshot 2024-09-12 at 2 20 36 PM

Can we consider to revert this part of change to align x-axis for both plots. This makes it easier to compare both methods visually. Let me know what you think...Thanks!

P.S. don't worry about the auto checks, I can fixes those formatting issues after we finalize the figures.

@whannah1
Copy link

Can we consider to revert this part of change to align x-axis for both plots. This makes it easier to compare both methods visually.

That seems like a good idea.

@chengzhuzhang
Copy link
Contributor

I slightly tweaked the plotting code, so that 3rd and 4th panels have same axis range, and with period of maximum power/variance labeled. The figure is shown below. @justin-richling @whannah1 and @jjbenedict let me know if you'd like any changes. Many thanks for this nice enhancement!

image

@jjbenedict
Copy link
Collaborator

@chengzhuzhang It's great to see the peak align so closely, thanks for adding this to E3SM Diags! My only comment: I think there's a typo in the wavelet y-axis label. And, variance should be m2 s-2...?

@chengzhuzhang
Copy link
Contributor

@jjbenedict thank you for catching the error. The units are now fixed with 3fd8aa5

@justin-richling the auto checks are also passed now with some formatting and after merging fixes from main. Appreciate your hard work with e3sm_diags!

@chengzhuzhang
Copy link
Contributor

@tomvothecoder as a heads-up. We will need to bring this PR to the cdat-migration branch.

@chengzhuzhang chengzhuzhang merged commit 2b303df into E3SM-Project:main Sep 17, 2024
4 checks passed
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

Successfully merging this pull request may close these issues.

4 participants