{lattice} subplots appear tilted

Question

I am trying to plot multiple .nc files in R using level plot. I have plotted the data but the plots are titled a bit it seems. I tried with different layouts but the layout is still titled to an angle. I also tried to change the height and width of the plot but still the problem is persisting. Any help with this problem would be appreciated. I have attached my sample plot here. I have tried it with the code below:

library(ncdf4)
library(raster)
library(rgdal)
library(rasterVis)
library(RColorBrewer)
library(colorRamps)
library(gridExtra)
library(terra)
path<- setwd('D:/DATA/2015/MIXING_RATIO/')

# Set the path to the shapefile
shapefile_path <- "C:/Users/IITM/Downloads/India_State/India_State/India_State.shp"

# Load the shapefile
shapefile <- readOGR(shapefile_path)
# Get a list of all netCDF files in the directory
nc_files <- list.files(path, pattern = ".nc", full.names = TRUE)
nc_files
raster_list <- list()
list()
# Loop through each netCDF file
for (nc_file in nc_files) {
  # Read the netCDF data
  data <- nc_open(nc_file)
  
  # Extract the necessary variables
  lon <- ncvar_get(data, "lon")
  lat <- ncvar_get(data, "lat")
  dust_MR <- ncvar_get(data, "DU002")
  
  # Plot the first level slice of the data
  dust1.slice1 <- dust_MR[, , 60]
  
  # Create a raster object
  r <- raster(t(dust1.slice1), xmn = min(lon), xmx = max(lon), ymn = min(lat), ymx = max(lat), crs = CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0"))
  
  # Flip the raster vertically
  r <- flip(r, direction = 'y')
  
  # Crop the raster based on the shapefile extent
  r_crop <- crop(r, extent(shapefile))
  
  # Mask the raster based on the shapefile geometry
  r_mask <- mask(x = r_crop, mask = shapefile)
  
  # Add the raster to the list
  raster_list[[nc_file]] <- r_mask
  
  # Close the netCDF file
  nc_close(data)
}


# Stack the rasters
raster_stack <- stack(raster_list)
crs(raster_stack)
DUST<- stack(mget(rep("raster_stack")))
names(DUST)
raster_names<- gsub('D..DATA.2015.MIXING_RATIO.MERRA2_400.inst3_3d_aer_Nv.|.SUB.nc'," ",names(DUST))
raster_names
raster_dates <- as.Date(trimws(raster_names), format = "%Y%m%d")
raster_dates
txtdates <- format(raster_dates, "%d%b%Y")
txtdates
substitute(paste(bold("txtdates")))
#plot(raster_stack, nc = 4)
library(latticeExtra)
v<- vect(shapefile)
llines.SpatVector <- function(x, ...) {
  xy <- crds(x, list=TRUE)
  names(xy) <- c("x", "y")
  lattice::llines(xy, ...)
}
lns <- as.lines(v)
#line_color <- "white"
# Create a single plot using levelplot
levelplot(raster_stack,layout=c(6,4),names.attr=txtdates,
          margin = FALSE,colorkey = list(space = "right"),
          xlab = substitute(paste(bold('Longitude'))),
          ylab = substitute(paste(bold('Latitude'))),
          main = "Dust Mixing Ratio (bin 2)(Kg/Kg)")+layer(llines(lns))

Thank You for your suggestions.