背景

我正在开发一个名为starBliss的自定义星图包。虽然我认为我成功地从我的question中找出了我需要的东西，但看起来有一些边缘情况，事情开始破裂。
有关更多背景信息，请点击此处查看GitHub问题

问题

我目前正在开发一个函数，该函数获取一个星座线shapefile，并基于经度和纬度值使用Lambert Azimuthal Equal Area投影将其重新投影，以创建给定位置在任何给定时间的夜空。
在当前starBliss包中可以看到一个成功的实施：

#devtools::install_github("benyamindsmith/starBliss")
library(ggplot2)
library(starBliss)
p<- plot_starmap(location= "Toronto, ON, Canada",
             date="2022-01-17",
             style="black",
             line1_text="Toronto",
             line2_text ="January 17th, 2023",
             line3_text="43.6532° N, 79.3832° W")
ggsave('toronto_black.png', plot = p, width = unit(10, 'in'), 
       height = unit(15, 'in'))

然而，这种方法确实会遇到一些问题，例如：
（图像被裁剪）

plot_starmap(
  location= "Caracas, Venezuela",
  date = as.Date("1991-03-17"),
  style = "black")

上面的代码创建了一些非对角线（用红色圈出）

一个可重现的示例

我组合了一个函数，它将获取星座线数据并将其转换为当前状态。当我使用该函数并将其绘制在ggplot2中时，问题仍然存在。

library(tidyverse)
library(sf)
library(tidygeocoder)
library(lubridate)  

custom_starmap <- function(location,
                           date){
  
  # Formatting Date properly
  date<- as.Date(date)
  # Formatted date
  dt<- lubridate::ymd(date)
  # Get Latitude and Longitude for ProjString
  
  # For Latitude
  
  suppressMessages(
    capture.output(
      lat <- tibble(singlelineaddress = location) %>%
        geocode(address=singlelineaddress,method = 'arcgis') %>% .[["lat"]]
    )
  )
  
  # Reference date used for calculating longitude
  ref_date <- paste0(year(dt),"01","01",sep="-") %>% ymd()
  # Resulting longitude
  lon <- (-as.numeric(difftime(ref_date,dt, units="days"))/365)*360
  
  
  # The CRS
  projString <- paste0("+proj=laea +x_0=0 +y_0=0 +lon_0=",round(lon,4), " +lat_0=", round(lat,4))
  
  
  # Data Transformation
  flip <- matrix(c(-1, 0, 0, 1), 2, 2)
  
  hemisphere <- st_sfc(st_point(c(lon, lat)), crs = 4326) %>%
    st_buffer(dist = 1e7) %>%
    st_transform(crs = projString)
  
  # Data source for constellation lines
  
  url1 <- "https://raw.githubusercontent.com/benyamindsmith/starBliss/main/data/constellations.lines.json"
  
  # Reading Data
  invisible(
    capture.output(
      constellation_lines_sf <- invisible(st_read(url1, stringsAsFactors = FALSE)) %>%
        st_wrap_dateline(options = c("WRAPDATELINE=YES", "DATELINEOFFSET=360")) %>%
        st_transform(crs = projString) %>%
        st_intersection(hemisphere) %>%
        filter(!is.na(st_is_valid(.))) %>%
        mutate(geometry = geometry * flip)
    )
  )
  
  
  st_crs(constellation_lines_sf) <- projString
  
  return(constellation_lines_sf)
}

# The data

df<-custom_starmap(location= "Caracas, Venezuela",
                   date = as.Date("1991-03-17"))

df

> Simple feature collection with 49 features and 2 fields
Geometry type: GEOMETRY
Dimension:     XY
Bounding box:  xmin: 8967611 ymin: -8898251 xmax: -8714977 ymax: 9004400
CRS:           +proj=laea +x_0=0 +y_0=0 +lon_0=73.9726 +lat_0=10.488
First 10 features:
    id rank                       geometry
1  And    1 MULTILINESTRING ((3542468 3...
2  Ant    3 LINESTRING (-6234955 -52010...
3  Aqr    2 MULTILINESTRING ((8967611 -...
4  Ari    1 LINESTRING (3098546 2071855...
5  Aur    1 MULTILINESTRING ((-1307725 ...
6  Cae    3 LINESTRING (557848.4 -59059...
7  Cam    2 MULTILINESTRING ((-24783.5 ...
8  Cnc    2 MULTILINESTRING ((-6264812 ...
9  CMa    1 MULTILINESTRING ((-2356827 ...
10 CMi    2 LINESTRING (-4432439 -32157..

当我绘制这些数据时，可以看到这些线：

df %>% ggplot()+geom_sf()

（为清楚起见，用红色圈出）

我该如何解决这个问题？我正在使用的CRS格式是否有问题？或者我需要裁剪行吗？