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

First commit of Density test #140

Merged
merged 7 commits into from
Aug 22, 2018
Merged

Conversation

fjaviersanchez
Copy link
Contributor

This PR tries to cover one of the tests that are in #123 and solves #134. I am selecting the galaxies from the DM object catalogs using their extendedness values (in particular I require extendedness==1). I have an example here: https://portal.nersc.gov/project/lsst/descqa/v2/?run=2018-08-14_10&test=DensityVersusExtinction&right=2018-08/2018-08-14_10/DensityVersusExtinction/dc2_coadd_run1.1p/dc2_coadd_run1.1p_density_vs_extinction.png

@rmjarvis Right now everything is very vanilla and I assumed that there's no depth variation but no cuts have been made to ensure this. Should I include something to do this?

Copy link
Member

@yymao yymao left a comment

Choose a reason for hiding this comment

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

@fjaviersanchez thanks for the PR. I made some minor suggestions/comments.


catalog_data = catalog_instance.get_quantities(['ra', 'dec', 'extendedness'])
galaxies = catalog_data['extendedness']==1
data_map = create_hp_map(catalog_data['ra'][galaxies], catalog_data['dec'][galaxies], self.nside)
Copy link
Member

Choose a reason for hiding this comment

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

Lines 34-36: (1) add has_quantities and skip test if not all needed quantities exist, (2) use filters instead of masking the data afterwards.

        if not catalog_instance.has_quantities(['ra', 'dec', 'extendedness']):
            return TestResult(skipped=True, summary='catalog does not have needed quantities')
        catalog_data = catalog_instance.get_quantities(['ra', 'dec', 'extendedness'], filters=['extendedness == 1'])
        data_map = create_hp_map(catalog_data['ra'], catalog_data['dec'], self.nside)

mask = data_map>0 # This is a good approximation if the pixels are big enough
xmin = np.percentile(self.validation_data[mask], 5)
xmax = np.percentile(self.validation_data[mask], 95)
data_map = data_map/(3600*hp.nside2pixarea(self.nside, degrees=True)) # To get the density in arcmin^-2
Copy link
Member

Choose a reason for hiding this comment

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

suggestion for Line 40

data_map /= (3600*hp.nside2pixarea(self.nside, degrees=True))

data_map = create_hp_map(catalog_data['ra'][galaxies], catalog_data['dec'][galaxies], self.nside)
mask = data_map>0 # This is a good approximation if the pixels are big enough
xmin = np.percentile(self.validation_data[mask], 5)
xmax = np.percentile(self.validation_data[mask], 95)
Copy link
Member

Choose a reason for hiding this comment

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

suggestion for Lines 38-39

xmin, xmax = np.percentile(self.validation_data[mask], [5, 95])

Auxiliary function to generate HEALPix maps from catalogs.
It reads the ra and dec in degrees and returns a HEALPix map
"""
pixnums = hp.ang2pix(nside, np.pi/2-np.radians(dec), np.radians(ra))
Copy link
Member

Choose a reason for hiding this comment

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

    pixnums = hp.ang2pix(nside, ra, dec, lonlat=True)

pixnums = hp.ang2pix(nside, np.pi/2-np.radians(dec), np.radians(ra))
return np.bincount(pixnums, minlength=hp.nside2npix(nside))

class DensityVersusExtinction(BaseValidationTest):
Copy link
Member

Choose a reason for hiding this comment

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

Will this class better be named as DensityVersusSkyPosition or something along the line? In particular, since validation_map_filename and xlabel are both configurable options, one can supply another validation map that is not related to extinction, right? Then it seems to me that the class name is too restrictive?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

You are totally right. Thanks!

Copy link
Member

@yymao yymao left a comment

Choose a reason for hiding this comment

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

@fjaviersanchez, thanks for the update. One more minor change!

pixnums = hp.ang2pix(nside, np.pi/2-np.radians(dec), np.radians(ra))
return np.bincount(pixnums, minlength=hp.nside2npix(nside))
pixnums = hp.ang2pix(nside, ra, dec, lonlat=True)
return 1.0*np.bincount(pixnums, minlength=hp.nside2npix(nside))
Copy link
Member

Choose a reason for hiding this comment

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

to make it more explicit about the intend

    return np.bincount(pixnums, minlength=hp.nside2npix(nside)).astype(np.float)

yymao
yymao previously approved these changes Aug 16, 2018
Copy link
Member

@yymao yymao left a comment

Choose a reason for hiding this comment

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

thanks, @fjaviersanchez. Should I wait for @rmjarvis to comment before merging?

@fjaviersanchez
Copy link
Contributor Author

Thank you @yymao. Yes, let's wait for @rmjarvis to see whether he has some suggestions for additions/refinements please.

rmjarvis
rmjarvis previously approved these changes Aug 16, 2018
Copy link

@rmjarvis rmjarvis left a comment

Choose a reason for hiding this comment

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

Looks fine. Just a couple comments.


class DensityVersusSkyPosition(BaseValidationTest):
"""
This tests checks the object density as a function

Choose a reason for hiding this comment

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

typo: test

This tests checks the object density as a function
of another map (e.g: extinction, airmass, etc)
"""
def __init__(self,**kwargs): # pylint: disable=W0231

Choose a reason for hiding this comment

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

Why use **kwargs here rather than list the kwargs? It looks like the are all required, so normally that would be good to list explicitly in the function declaration. Also worth mentioning what the parameters are in the doc string.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I did this to be able to "recycle" this script with a different configuration file (and to match the style of most of the tests in descqa). @yymao is there an easy way to use positional arguments and be able to read from different config files?

Copy link
Member

Choose a reason for hiding this comment

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

@fjaviersanchez I'm not quite sure what you mean "use positional arguments and be able to read from different config files." It is OK to use positional arguments in DESCQA tests; the only requirement is that you must allow additional keyword arguments at the end. See this test for an example. That being said, since these test classes are rarely being called by users in interactive environments, I'm ok with just using **kwargs. The more important thing is to specify the arguments in docstring which I see you already did.

@fjaviersanchez fjaviersanchez dismissed stale reviews from rmjarvis and yymao via 209805f August 22, 2018 02:21
Copy link
Member

@yymao yymao left a comment

Choose a reason for hiding this comment

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

@fjaviersanchez thanks for the updates. One more small issue. The class is now called DensityVersusSkyPosition as it can accept different validation data, but I think the config yaml file should be named something like DensityVersusExtinction.yaml because this particular config is calling the DensityVersusSkyPosition class to plot Density vs. Extinction, and hence the config itself refers to a specialized test. Does this make sense?

@yymao yymao merged commit c189758 into LSSTDESC:master Aug 22, 2018
@yymao yymao mentioned this pull request Aug 23, 2018
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

3 participants