How to creat bands difference image from MODIS using satpy

[Frank]

Dear All

Please give your's suggetions for how to creat bands difference image from MODIS using satpy. I got the following errors

>>> from pyhdf.SD import SD, SDC

>>>  

... import numpy as np

>>> import matplotlib.pyplot as plt

>>> import matplotlib as mpl

>>> import matplotlib.cm as cm

>>> import pprint

>>> from pyhdf.SD import SD, SDC

>>> from scipy.misc import bytescale

>>> import numpy as np

>>> import warnings

>>> from pyresample import load_area, save_quicklook, SwathDefinition

>>> from pyresample.kd_tree import resample_nearest

>>> from satpy import Scene, MultiScene

>>> import numpy as np

>>> import pandas as pd

>>> import matplotlib.pyplot as plt

>>> from pyresample import geometry

>>> from satpy import Scene

>>> from glob import glob

>>> from satpy import DatasetID

>>> import matplotlib.pyplot as plt

>>> import cartopy.crs as ccrs

>>> from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter

>>> from matplotlib.offsetbox import OffsetImage, AnnotationBbox

>>> import cartopy.crs as ccrs

>>> from PIL import Image

>>> from matplotlib.cbook import get_sample_data

>>> import cartopy.io.shapereader as shpreader

>>> import cartopy.feature as cf

>>> import cartopy.crs as ccrs

>>> from contextlib import suppress

>>> from mpl_toolkits.basemap import Basemap

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

  File "/home/vankayk/anaconda3/lib/python3.6/site-packages/mpl_toolkits/basemap/__init__.py", line 155, in <module>

    pyproj_datadir = os.environ['PROJ_LIB']

  File "/home/vankayk/anaconda3/lib/python3.6/os.py", line 669, in __getitem__

    raise KeyError(key) from None

KeyError: 'PROJ_LIB'

>>> from pyresample.geometry import AreaDefinition

>>> from satpy.config import check_satpy

>>> from satpy.writers import to_image

>>> from satpy.composites import PaletteCompositor

>>> import numpy as np

>>> from satpy.composites import BWCompositor

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

ImportError: cannot import name 'BWCompositor'

>>> from satpy.enhancements import colorize

>>> from satpy.writers import to_image

>>> from satpy import Scene, MultiScene

>>> from satpy.multiscene import timeseries

>>> from glob import glob

>>> from pyresample.geometry import AreaDefinition

>>> my_area = AreaDefinition(...)

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

TypeError: __init__() missing 6 required positional arguments: 'description', 'proj_id', 'projection', 'width', 'height', and 'area_extent'

>>> from satpy import Scene, MultiScene

>>> from glob import glob

>>> import folium

>>> from folium import plugins

>>> from scipy.ndimage import imread

>>> from cartopy.img_transform import warp_array

>>> import cartopy.crs as ccrs

>>>

>>> import folium

>>> from folium import plugins

>>> from scipy.ndimage import imread

>>> from cartopy.img_transform import warp_array

>>> import cartopy.crs as ccrs

>>> import cartopy.io.img_tiles as cimgt

>>> fnames = glob('MYD021KM.A2019286.1025.061.2019287145750.hdf')

>>>

>>>

>>> mod02_scn = Scene(reader='modis_l1b', filenames=fnames)

>>> mod02_scn.load(["3"],resolution=1000)

>>> mod02_scn.load(["31"],resolution=1000)

>>>

>>> mod02_scn1 = mod02_scn["3"] - mod02_scn["31"]

>>> mod02_scn2 = mod02_scn["3"] + mod02_scn["31"]

>>>

>>> mod02_scn3 = mod02_scn1 / mod02_scn2

>>>

>>> area_id = 'seviri_0deg'

>>> description = 'Seviri 0 Degree'

>>> proj_id = 'seviri_0deg'

>>> #x_size = 3712

... #y_size = 3712

... width = 800 # width of the result domain in pixels

>>> height = 800 # height of the result domain in pixels

>>>

>>> llx = 35E5 # projection x coordinate of lower left corner of lower left pixel

>>> lly = 10E5 # projection y coordinate of lower left corner of lower left pixel

>>> urx = 62E5 # projection x coordinate of upper right corner of upper right pixel

>>> ury = 34E5 # projection y coordinate of upper right corner of upper right pixel

>>>

>>> area_extent = (llx,lly,urx,ury)

>>> proj_dict = {'a': 6378169.0, 'b': 6378169.0,'units': 'm', 'lon_0': 0.0,'proj': 'eqc', 'lat_0': 0.0}

>>>

>>> area_def = geometry.AreaDefinition(area_id, description, proj_id, proj_dict, width, height, area_extent)

>>>

>>> local_scn = mod02_scn3.resample(area_def)

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

  File "/home/vankayk/anaconda3/lib/python3.6/site-packages/xarray/core/common.py", line 964, in resample

    indexer = either_dict_or_kwargs(indexer, indexer_kwargs, "resample")

  File "/home/vankayk/anaconda3/lib/python3.6/site-packages/xarray/core/utils.py", line 256, in either_dict_or_kwargs

    "the first argument to .%s must be a dictionary" % func_name

ValueError: the first argument to .resample must be a dictionary

>>> image = np.asarray(local_scn.transpose(0, 1)

... #image = np.transpose(np.matrix(image))

...

... image = np.interp(image, (np.nanpercentile(image,1), np.nanpercentile(image,99)), (0, 1))

  File "<stdin>", line 4

    image = np.interp(image, (np.nanpercentile(image,1), np.nanpercentile(image,99)), (0, 1))

        ^

SyntaxError: invalid syntax

>>> crs = local_scn.attrs["area"].to_cartopy_crs()

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

NameError: name 'local_scn' is not defined

>>> fig = plt.subplots(figsize=(16,12))

>>> fig = plt.subplots(figsize=(16,12))

>>> ax = plt.axes(projection=crs)

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

NameError: name 'crs' is not defined

>>>

>>> #ax.coastlines(resolution='50m', color='hotpink',linewidth=1)

...

>>>

>>> ax.imshow(image, transform=crs, extent=crs.bounds, origin="upper").add_to(map_hooray)

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

NameError: name 'ax' is not defined

>>> ax.add_feature(cartopy.feature.BORDERS, linestyle=':',linewidth=1)

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

NameError: name 'ax' is not defined

>>> ax.coastlines(resolution="10m",color="grey")

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

NameError: name 'ax' is not defined

>>> #ax.add_feature(cartopy.feature.LAND,facecolor='0.7',color="grey")

... ax.add_feature(cartopy.feature.LAND, zorder=12, edgecolor='#7f7f7f', facecolor='#B1B2B4')

Traceback (most recent call last):

  File "<stdin>", line 2, in <module>

NameError: name 'ax' is not defined

>>>

>>> ax.map_.add_children(folium.raster_layers.ImageOverlay(image,mercator_project=True, opacity=0.8,bounds =[[11, 38], [30, 48]]))

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

NameError: name 'ax' is not defined

>>> map_.save('test123.html')

Traceback (most recent call last):

  File "<stdin>", line 1, in <module>

NameError: name 'map_' is not defined

>>> plt.savefig('RAO123.png')

Thanks.

Frank

[David Hoese]

Hi Frank,

Is AI giving you this code to use? Where are you finding examples that have this type of code?

1. You are using Python 3.6. Satpy has not supported Python 3.6 for many years. Please create a new environment with Python 3.12 or newer.

2. I see mentions of basemap. This package is unmaintained and cartopy and matplotlib should be the only matplotlib-based plotting utilities you need.

3. Some of your errors are syntax errors or weird environment library loading errors or not passing arguments to things or using variables that don't exist. If you're following examples or being lead by AI then you need to find a complete example or figure out how to connect the different examples you're being given. This isn't something wrong with Satpy, but with the code you're being given.

As said in my private message to you, your messages seem like AI dumping errors from posts across mailing lists and forums. Please take time to formulate your responses with what you're trying to do, what isn't working, what your environment is, and what your code looks like. I won't be able to dedicate much more time to trying to decipher these messages.

Dave

--
You received this message because you are subscribed to the Google Groups "pytroll" group.
To unsubscribe from this group and stop receiving emails from it, send an email to pytroll+u...@googlegroups.com.
To view this discussion, visit https://groups.google.com/d/msgid/pytroll/c4eb5d21-4bbd-40ee-9d2d-e7e434ab9dccn%40googlegroups.com.

[Frank]

I apologize for that, I'm new to pytroll like I said on the email I sent to you, I will do better next time, I'm still learning pytroll and I use AI for guidance till I found this group