nina_99

Boa tarde, atualizei o link:https://drive.google.com/uc?export=download&id=1wykhadrVHP2I-nqlrfpNd3kQ0ErUoYxh

Obrigada, adicionei o link para baixar o arquivo.

Link do arquivo nc: https://drive.google.com/uc?export=download&id=1wykhadrVHP2I-nqlrfpNd3kQ0ErUoYxh import netCDF4 as nc import numpy as np import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeature # Open NetCDF file nc_file = '/content/drive/MyDrive/rs4/adaptor.mars.internal-1719616919.4771404-4624-16-f31563f3-4d1d-4422-820c-3ca7cd9b2e78.nc' dataset = nc.Dataset(nc_file, 'r') # Extract longitude, latitude, and time variables lon_var = dataset.variables['longitude'][:] lat_var = dataset.variables['latitude'][:] time_var = dataset.variables['time'] mwd_var = dataset.variables['mwd'] # Define the time index to plot time_idx = 1 # Adjust according to your specific time index # Extract the mwd data for the specified time index mwd_data = mwd_var[time_idx, :, :] # Check for missing values missing_value = mwd_var._FillValue if '_FillValue' in mwd_var.ncattrs() else None if missing_value is not None: print(f"Missing value representation: {missing_value}") # Replace missing values with NaN mwd_data = np.where(mwd_data == missing_value, np.nan, mwd_data) # Apply scale factor and add offset if needed scale_factor = mwd_var.scale_factor if 'scale_factor' in mwd_var.ncattrs() else 1.0 add_offset = mwd_var.add_offset if 'add_offset' in mwd_var.ncattrs() else 0.0 mwd_data = mwd_data * scale_factor + add_offset # Print the range of the data to verify print(f"MWD data range: min={np.nanmin(mwd_data)}, max={np.nanmax(mwd_data)}") # Ensure the data values are within a reasonable range (0 to 360 degrees) mwd_data = np.where((mwd_data >= 0) & (mwd_data <= 360), mwd_data, np.nan) # Create U and V components of the wave direction # Assuming MWD (mean wave direction) is in degrees from North mwd_rad = np.deg2rad(mwd_data) U = np.sin(mwd_rad) V = np.cos(mwd_rad) # Plot the wave direction as quivers on a global map plt.figure(figsize=(14, 7), dpi=150) ax = plt.axes(projection=ccrs.PlateCarree()) ax.set_global() ax.coastlines() ax.add_feature(cfeature.BORDERS, linestyle=':') ax.gridlines() # Plot the quivers quiver = ax.quiver(lon_var, lat_var, U, V, mwd_data, transform=ccrs.PlateCarree(), scale=100, cmap='viridis') plt.colorbar(quiver, ax=ax, label='MWD') plt.title(f'Mean Wave Direction (MWD) at time index {time_idx}') plt.xlabel('Longitude') plt.ylabel('Latitude') plt.show() # Close the NetCDF dataset dataset.close()

Boa tarde, a velocidade deveria ser um so resultado, no exemplo acima esta repetindo o mesmo valor 5 vezes dentro do colchete. Deveria sair desta forma: Initial guess V: [10.] Iteration 1, V: [20.], step: 5.0, sigma0_calc: [0.00436189] Iteration 2, V: [25.], step: 2.5, sigma0_calc: [0.02770691] Iteration 3, V: [27.5], step: 1.25, sigma0_calc: [0.0456733] Iteration 4, V: [28.75], step: 0.625, sigma0_calc: [0.05417324] Iteration 5, V: [29.375], step: 0.3125, sigma0_calc: [0.05813452] Final V after iterations: [29.9609375] Wind Speed (m/s): [29.9609375]

Pode me passar um link de como editar o codigo neste formato?

Obrigada, adicionei as tags.

Olá, estou tendo problemas para executar a função cmod5n_inverse, inseri os parâmetros phi(direção do vento-ângulo azimutal), e o ângulo de incidência e sigma0 extraídos do SNAP de uma imagem do Sentinel-1, usando regressão linear. A função é deveria retornar V, a velocidade do vento em metros por segundo, porém os resultados parecem se repetir. Assim: Suposição inicial V: [10. 10. 10. 10. 10.] Iteração 1, V: [20. 20. 20. 20. 20.], passo: 5.0, sigma0_calc: [0.00342216 0.00331692 0.00299704 0.00234962 0.00433662] Iteração 2, V: [25. 25. 25. 25. 25.], passo: 2.5, sigma0_calc: [0.02235289 0.02171421 0.01968306 0.01462111 0.02769178] Iteração 3, V: [27,5 27,5 27,5 27,5 27,5], passo: 1,25, sigma0_calc: [0,03657467 0,03551247 0,0321561 0,0236992 0,04577651] Iteração 4, V: [28,75 28,75 28,75 28,75 28,75], passo: 0,625, sigma0_calc: [0,0431854 0,04192145 0,0379445 0,0279204 0,05434188] Final V após iterações: [28,75 28,75 28,75 28,75 28,75] Velocidade do vento (m/s): [28,75 28,75 28,75 28,75 28,75] Alguém pode me informar qual é o possível erro? Aqui está o seguinte código usado para obter os resultados: sigma0_obs=[4.32674215,3.72549955,3.25660927,2.93046451,4.58418151] incidence=[64.1597756,64.8628612,67.0776741,72.8808313,57.7738178] phi=[-261.571387,-262.019219,-264.074833,-267.473085,-270.674255] sigma0_obs = np.array(sigma0_obs) incidence=np.array(incidence) phi=np.array(phi) import numpy as np import warnings warnings.simplefilter("ignore", RuntimeWarning) def cmod5n_forward(v,phi,theta) : '''! --------- ! cmod5n_forward(v, phi, theta) ! inputs: ! v in [m/s] wind velocity (always >= 0) ! phi in [deg] angle between azimuth and wind direction ! (= D - AZM) ! theta in [deg] incidence angle ! output: ! CMOD5_N NORMALIZED BACKSCATTER (LINEAR) ! ! All inputs must be Numpy arrays of equal sizes ! ! A. STOFFELEN MAY 1991 ECMWF CMOD4 ! A. STOFFELEN, S. DE HAAN DEC 2001 KNMI CMOD5 PROTOTYPE ! H. HERSBACH JUNE 2002 ECMWF COMPLETE REVISION ! J. de Kloe JULI 2003 KNMI, rewritten in fortan90 ! A. Verhoef JAN 2008 KNMI, CMOD5 for neutral winds ! K.F.Dagestad OCT 2011 NERSC, Vectorized Python version !--------------------------------------------------------------------- ''' from numpy import cos, exp, tanh, array DTOR = 57.29577951 THETM = 40. THETHR = 25. ZPOW = 1.6 # NB: 0 added as first element below, to avoid switching from 1-indexing to 0-indexing C = [0, -0.6878, -0.7957, 0.3380, -0.1728, 0.0000, 0.0040, 0.1103, 0.0159, 6.7329, 2.7713, -2.2885, 0.4971, -0.7250, 0.0450, 0.0066, 0.3222, 0.0120, 22.7000, 2.0813, 3.0000, 8.3659, -3.3428, 1.3236, 6.2437, 2.3893, 0.3249, 4.1590, 1.6930] Y0 = C[19] PN = C[20] A = C[19]-(C[19]-1)/C[20] B = 1./(C[20]*(C[19]-1.)**(3-1)) # ! ANGLES FI=phi/DTOR CSFI = cos(FI) CS2FI= 2.00 * CSFI * CSFI - 1.00 X = (theta - THETM) / THETHR XX = X*X # ! B0: FUNCTION OF WIND SPEED AND INCIDENCE ANGLE A0 =C[ 1]+C[ 2]*X+C[ 3]*XX+C[ 4]*X*XX A1 =C[ 5]+C[ 6]*X A2 =C[ 7]+C[ 8]*X GAM=C[ 9]+C[10]*X+C[11]*XX S0 =C[12]+C[13]*X # V is missing! Using V=v as substitute, this is apparently correct V=v S = A2*V S_vec = S.copy() SlS0 = [S_vec<S0] S_vec[SlS0]=S0[SlS0] A3=1./(1.+exp(-S_vec)) SlS0 = (S<S0) A3[SlS0]=A3[SlS0]*(S[SlS0]/S0[SlS0])**( S0[SlS0]*(1.- A3[SlS0])) #A3=A3*(S/S0)**( S0*(1.- A3)) B0=(A3**GAM)*10.**(A0+A1*V) # ! B1: FUNCTION OF WIND SPEED AND INCIDENCE ANGLE B1 = C[15]*V*(0.5+X-tanh(4.*(X+C[16]+C[17]*V))) B1 = C[14]*(1.+X)- B1 B1 = B1/(exp( 0.34*(V-C[18]) )+1.) # ! B2: FUNCTION OF WIND SPEED AND INCIDENCE ANGLE V0 = C[21] + C[22]*X + C[23]*XX D1 = C[24] + C[25]*X + C[26]*XX D2 = C[27] + C[28]*X V2 = (V/V0+1.) V2ltY0 = V2<Y0 V2[V2ltY0] = A+B*(V2[V2ltY0]-1.)**PN B2 = (-D1+D2*V2)*exp(-V2) # ! CMOD5_N: COMBINE THE THREE FOURIER TERMS CMOD5_N = B0*(1.0+B1*CSFI+B2*CS2FI)**ZPOW return CMOD5_N def cmod5n_inverse(sigma0_obs, phi, incidence, iterations=10) : from numpy import ones, array # First guess wind speed V = array([10.]) * ones(sigma0_obs.shape) step = 10. # Debugging prints print(f"Initial guess V: {V}") # Iterating until error is smaller than threshold for iterno in range(1, iterations): sigma0_calc = cmod5n_forward(V, phi, incidence) ind = sigma0_calc - sigma0_obs > 0 V = V + step V[ind] = V[ind] - 2 * step step = step / 2 # Debugging prints print(f"Iteration {iterno}, V: {V}, step: {step}, sigma0_calc: {sigma0_calc}") # Debugging prints print(f"Final V after iterations: {V}") return V # Call cmod5n_inverse to get the wind speed (V) for the given observation V = cmod5n_inverse(sigma0_obs,phi,incidence, iterations=5) # Print the resulting wind speed (V) print("Wind Speed (m/s):", V)