[yt-users] cut_region bafflement

Kacper Kowalik xarthisius.kk at gmail.com
Thu Oct 8 09:26:08 PDT 2015


On 10/07/2015 07:24 PM, Stephanie Tonnesen wrote:
> Hi Nathan and all,
> 
> I have pasted a code that shows my problem using the isolated galaxy
> output.  You will have to change the line of code that reads in the
> isolated galaxy output, but otherwise it should run for anyone.  See the
> madness yourself!!!
> 
> http://paste.yt-project.org/show/5943/
> 
> Thanks again,
> Stephanie

Hi Stephanie,
while we investigate possible issue with cut_region chaining could you
try using single 'cut_region' call? You can create single mask using
binary operators e.g.:

r4 = ds.cut_region(["
   (obj['Metal_Density']/obj['density'] > 0.33) &
   (obj['x'].in_units('kpc') < 2.0) &
   (obj['z'].in_units('kpc') >= 0.0)
"])

Cheers,
Kacper

> --
> Dr. Stephanie Tonnesen
> Alvin E. Nashman Postdoctoral Fellow
> Carnegie Observatories, Pasadena, CA
> stonnes at gmail.com
> 
> On Wed, Oct 7, 2015 at 4:33 PM, Stephanie Tonnesen <stonnes at gmail.com>
> wrote:
> 
>> Hi Nathan,
>>
>> Thanks for the quick response!  I will see if I can reproduce this on a
>> test dataset.
>>
>> Different is the confusing part to me.  These all should be selecting
>> exactly the same region.  If P3 and P4 are trying to select a region from
>> ds.all_data() instead of ds.disk() then all three should be different.
>> Another strange thing, if I select x > 65 kpc (only the top half of the
>> cylinder, then the volume is 3.7384e-4, half of the values that don't seem
>> right to me!
>>
>> If I do
>> r4c = ds.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>> then the volume doesnt change, is 7.4768e-4
>>
>> If I do:
>> r4 = ds.all_data()
>> r4c = r4.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>> then the volume is, as it should be, 0.3
>>
>> *I think that something going wrong when I have more than one cut_region
>> command.*  If I do:
>>
>> r3 = ds.disk([0.5,0.5,0.5],[-1,0,0],(4,"kpc"),(2,"kpc"))
>> r3c = r3.cut_region(["(obj['x'].in_units('kpc') > 60.0)"])
>>  or
>> r2 = ds.disk([0.5,0.5,0.5],[-1,0,0],(4,"kpc"),(2,"kpc"))
>> r2c = r2.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>>
>> then the volume is 9.346e-5, as it should be.
>>
>> But when I put those together:
>>
>> r3 = ds.disk([0.5,0.5,0.5],[-1,0,0],(4,"kpc"),(2,"kpc"))
>> r3r = r3.cut_region(["(obj['x'].in_units('kpc') > 60.0)"])
>> r3c = r3r.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>>
>> the volume is 7.4768e-4 again!
>>
>> And when I do this:
>> r4 = ds.all_data()
>> r4r = r4.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>> r4c = r4r.cut_region(["(obj['z'].in_units('kpc') > 0.0)"])
>>
>> I get 2.4467, which is more than my box size.
>>
>> --
>> Dr. Stephanie Tonnesen
>> Alvin E. Nashman Postdoctoral Fellow
>> Carnegie Observatories, Pasadena, CA
>> stonnes at gmail.com
>>
>> On Wed, Oct 7, 2015 at 3:30 PM, Nathan Goldbaum <nathan12343 at gmail.com>
>> wrote:
>>
>>> Hi Stephanie,
>>>
>>> I'm not sure off-hand why this is happening.
>>>
>>> It may help one of us to be able to reproduce and understand what you're
>>> seeing if you can phrase your question in the form of an example that one
>>> of can run locally. This would make use of one of the test datasets on
>>> yt-project.org/data, and a full script that we can run locally pasted to
>>> paste.yt-project.org/data
>>>
>>> Finally, just to clarify, the confusion is that the volumes of P3 and P4
>>> are *bigger* than P2 and identical?
>>>
>>> It might be the case that the derived quantity infrastructure is somehow
>>> selecting the data as if it were selecting from ds.all_data() instead of
>>> your disk selector.  Can you try computing what the total volume of cells
>>> in a data object constructed by applying cut_regions on ds.all_data()
>>> instead of ds.disk()?
>>>
>>> -Nathan
>>>
>>> On Wed, Oct 7, 2015 at 5:21 PM, Stephanie Tonnesen <stonnes at gmail.com>
>>> wrote:
>>>
>>>> Dear yt-users,
>>>>
>>>> I am running
>>>>
>>>> Version = 3.2-dev
>>>> Changeset = dc2467c4eae7 (yt) tip
>>>>
>>>> I have been fiddling around with cut_region in a disk region and have
>>>> gotten really confused.  I was initially trying to see whether defining a
>>>> disk region such that the cylinder height was on the x-axis changed how I
>>>> need to make a cut in the cylinder height.  Anyway, as I was fiddling with
>>>> these cuts I found that if I made a "spatial" cut (rather than a different
>>>> variable cut) I was getting different volumes.  Let me show you:
>>>>
>>>> small bit of background:  the whole box is 130 kpc across, so 0.5 0.5
>>>> 0.5 is 65 65 65 kpc.  Also, I am cutting and pasting bits here so the
>>>> indentation is off.
>>>>
>>>> ds =
>>>> yt.load("/data/stonnes/dwarf_nowind_lowmassgal_bigtracer/DD"+loop[i]+"/DD44sMBhn25now"+loop[i])
>>>>
>>>> r2 = ds.disk([0.5,0.5,0.5],[-1,0,0],(4,"kpc"),(2,"kpc"))
>>>>     r2c = r2.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>>>>
>>>>     r3 = ds.disk([0.5,0.5,0.5],[-1,0,0],(4,"kpc"),(2,"kpc"))
>>>>     r3r = r3.cut_region(["(obj['x'].in_units('kpc') > 63.0)"])
>>>>     r3r = r3r.cut_region(["(obj['x'].in_units('kpc') < 67.0)"])
>>>>     r3c = r3r.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>>>>
>>>>     r4 = ds.disk([0.5,0.5,0.5],[-1,0,0],(4,"kpc"),(2,"kpc"))
>>>>     r4r = r4.cut_region(["(obj['x'].in_units('kpc') > 61.)"])
>>>>     r4c = r4r.cut_region(["obj['Metal_Density']/obj['density'] > 0.33"])
>>>>
>>>>  P2 =
>>>> np.append(P2,r2c.quantities.total_quantity([("index","cell_volume")]))
>>>>  P3 =
>>>> np.append(P3,r3c.quantities.total_quantity([("index","cell_volume")]))
>>>>  P4 =
>>>> np.append(P4,r4c.quantities.total_quantity([("index","cell_volume")]))
>>>>
>>>> P2 volume is 9.346e-5 : within 2% of my calculation for the volume I
>>>> should be getting
>>>> P3 volume is 7.4768e-4
>>>> P4 volume is 7.4768e-4
>>>>
>>>> When I print out the x- y- and z- maxes and mins, they are the same for
>>>> each of the selected regions.
>>>> (By doing this:
>>>> m2 = np.append(m2,r2c.quantities.extrema([("index","z")])[0])
>>>> rho2 = np.append(rho2,r2c.quantities.extrema([("index","z")])[1])
>>>>
>>>> )
>>>>
>>>> Can anyone help me figure out/tell me what is going on here?
>>>>
>>>> Thanks!
>>>> Stephanie
>>>>
>>>> --
>>>> Dr. Stephanie Tonnesen
>>>> Alvin E. Nashman Postdoctoral Fellow
>>>> Carnegie Observatories, Pasadena, CA
>>>> stonnes at gmail.com
>>>>
>>>> _______________________________________________
>>>> yt-users mailing list
>>>> yt-users at lists.spacepope.org
>>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>>
>>>>
>>>
>>> _______________________________________________
>>> yt-users mailing list
>>> yt-users at lists.spacepope.org
>>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>>
>>>
>>
> 
> 
> 
> _______________________________________________
> yt-users mailing list
> yt-users at lists.spacepope.org
> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
> 


-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 819 bytes
Desc: OpenPGP digital signature
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20151008/ad506efd/attachment.sig>


More information about the yt-users mailing list