[yt-users] Problem with making SlicePlot and ProjectionPlot from new derived field
tazkera haque
h.tazkera at gmail.com
Thu Oct 5 17:45:12 PDT 2017
HI Nathan and Bili,
Thanks very much for the help. This makes perfect sense now. I was not
aware of the difference of 'gas' and 'Gas' fields. Also the unit override
function is very helpful in this case.
Best
Tazkera
On Thu, Oct 5, 2017 at 7:40 PM, Bili Dong - Gmail <qobilidop at gmail.com>
wrote:
> Hi Tazkera,
>
> Now I understand where the error comes from. It's not from the unit
> system but from the difference between 'Gas' and 'gas' fields. Let me show
> you the code that works and explain later. The corrected code looks like
> this:
>
> ```
> dd = ds.all_data()
>
> # derive a new domain width with conversion factor
> (6Mpc*current_time/hubble_constant)
> def _NewDomainWidth(field, data):
> return ds.domain_width*((3.086e18*6e6*0.909090893796)/(0.7))
>
> ds.add_field(('gas', 'PhysicalDomainWidth'), function=_NewDomainWidth,
> units="code_length")
>
> # derive a new domain center with conversion factor
> (6Mpc*current_time/hubble_constant)
> def _NewDomainCenter(field, data):
> return ds.domain_center*((3.086e18*6e6*0.909090893796)/(0.7))
>
> ds.add_field(('gas', 'PhysicalDomainCenter'), function=_NewDomainCenter,
> units="code_length")
>
> # derive a new electron pressure field
> def _ElectronPressure(field, data):
> return data['Gas','Density']*(5.49e-7)*data['Gas','Temperature']
>
> ds.add_field(('Gas', 'ElectronPressure'), function=_ElectronPressure,
> units="(code_mass*K)/(code_length**3)",
> particle_type=True)
> f_smooth = ds.add_smoothed_particle_field(('Gas', 'ElectronPressure'))
> # f_smooth equals ('deposit', 'Gas_smoothed_ElectronPressure')
>
>
> # integration of pressure field along line of sight
> yt.ProjectionPlot(ds, 'z', f_smooth, width = dd['gas',
> 'PhysicalDomainWidth'],
> center = dd['gas', 'PhysicalDomainCenter'] )
> ```
>
> First, let's go back to the error message. Where do the magic numbers (39678,)
> (54608,) come from? Print out dd['Gas', 'Density'].shape and dd['gas',
> 'density'].shape and you'll see. The issue with your original code is
> that you tried to create a 'gas' field from a 'Gas' field but their
> dimension doesn't match. Here 'Gas' field is a particle field retrieved
> from your snapshot directly, but 'gas' field is an SPH field calculated
> from the particle field. They have different dimensions. More detailed
> documentation could be found here:
>
> http://yt-project.org/docs/dev/analyzing/fields.html#sph-fields
>
> That's for solving the error. As a side note, I do not recommend defining
> fields like your dd['gas', 'PhysicalDomainCenter'], though it works (to
> my surprise). Because usually we assume the same field type means the
> same dimension (as in the first number of shape). Your case is more
> suitable to define a new variable.
>
> Hope it helps.
>
> Best,
> Bili
>
>
> On Thu, Oct 5, 2017 at 3:54 PM, Bili Dong - Gmail <qobilidop at gmail.com>
> wrote:
>
>>
>> On Thu, Oct 5, 2017 at 3:48 PM, Bili Dong - Gmail <qobilidop at gmail.com>
>> wrote:
>>
>>> Hi Tazkera,
>>>
>>> I think yt should be able to understand your code units and do the unit
>>> conversion for you, or easily configured to interpret your code units. See
>>> an example here
>>> <http://yt-project.org/docs/dev/examining/loading_data.html#tipsy-data>.
>>> This I feel is the better solution.
>>>
>>> But if you want to the conversion by hand, a quick fix to your code is
>>> to define new variables for domain width and domain center:
>>>
>>> ```
>>> new_domain_width = ds.domain_width*((3.086e18*6e6*0.909090893796)/(0.7))
>>> new_domain_center = ds.domain_center*((3.086e18*
>>> 6e6*0.909090893796)/(0.7))
>>> yt.ProjectionPlot(ds, 'z', ('gas', 'ElectronPressure'), width=
>>> new_domain_width[0],
>>> center=new_domain_center)
>>> ```
>>>
>>> The issue with your code is that ds.domain_width and ds.domain_center
>>> are not fields themselves. If you try to add new fields for those constant
>>> values, you'll get a long list of those values as you could read from the
>>> error message. You could print out dd['gas', 'PhysicalDomainWidth'] or dd['gas',
>>> 'PhysicalDomainCenter'] in your code to see what they look like. That
>>> might help you understand what was wrong.
>>>
>>
>> Sorry I was wrong. This is not true.
>>
>>>
>>> Best,
>>> Bili
>>>
>>>
>>> On Thu, Oct 5, 2017 at 2:59 PM, tazkera haque <h.tazkera at gmail.com>
>>> wrote:
>>>
>>>> Dear yt people,
>>>>
>>>> I am working on a tipsy output file from a n-body simulation. My file
>>>> contains some scaled units for fields like density, coordinates, mass etc.
>>>> To make a physical sense of the data, I had to convert those fields units
>>>> to some physical units (c.g.s) and I derived new fields having cgs unit.
>>>> for example:
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> *dd = ds.all_data()# derive a new domain width with conversion factor
>>>> (6Mpc*current_time/hubble_constant) def _NewDomainWidth(field, data):
>>>> return
>>>> ds.domain_width*((3.086e18*6e6*0.909090893796)/(0.7))ds.add_field(('gas',
>>>> 'PhysicalDomainWidth'), function=_NewDomainWidth, units="code_length")*
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> *# derive a new domain center with conversion factor
>>>> (6Mpc*current_time/hubble_constant) def _NewDomainCenter(field, data):
>>>> return
>>>> ds.domain_center*((3.086e18*6e6*0.909090893796)/(0.7))ds.add_field(('gas',
>>>> 'PhysicalDomainCenter'), function=_NewDomainCenter, units="code_length")*
>>>>
>>>> note that here "code_length" is actually in "cm". but I could not just
>>>> do dd["Gas",'Coordinates'].in_units('cm'), because this conversion
>>>> does not take into account some constant corrections and returns the wrong
>>>> length magnitude.
>>>>
>>>> My goal was to make sliceplot and projectionplot of a newly derived
>>>> field :
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> *# derive a new electron pressure fielddef _ElectronPressure(field,
>>>> data): return
>>>> data['Gas','Density']*(5.49e-7)*data['Gas','Temperature']ds.add_field(('gas',
>>>> 'ElectronPressure'), function=_ElectronPressure,
>>>> units="(code_mass*K)/(code_length**3)")*
>>>>
>>>> *# integration of pressure field along line of sight*
>>>>
>>>> *yt.ProjectionPlot(ds, 'z', ('gas', 'ElectronPressure'), width =
>>>> dd['gas', 'PhysicalDomainWidth'], center = dd['gas',
>>>> 'PhysicalDomainCenter'] )*
>>>>
>>>>
>>>> *and it is returning me *
>>>>
>>>> ValueError Traceback (most recent call last)<ipython-input-8-fa351b6c09d1> in <module>() 1 yt.ProjectionPlot(ds, 'z', ('gas', 'ElectronPressure'), width = dd['gas', 'PhysicalDomainWidth'], ----> 2 center = dd['gas', 'PhysicalDomainCenter'] )
>>>>
>>>> ValueError: operands could not be broadcast together with shapes (39678,) (54608,)
>>>>
>>>> Now I cannot use width =ds.domain_width and center = ds.domain_center for the same unit error. I was wondering how I can make the projectionplot for my newly derived field and newly derived domain width and domain center.
>>>>
>>>>
>>>>
>>>> Your suggestion is most welcome
>>>>
>>>>
>>>> Best
>>>>
>>>> Tazkera
>>>>
>>>>
>>>> _______________________________________________
>>>> 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 --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20171005/173aae73/attachment-0001.html>
More information about the yt-users
mailing list