How to make the result of dates in the loop correlative?

I made a code to calculate daily data from MODIS 8-day data, as it is known the last data of the year is always divided by 5 or 6 and not by 8 depending on whether it is a leap year or not so I added it, the problem is the result because I have repeated dates and other missing ones, for example September 20 is repeated twice but I don't have 21. How could I solve this?

So that the code is well understood, the first thing I did was change the crs of my rasters one by one because I had errors applying it directly to the stack at the end, this worked fine, and then with that list with my reprojected rasters I worked with the loop I have doubts

Here is my code

library(raster)
library(terra)
library(lubridate)

archivos_ET <- list.files(path = "/", pattern = "ET_\\d{4}_\\d{3}", full.names = TRUE)

#Save all raster et files in a list called et_proy

# Create a list to store the daily rasters
et_dia <- list()

for (i in 1:length(et_proy)) {
  # # Load the raster
  raster_actual <- et_proy[[i]]
  # Get the year and day of the current file name
  partes_nombre <- strsplit(basename(archivos_ET[i]), "_")[[1]]
  ano <- as.integer(partes_nombre[2])
  dia <- as.integer(gsub(".tif", "", partes_nombre[3]))
  # Check if the year is a leap year
  es_bisiesto <- ifelse(ano %% 4 == 0 & (ano %% 100 != 0 | ano %% 400 == 0), TRUE, FALSE)
  # Determine the divisor according to whether it is a leap year and the day is 361 (it is the last day of the modis products in a year)
  divisor <- ifelse(dia == 361, ifelse(es_bisiesto, 6, 5), 8)
  # Loop to split the current raster into daily rasters and store them
  for (j in 1:(divisor-1)) {
    # Calcular la fecha del raster diario
    fecha <- as.Date(paste(ano, "-01-01", sep = "")) + dia + j - 1
    # Calculate date from daily raster
    raster_diario <- raster_actual / divisor
    # Assign date to daily raster
    names(raster_diario) <- format(fecha, "%Y-%m-%d")
    # Add the daily raster to the list
    et_dia[[length(et_dia) + 1]] <- raster_diario
  }
  # Calculate the date of the last raster
  fecha_ultimo <- as.Date(paste(ano, "-01-01", sep = "")) + dia + divisor - 2
  # Calculate the last raster
  raster_ultimo <- raster_actual / divisor
  # Assign date to last raster
  names(raster_ultimo) <- format(fecha_ultimo, "%Y-%m-%d")
  # Add the last raster to the list
  et_dia[[length(et_dia) + 1]] <- raster_ultimo
}

And here is an example of my result

[[993]]
class      : RasterLayer 
dimensions : 33, 56, 1848  (nrow, ncol, ncell)
resolution : 346, 463  (x, y)
extent     : 602219.9, 621595.9, 5353158, 5368437  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=18 +south +datum=WGS84 +units=m +no_defs 
source     : memory
names      : X2016.09.20 
values     : 0.6063787, 2.016294  (min, max)


[[994]]
class      : RasterLayer 
dimensions : 33, 56, 1848  (nrow, ncol, ncell)
resolution : 346, 463  (x, y)
extent     : 602219.9, 621595.9, 5353158, 5368437  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=18 +south +datum=WGS84 +units=m +no_defs 
source     : memory
names      : X2016.09.20 
values     : 0.6063787, 2.016294  (min, max)


[[995]]
class      : RasterLayer 
dimensions : 33, 56, 1848  (nrow, ncol, ncell)
resolution : 346, 463  (x, y)
extent     : 602219.9, 621595.9, 5353158, 5368437  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=18 +south +datum=WGS84 +units=m +no_defs 
source     : memory
names      : X2016.09.22 
values     : 1.3125, 2.177313  (min, max)


[[996]]
class      : RasterLayer 
dimensions : 33, 56, 1848  (nrow, ncol, ncell)
resolution : 346, 463  (x, y)
extent     : 602219.9, 621595.9, 5353158, 5368437  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=18 +south +datum=WGS84 +units=m +no_defs 
source     : memory
names      : X2016.09.23 
values     : 1.3125, 2.177313  (min, max)

I also realized that the first data of my series of years, that is, 2014-01-01, is not there either.

Here is a general solution.

Example data at 8-day time step

library(terra)
r <- rast(system.file("ex/logo.tif", package="terra"))   
s <- c(r/c(1:3), r/c(4:6))
time(s) <- as.Date("2001-01-03") + seq(0, 40, 8)
time(s)
#[1] "2001-01-03" "2001-01-11" "2001-01-19" "2001-01-27" "2001-02-04" "2001-02-12"

Add one layer to make the smallest time-step one day

x <- s[[1]]
time(x) <- time(x)-1
x <- c(x, s)

Now use fillTime to add (empty) layers inbetween the 8 days intervals, and remove the added layer.

sx <- fillTime(x)[[-1]]
time(sx)
# [1] "2001-01-03" "2001-01-04" "2001-01-05" "2001-01-06" "2001-01-07" "2001-01-08" "2001-01-09" "2001-01-10" "2001-01-11" "2001-01-12" "2001-01-13"
#[12] "2001-01-14" "2001-01-15" "2001-01-16" "2001-01-17" "2001-01-18" "2001-01-19" "2001-01-20" "2001-01-21" "2001-01-22" "2001-01-23" "2001-01-24"
#[23] "2001-01-25" "2001-01-26" "2001-01-27" "2001-01-28" "2001-01-29" "2001-01-30" "2001-01-31" "2001-02-01" "2001-02-02" "2001-02-03" "2001-02-04"
#[34] "2001-02-05" "2001-02-06" "2001-02-07" "2001-02-08" "2001-02-09" "2001-02-10" "2001-02-11" "2001-02-12" "2001-02-13" "2001-02-14" "2001-02-15"
#[45] "2001-02-16" "2001-02-17" "2001-02-18" "2001-02-19" "2001-02-20" "2001-02-21" "2001-02-22" "2001-02-23" "2001-02-24" "2001-02-25" "2001-02-26"
#[56] "2001-02-27" "2001-02-28" "2001-03-01" "2001-03-02" "2001-03-03" "2001-03-04" "2001-03-05" "2001-03-06" "2001-03-07" "2001-03-08"

Use approximate to estimate missing values by linear interpolation.

a <- approximate(sx)

plot(time(a), minmax(a)[2,])

A more "manual" approach:

stm <- time(s)
tms <- seq(stm[1], stm[length(stm)], 1)
r <- rast(s, nlyr=length(tms))
time(r) <- tms

i <- match(stm, tms)
r[[i]] <- s

And take it from there.

I think I have solved my question, the solution was to have created an initial date and to have reviewed where the name was assigned.

et_dia <- list()
fecha_actual <- as.Date("2013-12-31")

for (i in 1:length(et_proy)) {
  raster_actual <- et_proy[[i]]
  partes_nombre <- strsplit(basename(archivos_ET[i]), "_")[[1]]
  ano <- as.integer(partes_nombre[2])
  dia <- as.integer(substr(partes_nombre[3], 1, 3))
  es_bisiesto <- ifelse(ano %% 4 == 0 & (ano %% 100 != 0 | ano %% 400 == 0), TRUE, FALSE)
  divisor <- ifelse(dia == 361, ifelse(es_bisiesto, 6, 5), 8)
  for (j in 1:(divisor)) {
    fecha_actual <- fecha_actual + 1
    raster_diario <- raster_actual / divisor
    names(raster_diario) <- format(fecha_actual, "%Y-%m-%d")
    et_dia[[length(et_dia) + 1]] <- raster_diario
  }
}