Using python library Rasterio to create a subset of a TIFF image and then display it and save it?












2















I have two Raster images, one from band 4 with a B4 at the end and another from band 5 with B5 at the end. I want to subset the B5 raster to 800x600, then display it and save it as a GeoTiff. Then I want to compute the NDVI (I assume I'll need both the B4 and B5 to do this, but not sure). Then I want to display the NDVI subset of the B5 raster. Display it and save it as a GeoTiff.



How would I create something like a 800 x 600 pixel subset of a TIFF raster image? I want to also take that TIFF and generate an NDVI image for that subset.



NOTE: I am working with a Landsat image. The image has B5 at the end of the file title.



What I've done so far:



import rasterio
from rasterio.windows import Window
import matplotlib.pyplot plt # for later use

with rasterio.open('MyRasterImage.tif') as src:
w = src.read(1, window=Window(0, 0, 800, 600))


I want to display this using Spyder or Jupyter notebooks. So I thought to use matplotlib and did the flowing code:



# Plot
plt.imshow(w)
plt.show()


Doing this generates a 800x600 matplotlib window, but it's all purple, not sure why its producing this.



Now I want to be able to display this 800x600 image. Then after that I want preform an NDVI on that subset 800x600 image. Then display the subset 800x600 image with NDVI showing.



I know the formuala is: NDVI = (NIR - red) / (NIR + red)



But how do I extract out NIR and red here from this single Landsat image?



My attempt at that:



band1 = dataset.read(1)
band2 = dataset.read(2)
band3 = dataset.read(3)
print(band[2])


When I run that code for the bands I get the error:



rasterio indexerror: band index 2 out of range (not in (1,))


When I run this code:



print(w.count)


It returns '1'.



So this means that the Landsat image only has one band? But in order to do NDVI don't I need 3 bands?



I am thinking of writing some code like this to get the NDVI from that raster. But not sure how to go about extracting out the bands:



# We handle the connections with "with"
with rasterio.open(bands[0]) as src:
b3 = src.read(1)

with rasterio.open(bands[1]) as src:
b4 = src.read(1)

# Allow division by zero
numpy.seterr(divide='ignore', invalid='ignore')

# Calculate NDVI
ndvi = (b4.astype(float) - b3.astype(float)) / (b4 + b3)


This code doesn't work because bands isn't define as anything so I don't know how to define bands to get the NDVI.



After this I am not sure how to go about both displaying and saving the rendered image.










share|improve this question





























    2















    I have two Raster images, one from band 4 with a B4 at the end and another from band 5 with B5 at the end. I want to subset the B5 raster to 800x600, then display it and save it as a GeoTiff. Then I want to compute the NDVI (I assume I'll need both the B4 and B5 to do this, but not sure). Then I want to display the NDVI subset of the B5 raster. Display it and save it as a GeoTiff.



    How would I create something like a 800 x 600 pixel subset of a TIFF raster image? I want to also take that TIFF and generate an NDVI image for that subset.



    NOTE: I am working with a Landsat image. The image has B5 at the end of the file title.



    What I've done so far:



    import rasterio
    from rasterio.windows import Window
    import matplotlib.pyplot plt # for later use

    with rasterio.open('MyRasterImage.tif') as src:
    w = src.read(1, window=Window(0, 0, 800, 600))


    I want to display this using Spyder or Jupyter notebooks. So I thought to use matplotlib and did the flowing code:



    # Plot
    plt.imshow(w)
    plt.show()


    Doing this generates a 800x600 matplotlib window, but it's all purple, not sure why its producing this.



    Now I want to be able to display this 800x600 image. Then after that I want preform an NDVI on that subset 800x600 image. Then display the subset 800x600 image with NDVI showing.



    I know the formuala is: NDVI = (NIR - red) / (NIR + red)



    But how do I extract out NIR and red here from this single Landsat image?



    My attempt at that:



    band1 = dataset.read(1)
    band2 = dataset.read(2)
    band3 = dataset.read(3)
    print(band[2])


    When I run that code for the bands I get the error:



    rasterio indexerror: band index 2 out of range (not in (1,))


    When I run this code:



    print(w.count)


    It returns '1'.



    So this means that the Landsat image only has one band? But in order to do NDVI don't I need 3 bands?



    I am thinking of writing some code like this to get the NDVI from that raster. But not sure how to go about extracting out the bands:



    # We handle the connections with "with"
    with rasterio.open(bands[0]) as src:
    b3 = src.read(1)

    with rasterio.open(bands[1]) as src:
    b4 = src.read(1)

    # Allow division by zero
    numpy.seterr(divide='ignore', invalid='ignore')

    # Calculate NDVI
    ndvi = (b4.astype(float) - b3.astype(float)) / (b4 + b3)


    This code doesn't work because bands isn't define as anything so I don't know how to define bands to get the NDVI.



    After this I am not sure how to go about both displaying and saving the rendered image.










    share|improve this question



























      2












      2








      2








      I have two Raster images, one from band 4 with a B4 at the end and another from band 5 with B5 at the end. I want to subset the B5 raster to 800x600, then display it and save it as a GeoTiff. Then I want to compute the NDVI (I assume I'll need both the B4 and B5 to do this, but not sure). Then I want to display the NDVI subset of the B5 raster. Display it and save it as a GeoTiff.



      How would I create something like a 800 x 600 pixel subset of a TIFF raster image? I want to also take that TIFF and generate an NDVI image for that subset.



      NOTE: I am working with a Landsat image. The image has B5 at the end of the file title.



      What I've done so far:



      import rasterio
      from rasterio.windows import Window
      import matplotlib.pyplot plt # for later use

      with rasterio.open('MyRasterImage.tif') as src:
      w = src.read(1, window=Window(0, 0, 800, 600))


      I want to display this using Spyder or Jupyter notebooks. So I thought to use matplotlib and did the flowing code:



      # Plot
      plt.imshow(w)
      plt.show()


      Doing this generates a 800x600 matplotlib window, but it's all purple, not sure why its producing this.



      Now I want to be able to display this 800x600 image. Then after that I want preform an NDVI on that subset 800x600 image. Then display the subset 800x600 image with NDVI showing.



      I know the formuala is: NDVI = (NIR - red) / (NIR + red)



      But how do I extract out NIR and red here from this single Landsat image?



      My attempt at that:



      band1 = dataset.read(1)
      band2 = dataset.read(2)
      band3 = dataset.read(3)
      print(band[2])


      When I run that code for the bands I get the error:



      rasterio indexerror: band index 2 out of range (not in (1,))


      When I run this code:



      print(w.count)


      It returns '1'.



      So this means that the Landsat image only has one band? But in order to do NDVI don't I need 3 bands?



      I am thinking of writing some code like this to get the NDVI from that raster. But not sure how to go about extracting out the bands:



      # We handle the connections with "with"
      with rasterio.open(bands[0]) as src:
      b3 = src.read(1)

      with rasterio.open(bands[1]) as src:
      b4 = src.read(1)

      # Allow division by zero
      numpy.seterr(divide='ignore', invalid='ignore')

      # Calculate NDVI
      ndvi = (b4.astype(float) - b3.astype(float)) / (b4 + b3)


      This code doesn't work because bands isn't define as anything so I don't know how to define bands to get the NDVI.



      After this I am not sure how to go about both displaying and saving the rendered image.










      share|improve this question
















      I have two Raster images, one from band 4 with a B4 at the end and another from band 5 with B5 at the end. I want to subset the B5 raster to 800x600, then display it and save it as a GeoTiff. Then I want to compute the NDVI (I assume I'll need both the B4 and B5 to do this, but not sure). Then I want to display the NDVI subset of the B5 raster. Display it and save it as a GeoTiff.



      How would I create something like a 800 x 600 pixel subset of a TIFF raster image? I want to also take that TIFF and generate an NDVI image for that subset.



      NOTE: I am working with a Landsat image. The image has B5 at the end of the file title.



      What I've done so far:



      import rasterio
      from rasterio.windows import Window
      import matplotlib.pyplot plt # for later use

      with rasterio.open('MyRasterImage.tif') as src:
      w = src.read(1, window=Window(0, 0, 800, 600))


      I want to display this using Spyder or Jupyter notebooks. So I thought to use matplotlib and did the flowing code:



      # Plot
      plt.imshow(w)
      plt.show()


      Doing this generates a 800x600 matplotlib window, but it's all purple, not sure why its producing this.



      Now I want to be able to display this 800x600 image. Then after that I want preform an NDVI on that subset 800x600 image. Then display the subset 800x600 image with NDVI showing.



      I know the formuala is: NDVI = (NIR - red) / (NIR + red)



      But how do I extract out NIR and red here from this single Landsat image?



      My attempt at that:



      band1 = dataset.read(1)
      band2 = dataset.read(2)
      band3 = dataset.read(3)
      print(band[2])


      When I run that code for the bands I get the error:



      rasterio indexerror: band index 2 out of range (not in (1,))


      When I run this code:



      print(w.count)


      It returns '1'.



      So this means that the Landsat image only has one band? But in order to do NDVI don't I need 3 bands?



      I am thinking of writing some code like this to get the NDVI from that raster. But not sure how to go about extracting out the bands:



      # We handle the connections with "with"
      with rasterio.open(bands[0]) as src:
      b3 = src.read(1)

      with rasterio.open(bands[1]) as src:
      b4 = src.read(1)

      # Allow division by zero
      numpy.seterr(divide='ignore', invalid='ignore')

      # Calculate NDVI
      ndvi = (b4.astype(float) - b3.astype(float)) / (b4 + b3)


      This code doesn't work because bands isn't define as anything so I don't know how to define bands to get the NDVI.



      After this I am not sure how to go about both displaying and saving the rendered image.







      python rasterio






      share|improve this question















      share|improve this question













      share|improve this question




      share|improve this question








      edited Nov 18 '18 at 16:26







      yuen2

















      asked Nov 18 '18 at 16:16









      yuen2yuen2

      194




      194
























          0






          active

          oldest

          votes











          Your Answer






          StackExchange.ifUsing("editor", function () {
          StackExchange.using("externalEditor", function () {
          StackExchange.using("snippets", function () {
          StackExchange.snippets.init();
          });
          });
          }, "code-snippets");

          StackExchange.ready(function() {
          var channelOptions = {
          tags: "".split(" "),
          id: "1"
          };
          initTagRenderer("".split(" "), "".split(" "), channelOptions);

          StackExchange.using("externalEditor", function() {
          // Have to fire editor after snippets, if snippets enabled
          if (StackExchange.settings.snippets.snippetsEnabled) {
          StackExchange.using("snippets", function() {
          createEditor();
          });
          }
          else {
          createEditor();
          }
          });

          function createEditor() {
          StackExchange.prepareEditor({
          heartbeatType: 'answer',
          autoActivateHeartbeat: false,
          convertImagesToLinks: true,
          noModals: true,
          showLowRepImageUploadWarning: true,
          reputationToPostImages: 10,
          bindNavPrevention: true,
          postfix: "",
          imageUploader: {
          brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
          contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
          allowUrls: true
          },
          onDemand: true,
          discardSelector: ".discard-answer"
          ,immediatelyShowMarkdownHelp:true
          });


          }
          });














          draft saved

          draft discarded


















          StackExchange.ready(
          function () {
          StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fstackoverflow.com%2fquestions%2f53362947%2fusing-python-library-rasterio-to-create-a-subset-of-a-tiff-image-and-then-displa%23new-answer', 'question_page');
          }
          );

          Post as a guest















          Required, but never shown

























          0






          active

          oldest

          votes








          0






          active

          oldest

          votes









          active

          oldest

          votes






          active

          oldest

          votes
















          draft saved

          draft discarded




















































          Thanks for contributing an answer to Stack Overflow!


          • Please be sure to answer the question. Provide details and share your research!

          But avoid



          • Asking for help, clarification, or responding to other answers.

          • Making statements based on opinion; back them up with references or personal experience.


          To learn more, see our tips on writing great answers.




          draft saved


          draft discarded














          StackExchange.ready(
          function () {
          StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fstackoverflow.com%2fquestions%2f53362947%2fusing-python-library-rasterio-to-create-a-subset-of-a-tiff-image-and-then-displa%23new-answer', 'question_page');
          }
          );

          Post as a guest















          Required, but never shown





















































          Required, but never shown














          Required, but never shown












          Required, but never shown







          Required, but never shown

































          Required, but never shown














          Required, but never shown












          Required, but never shown







          Required, but never shown







          這個網誌中的熱門文章

          Xamarin.form Move up view when keyboard appear

          Post-Redirect-Get with Spring WebFlux and Thymeleaf

          Anylogic : not able to use stopDelay()