How to calculate mean for each band in raster using GDAL and Python 3

Question

I have a raster with 6 bands and I need to calculate the mean for all bands in the raster without value 0

I have this script but it works for one band only:

import os
from osgeo import gdal
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import ogr 

raster = gdal.Open("D:/script/NDVI2015.tif")
bands = raster.RasterCount
bands
6

geotransform = raster.GetGeoTransform()
proj= raster.GetProjection()
data = raster.GetRasterBand(2).ReadAsArray().astype('float')
meandata(2) = np.mean(data[data != 0])#calcul mean without value 0

I'm very new with Python and GDAL.

Answer 1

Building on @Luke's comment, incorporate the band number into a for loop using range():

import gdal
import numpy as np

raster = gdal.Open(r'/path/to/your/image.tif')
bands = raster.RasterCount

for band in range(1, bands+1):
    data = raster.GetRasterBand(band).ReadAsArray().astype('float')
    mean = np.mean(data[data != 0]) #calculate mean without value 0
    print("Band %s: Mean = %s" % (band, round(mean, 2)))